the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Robust 4D climateoptimal flight planning in structured airspace using parallelized simulation on GPUs: ROOST V1.0
Abolfazl Simorgh
Manuel Soler
Daniel GonzálezArribas
Florian Linke
Benjamin Lührs
Maximilian M. Meuser
Simone Dietmüller
Sigrun Matthes
Hiroshi Yamashita
Feijia Yin
Federica Castino
Volker Grewe
Sabine Baumann
The climate impact of nonCO_{2} emissions, which are responsible for twothirds of aviation radiative forcing, highly depends on the atmospheric chemistry and weather conditions. Hence, by planning aircraft trajectories to reroute areas where the nonCO_{2} climate impacts are strongly enhanced, called climatesensitive regions, there is a potential to reduce aviationinduced nonCO_{2} climate effects. Weather forecast is inevitably uncertain, which can lead to unreliable determination of climatesensitive regions and aircraft dynamical behavior and, consequently, inefficient trajectories. In this study, we propose robust climateoptimal aircraft trajectory planning within the currently structured airspace considering uncertainties in standard weather forecasts. The ensemble prediction system is employed to characterize uncertainty in the weather forecast, and climatesensitive regions are quantified using the prototype algorithmic climate change functions. As the optimization problem is constrained by the structure of airspace, it is associated with hybrid decision spaces. To account for discrete and continuous decision variables in an integrated and more efficient manner, the optimization is conducted on the space of probability distributions defined over flight plans instead of directly searching for the optimal profile. A heuristic algorithm based on the augmented random search is employed and implemented on graphics processing units to solve the proposed stochastic optimization computationally fast. An opensource Python library called ROOST (V1.0) is developed based on the aircraft trajectory optimization technique. The effectiveness of our proposed strategy to plan robust climateoptimal trajectories within the structured airspace is analyzed through two scenarios: a scenario with a large contrail climate impact and a scenario with no formation of persistent contrails. It is shown that, for a nighttime flight from Frankfurt to Kyiv, a 55 % reduction in climate impact can be achieved at the expense of a 4 % increase in the operating cost.
The aviation industry has experienced strong growth in recent years (Lee et al., 2021). Air traffic is estimated to grow at a 4.3 % annual rate over the next 20 years (Scherer, 2019). Aviation's contribution to global warming through CO_{2} and nonCO_{2} emissions is currently responsible for 3.5 % of total anthropogenic radiative forcing (RF) (Lee et al., 2021). The nonCO_{2} climate impact includes changes in ozone and methane concentrations induced by nitrogen oxides (NO_{x}), water vapor (H_{2}O), hydrocarbons (HC), carbon monoxide (CO), sulfur oxides (SO_{x}), and increased cloudiness due to persistent contrail formation. Accounting for the growth rate of air traffic, a critical increase in its associated climate impact is foreseen.
Mitigating the climate impact associated with aviationinduced CO_{2} emissions requires progressing along with technical aspects, such as moving toward more efficient aircraft and using novel propulsion as well as sustainable aviation fuel. Technological improvements can only be moderately incorporated into the existing aircraft fleet. This is due to the aircraft's long life service and long phases of development, production, and certification. On the other hand, the climate impact caused by nonCO_{2} emissions, which is about twice as large as that from CO_{2} emissions (Lee et al., 2021), reveals a high spatial and temporal dependency (Grewe et al., 2014b). Such dependencies provide the possibility to mitigate their climate effects through operational strategies, particularly aircraft trajectory optimization to avoid areas sensitive to aircraft emissions, called climatesensitive regions (e.g., Simorgh et al., 2022).
Numerous studies have been proposed to reduce the climate impacts of nonCO_{2} emissions by changing aircraft maneuvers to avoid climatesensitive regions. These studies differ mainly in (1) how the climatesensitive areas are defined and (2) how climatefriendly trajectories are determined. The first attempts to consider climate hotspots were based on areas sensitive to the formation of persistent contrails (see Gierens et al., 2008). In order to provide information on the spatial and temporal dependency of nonCO_{2} effects, climate change functions (CCFs) were developed. These CCFs provide the climate impact of aviation emissions per flown kilometer and per emitted mass of the species as fivedimensional datasets (i.e., longitude, latitude, altitude, time, type of emission) (Matthes et al., 2012; Frömming et al., 2013; Grewe et al., 2014b). Due to their computational complexity, CCFs were not suited for realtime operations. Therefore, the socalled algorithmic climate change functions (aCCFs) were developed. The aCCFs provide a very fast computation of the individual nonCO_{2} climate impact, as they are based on mathematical formulas, which only need relevant local meteorological input parameters (e.g., van Manen and Grewe, 2019). The aCCFs are wellsuited for trajectory optimization tools due to their computational efficiency (Matthes et al., 2017). An enhanced and consistent set (with respect to emission scenario, metrics, etc.) of aCCFs (aCCFV1.0A) has recently been developed and introduced within the EU project FlyATM4E (see Yin et al., 2022; Dietmüller et al., 2022; Matthes et al., 2023).
As for climateoptimal trajectory planning methods, various strategies ranging from mathematical programming (e.g., Campbell et al., 2008) to metaheuristic (e.g., Yin et al., 2018a; Yamashita et al., 2020), indirect optimal control (e.g., Sridhar et al., 2011), and direct optimal control methods (e.g., Niklaß et al., 2017; Lührs et al., 2021, 2016; Matthes et al., 2020) have been adopted. For instance, the direct optimal control approach has been employed by Hartjes et al. (2016) to minimize flight time (or distance flown) in areas sensitive to persistent contrail formation, by Lührs et al. (2021) to minimize average temperature response over the next 20 years (ATR20) associated with nonCO_{2} emissions, and by Vitali et al. (2021) to minimize the global warming potential (GWP) of NO_{x}, H_{2}O, soot, SO_{2}, and contrails. Using aCCFs to quantify climate impacts, Yamashita et al. employed a genetic algorithm to determine climateoptimal aircraft trajectories (Yamashita et al., 2020, 2021). Regarding the optimization methodologies, the mathematical programming methods only apply to simplified aircraft trajectory optimization problems (e.g., in the study conducted by Campbell et al., 2008, the aircraft's dynamic behavior is represented with a linearized model). The metaheuristic methods (e.g., genetic algorithm) require very fast aircraft trajectory prediction in order to find an optimal solution with a large number of iterations; thus, the flight planning problem is usually approximated with a simplified but representative enough problem (e.g., in Yamashita et al., 2020, the optimization is defined with 11 decision variables to characterize lateral path and flight altitude, and the speed profile is considered constant). Finally, with the optimal control methods, the capability to model more accurate aircraft trajectory optimization problems is provided since the problem is represented as a dynamic optimization problem. Nevertheless, there are some drawbacks associated with addressing the formulated problem. The dynamic programming method (as an optimal control approach) results in the “curse of dimensionality” for complex problems (e.g., a full 4D aircraft trajectory optimization problem). Regarding the indirect optimal control approach, deriving analytical solutions using Pontryagin's maximum principle is daunting, especially for problems with singularities (e.g., only a 2D trajectory optimization problem has been addressed in the literature; Sridhar et al., 2011). The direct optimal control approach, despite being very flexible in modeling aircraft trajectory optimization problems (e.g., considering a full 4D dynamical model with nonlinear path and boundary constraints; Lührs et al., 2021), has a high sensitivity to initial conditions, and thus local optimality is its main drawback. In addition, considering the airspace structure with indirect and direct optimal control methods is not straightforward. Interested readers are referred to Simorgh et al. (2022) for our recent comprehensive survey on climateoptimal aircraft trajectory planning, reviewing both the approaches to model climatesensitive regions and trajectory planning methods. A classification of the most recent studies in this field is provided in Table 1.
To quantify the nonCO_{2} climate effects, specific weather variables are required. In the case of aCCFs, variables such as temperature (T), potential vorticity (PV), geopotential (Φ), relative humidity over ice (r_{hum}), and outgoing longwave radiation (OLR) are needed. These variables are obtained from standard weather forecasts. Several factors, including incomplete understanding of the state of the atmosphere, computational complexity, and nonlinear and sometimes chaotic dynamics, affect the accuracy of weather predictions, implying that the weather forecast is inevitably uncertain (WMO, 2012). These weatherforecastrelated uncertainties in the aCCFs and also in aircraft dynamical behavior (e.g., uncertainty in wind and temperature), if not accounted for within aircraft trajectory planning, can lead to inefficient trajectories. Previous research in the field of climateoptimal aircraft trajectory planning has been conducted in a deterministic manner, neglecting the inclusion of any sources of uncertainty (see Table 1) (Simorgh et al., 2022). A first step in managing and integrating meteorological uncertainties into aircraft path planning is to obtain reliable weather forecasts that can predict probable variations in meteorological conditions. To characterize weather forecast uncertainties, probabilistic weather forecasting (PWF) is typically used (AMSCouncil, 2008). Stateoftheart probabilistic weather forecasting is obtained from the ensemble prediction system (EPS), which provides N_{EPS} possible realizations of meteorological conditions called ensemble members (Bauer et al., 2015).
When accounting for the ensemble weather forecast in solving aircraft trajectory optimization, the computational time is an important issue that arises in addition to the capability to consider such uncertainties. This is due to the fact that, instead of taking one weather forecast and solving the trajectory optimization in a deterministic manner, the optimizer should be capable of considering N_{EPS} (e.g., N_{EPS}=50) different forecasts. Several studies have been proposed in the literature to determine robust aircraft trajectories in the presence of meteorological uncertainty quantified using EPS weather forecast (though not considering climate impact; Simorgh et al., 2022). However, these studies suffer mainly from computational perspectives (i.e., the computational time of the optimizer when considering weather uncertainty) and some restrictive assumptions (e.g., inaccurate modeling of the aircraft dynamical model or ignoring important decision variables such as the flight altitude) (Simorgh et al., 2022, Sect. 5.3). For instance, in the study conducted by GonzálezArribas et al. (2018), a robust aircraft trajectory optimization problem is solved using the optimal control approach within fully freerouting airspace. In this approach, the effects of uncertainty are included in the optimization problem by expanding the dynamical model of the aircraft (almost) linearly to the number of ensemble members, resulting in a larger dimensional deterministic optimization problem. Thus, it requires more computational time compared to the deterministic flight planning problem. Also, the flight planning problem is limited to optimizing only the lateral path. As for the structured airspace, Franco Espín et al. (2018) proposed uncertain flight plan optimization using mixedinteger linear programming (MILP). Similarly, the optimization in this study is performed in 2D airspace. In addition, the fuel burn and its associated nonlinearities are ignored, and it has a sizable computational cost of several minutes. In this respect, developing efficient trajectory optimization solvers capable of delivering robust climateoptimal trajectories with a computational time compatible with operations has been identified as a scientific gap (see Simorgh et al., 2022).
The focus of recent studies has been restricted to planning climateoptimal trajectories considering the concept of future freeroute airspace (see last column Table 1), which is thus not applicable for the structured airspace of today. The trajectory optimization problem constrained by the structure of airspace results in hybrid decision spaces (e.g., route and flight level are discrete, and speed schedule is continuous) (GonzálezArribas et al., 2023). The trajectory optimization problem with the combination of discrete and continuous decision variables is one of the most challenging optimization problems, typically solved using mixedinteger nonlinear programming with intensive computational cost (see, e.g., Bonami et al., 2013).
Soler et al. (2014)Grewe et al. (2014a)Hartjes et al. (2016)Lührs et al. (2016)Lim et al. (2017)Matthes et al. (2017)Niklaß et al. (2017)Yin et al. (2018b)Yin et al. (2018a)Niklaß et al. (2019)Yin et al. (2022)Yamashita et al. (2020)Matthes et al. (2020)Lührs et al. (2021)Yamashita et al. (2021)Drawing upon the brief literature review and the presented open problems, we aim to address the problem of determining robust climateoptimal aircraft trajectories within the structured airspace in this study. Our main contributions are summarized as follows: (1) full 4D climateoptimal trajectory planning within the currently structured airspace, (2) accounting for uncertain meteorological conditions and uncertainty associated with initial flight conditions such as departure time and aircraft initial mass, and (3) determining the optimized trajectory computationally very fast. The uncertainty in weather forecast is characterized using the ensemble prediction system, and aviation's climate impacts are quantified by employing the latest version of aCCFs (V1.0A). The concept of robustness that we refer to is the determination of the aircraft trajectory considering all possible realizations of meteorological variables provided using an EPS weather forecast. In other words, instead of planning a trajectory based on one forecast in a deterministic manner, we aim to determine a trajectory that is optimal considering the overall performance obtained from using different members of an ensemble weather forecast. In this respect, from the operational point of view, the optimized trajectory is tracked as determined, and the effects of meteorological uncertainties are reflected in the flight performance variables such as flight time, fuel burn, and climate impact. Mathematically, the perturbations due to the meteorological uncertainty are considered in the dynamical model of aircraft, and the proposed trajectory optimization is generic in terms of the objective function, in which a wide range of objectives, such as flight time, fuel consumption, emissions, and climate impact, and different statistics including expected values and variance of the performance variables, can be considered. Such flexibility in defining the cost function allows for solving a multiobjective optimization problem. Moreover, by penalizing the mean and variance of the objectives, the effects of uncertainty on flight variables can be controlled. In this study, the flight planning objective is a weighted sum of the simple operating cost (as a function of flight time and fuel consumption) and climate impact. The focus is restricted to optimizing the expected performance since, as will be shown in the simulation results, minimizing the averaged performance leads to reducing the uncertainty ranges for the considered case studies.
We employ the probabilistic flight planning method firstly developed by GonzálezArribas et al. (2023) to determine robust climateoptimal trajectories for three phases: climb, cruise, and descent. In this approach, to account for discrete and continuous decision variables in an integrated manner, the optimization is carried out on the space of probability distributions defined over flight plans instead of directly searching for the optimal profile. Then, the probability distribution over flight plans is parameterized, allowing the generation of multiple flight plans stochastically. The augmented random search algorithm is employed and implemented on graphics processing units (GPUs) to deliver a nearly optimal solution to the resulting stochastic optimization in seconds. We have developed an opensource Python library called ROOST V1.0 (Robust Optimization of Structured Trajectories) based on the proposed robust aircraft trajectory optimization technique, which is currently available via the DOI https://doi.org/10.5281/zenodo.7121862 (Simorgh, 2022). ROOST is a tool that efficiently uses the information provided by the prototype aCCFs (implemented in our recently developed Python library CLIMaCCF V1.0 (DOI: https://doi.org/10.5281/zenodo.6977273; Dietmüller et al., 2022) for planning climateoptimal trajectories accounting for the operational constraints and uncertainty. It should be noted that CLIMaCCF V1.0 is the first release of the Python library, which includes both versions of aCCFs, i.e., V1.0 and V1.0A. For performing aircraft trajectory optimization in this study, we use V1.0A aCCFs. Users need to input the required weather variables, route graphs, and aircraft type specifications to start working with the library. In addition, users should assign values to the weighting parameters associated with different objectives in the objective functions. This paper is mainly devoted to the climateoptimal aircraft trajectory planning algorithm implemented in the library and its application to optimize different case studies. Instructions to get started with the library can be found in the repository of ROOST.
The paper is arranged as follows. The robust climateoptimal aircraft trajectory planning problem is stated and formulated in Sect. 2 and solved in Sect. 3. The potential of our flight planning algorithm to reduce aviationinduced climate impacts under uncertain meteorological conditions is explored in Sect. 4. Finally, the discussion on the obtained results and concluding remarks are presented in Sect. 5.
Aviationinduced nonCO_{2} climate effects have a direct dependency on atmospheric conditions at the location and time of emissions. In this respect, generating flight plans that could potentially minimize emissions in areas with high sensitivity to aircraft emissions in terms of climate change is a step towards climatefriendly air transport. To enable such a flight planning strategy, information on climatesensitive areas is required first. Once climate information is available, it needs to be included in flight planning tools to generate ecoefficient trajectories. It is worth mentioning that to generate reliable trajectories, different sources of uncertainty need to be identified, understood, and quantified.
The focus of this study is on determining robust flight plans taking into account meteorological uncertainty and uncertainty associated with the initial flight conditions. In this respect, we start by presenting a general formulation of the dynamic optimization problem with uncertainty in Sect. 2.1, which is then employed to model the proposed robust flight planning (in Sect. 2.2). In Sect. 2.3, the effects of meteorological uncertainty on the efficiency of climateoptimal trajectories are discussed, and the motivation to solve robust trajectory optimization is provided.
2.1 General formulation of dynamic optimization problem with uncertainty
We present a general formulation of a dynamic optimization problem with the inclusion of uncertainty effects, focusing mainly on the dynamical model and objective function.
Let us consider a class of dynamical systems with uncertainty as follows:
where $\mathit{u}\in {\mathbb{R}}^{{n}_{u}}$, $\mathit{x}\in {\mathbb{R}}^{{n}_{x}}$, and $\mathit{z}\in {\mathbb{R}}^{{n}_{z}}$ are the vectors of control inputs as well as states and algebraic variables, respectively, and f is a vector field, mapping $\mathbb{R}\times {\mathbb{R}}^{{n}_{z}}\times {\mathbb{R}}^{{n}_{x}}\times {\mathbb{R}}^{{n}_{u}}\times {\mathbb{R}}^{{n}_{\mathit{\zeta}}}\to {\mathbb{R}}^{{n}_{x}}$. The uncertain parameters (denoted with $\mathit{\zeta}\in {\mathbb{R}}^{{n}_{\mathit{\zeta}}}$) are considered continuous random variables and assumed to have known probability distribution functions $\mathit{\zeta}(\cdot ):\mathrm{\Omega}\to {\mathbb{R}}^{{n}_{\mathit{\zeta}}}$, where Ω is a space of possible outcomes. The random variables take different values depending on the probable outcomes (i.e., ζ(ω) for ω∈Ω). The nonlinear function f(⋅) is assumed to be a measurable function in ζ. To emphasize the effects of random variables on the system's trajectories, all the variables were denoted with dependency on possible outcomes (e.g., x(t,ω)). For the sake of improving clarity, the abbreviated notation (e.g., x(t)) will be used in the following.
A general form of the cost functional considered for dynamic optimization problems with uncertainty is
where $\mathcal{M}:\mathbb{R}\times {\mathbb{R}}^{{n}_{x}}\times \mathbb{R}\times {\mathbb{R}}^{{n}_{x}}\to \mathbb{R}$ and $\mathcal{L}:\mathbb{R}\times {\mathbb{R}}^{{n}_{x}}\times {\mathbb{R}}^{{n}_{u}}\times {\mathbb{R}}^{{n}_{z}}\times {\mathbb{R}}^{{n}_{\mathit{\zeta}}}\to \mathbb{R}$ are the Mayer and Lagrangian terms, respectively. Notice that the objective function and also dynamical model are similar to what is usually considered within the context of optimal control theory. Although the objective function is written in terms of using the expectation operator, other statistics can be evaluated under this formulation. For instance, one can include the variance of a function A(ζ) as $\mathbb{V}\mathit{\left\{}A\mathit{\right\}}=\mathbb{E}\mathit{\left\{}\right(A\mathbb{E}\mathit{\left\{}A\mathit{\right\}}{)}^{\mathrm{2}}\mathit{\}}=\mathbb{E}\mathit{\{}{A}^{\mathrm{2}}\mathit{\}}\mathbb{E}\mathit{\{}A{\mathit{\}}}^{\mathrm{2}}$, where 𝕍{⋅} is the variance operator. Depending on the benchmark problem, several constraints such as equality and inequality path and boundary constraints are considered (see Simorgh et al., 2022).
The objective of the stated dynamic optimization problem is to find a feasible control policy (u∈𝒰) to minimize the performance index (Eq. 2) respecting a set of constraints, including dynamical constraints (Eq. 1). Depending on the benchmark problem, reformulations and approximations are normally made to address the required performance, such as computational complexity. In this study, we will slightly reformulate the optimization for the proposed path planning problem.
2.2 Robust climateoptimal flight planning problem formulation
The definition of the aircraft trajectory optimization problem mainly requires the aircraft dynamical model, flight planning objectives, and physical and operational constraints. To consider climate impact within aircraft trajectory planning, information on the climate impacts of CO_{2} and nonCO_{2} emissions is necessary and needs to be included in the objective function. In the following, the modeling of the mentioned elements is briefly presented.
2.2.1 Aircraft dynamical model
To determine reliable aircraft trajectories, accurate aircraft dynamical models are necessary. In this work, the pointmass model with the following equations of motion is used to represent the aircraft's dynamical behavior, as is usually considered within air traffic management studies.
Here, ϕ is the latitude, λ is the longitude, v is the true airspeed, h is the altitude, m is the mass, C_{T} is the thrust coefficient, γ is the climb angle, χ is the heading angle, ${C}_{L}\left(\mathit{\gamma}\right)=\left(\mathrm{2}mg\mathrm{cos}\mathit{\gamma}\right)/\left(\mathit{\rho}{v}^{\mathrm{2}}S\right)$, and (w_{x},w_{y}) are wind components. The Earth's ellipsoid radii of curvature in the meridian and the prime vertical are denoted with R_{M} and R_{N}, respectively, T is the magnitude thrust force, and the drag force is denoted with D. g is the Earth's gravity, FF is the fuel flow, and S is the wetted surface of the aircraft. The BADA4.2 model (Gallo et al., 2006) is employed to provide the aerodynamic and propulsive performance of the aircraft. Notice that we used the notation presented in Sect. 2.1 (i.e., ω) to represent the uncertainty in the temperature and components of wind. The realization of the state and control variables depends on uncertainty, e.g., for mass, $m=m(t,\mathit{\omega})$; in the interest of notational clarity, we use the abbreviated notation m.
Structured airspace
As trajectory optimization in this study is performed within the structured airspace, the evolutions of aircraft states are constrained. In the following, we briefly present our proposed modeling of airspace structure and flight plan.
We consider a directed acyclic graph $G=(V,E)$ to model the airspace with V as the navigation waypoints and e∈E as the airway edges connecting waypoints. The trajectory is assumed to start at the end of the standard instrument departure procedure (SID) and end at the beginning of the standard instrument arrival procedure (STAR) to the destination airport, denoted as o∈V and d∈V, respectively. We define the flight plan ℱ with a tuple $(\mathcal{R},\stackrel{\mathrm{\u203e}}{\text{FL}},\stackrel{\mathrm{\u203e}}{\text{M}},\mathcal{C},\mathcal{D},{d}_{\mathcal{D}})$. In the flight plan (ℱ), the route (or lateral path) denoted as ℛ includes a sequence of waypoints, i.e., $\mathcal{R}:=({r}_{\mathrm{0}},{r}_{\mathrm{1}},\mathrm{\cdots},{r}_{{n}_{r}})$. The vertical profile of the cruise, $\stackrel{\mathrm{\u203e}}{\text{FL}}$, is composed of an ordered sequence of tuples of the form (r_{k},FL_{k}), indicating that, if the aircraft is in the cruise phase, the flight level will be changed to FL_{k} when the waypoint r_{k} is reached (see Fig. 1). The Mach schedule $\stackrel{\mathrm{\u203e}}{\text{M}}:=({\text{M}}_{\mathrm{0}},{\text{M}}_{\mathrm{1}},\mathrm{\cdots},{\text{M}}_{{n}_{r}})$ indicates the target Mach number M_{k} at waypoint r_{k} during the cruise phase. Climb and descent profiles $\mathcal{C},\mathcal{D}:\mathbb{R}\to \mathbb{R}$ are represented by continuous and piecewisedifferentiable functions mapping the altitude to the target airspeed during the climb and descent phases, respectively. Finally, a scalar variable d_{𝒟} shows the distance to go to the destination node at which the aircraft should end the cruise and start the descent phase.
2.2.2 Objective function
The goals of the aircraft trajectory optimization problem are interpreted mathematically and defined as an objective function to be minimized (or maximized). In addition to the climate impact, the operating cost is a crucial aspect that needs to be considered as it is one of the main interests of airliners. Generally, there is a tradeoff between the operating cost and climate impacts. This is due to the fact that rerouting areas sensitive to climate increases operational costs as the aircraft tends to fly longer routes (Niklaß et al., 2021). To consider both objectives in the trajectory optimization, we define the following objective function:
where ψ_{CST} and ψ_{CLM} are weighting parameters penalizing cost and climate impact, respectively. In the following, the proposed modeling of these two objectives is presented.
Operating cost
There are various approaches to account for operating costs within aircraft path planning. Flight time and/or fuel are common objectives. However, more realistic cost metrics exist, which include additional costs such as flight crew, cabin crew, and landing fees. Interested readers are referred to Simorgh et al. (2022) for a classification of these metrics.
In this study, we use simple operating cost (SOC) as a metric expressing cost in USD with linear relation to flight time and fuel consumption:
where t_{0} and t_{f} are the initial flight time and final flight time weighted by ψ_{t}=0.75 USD s^{−1}, and m_{0} and m_{f} are the initial mass and final mass weighted by ψ_{m}=0.51 [USD kg^{−1}]. The expectation operator is used here because the meteorological uncertainty (included in the aircraft dynamical model) and also uncertainty due to initial conditions directly affect the flight time and fuel consumption.
In spite of considering only flight time and fuel consumption to represent the operating cost, it was reported in Table 4 of Yamashita et al. (2020) that employing SOC and a more comprehensive metric such as cash operating cost within trajectory optimization delivered almost similar results. This is because these two metrics, in the end, consider flight time and fuel consumption as inputs to estimate the operating cost.
Climate effects induced by aviation
Numerous approaches have been proposed in the literature to consider climate impact within aircraft trajectory planning strategies (see Table 1 and Fig. 2 of Simorgh et al., 2022). In this study, the climate impact of aircraft operations is modeled using the prototype aCCFs. These aCCFs take as input specific meteorological variables (e.g., temperature and relative humidity) and estimate the climate impact in terms of average temperature change. The suitability of aCCFs for climateoptimal trajectory planning can be justified as follows:

aCCFs account for the temporal and spatial dependency of climate impacts associated with nonCO_{2} species, including ozone and methane, resulting from NO_{x} emissions, water vapor emissions, and persistent contrails;

aCCFs estimate the climate impact associated with aircraft emissions computationally in real time, making it wellsuited for climateoptimal trajectory planning; and

aCCFs directly quantify climate impacts in average temperature change.
The aCCF V1.0 reported by Yin et al. (2022) is the first complete and consistent set of prototype aCCFs providing spatially and temporally dependent nonCO_{2} climate effects in terms of average temperature response over the next 20 years for a pulse emission scenario (PATR20). In this study, we use the latest version of aCCFs, i.e., aCCF V1.0A developed by Matthes et al. (2023) within the EU project FlyATM4E, which has been calibrated to the stateoftheart climate response model AirClim (Dahlmann et al., 2016).
The V1.0A aCCFs are multiplied by some factors in order to be more suitable for the flight planning application.

Metric conversion. The selection of a suitable metric depends on the question to be answered (see Grewe and Dahlmann, 2015, for more details); therefore, different questions require the use of different metrics. The PATR20 metric was selected as a metric for the aCCFs to assess the impact of a “simple” pulse emission. However, factors are available (see Dietmüller et al., 2022) to convert PATR20 to other metrics: for example, assuming the future emission scenario or longer time horizons. This way, one can select the emission scenario and time horizon that are bestsuited for the question. In this study, the FATR20 metric is used to assess the climate effect reduction obtained by steadily applying a specific routing strategy under the assumption of a future businessasusual emission scenario. To this end, the aCCFs based on the pulse emission scenario are converted to the future scenario (FATR20) using values reported in Table 3 of Dietmüller et al. (2022). These conversion factors were derived by simulations with the climate response model AirClim (Dahlmann et al., 2016): one simulation with pulse emission and one with the future emission scenario. By comparing the two simulations, the factors can be derived.

Efficacy. Efficacies were introduced to take into account the fact that the radiative forcing of some nonCO_{2} forcing agents (e.g., ozone, methane, contrails) is more or less effective in changing the global mean nearsurface temperature per unit forcing compared to the response of CO_{2} forcing (Hansen et al., 2005; Ponater et al., 2007). Dietmüller et al. (2022) have summarized the efficacy parameters reported by Lee et al. (2021). For a detailed explanation of efficacy, the reader is referred to the stateoftheart literature (e.g., Ponater et al., 2007; Rap et al., 2010; Bickel et al., 2020).
For the sake of compactness of notation, FATR20 (from aCCF V1.0A) with efficacy parameters is replaced with ATR in the following.
For the aCCFs of (daytime and nighttime) contrails, the ice supersaturation is applied using temperature and relative humidity over ice in order to predict regions where persistent contrails are expected to form, called persistent contrail formation areas (PCFAs) (Schmidt, 1941; Appleman, 1953). To represent the climate impacts using aCCFs in the average temperature change (i.e., [K]), fuel consumption rate, NO_{x} emissions, and distance flown through PCFAs are required.
The geographical aCCF pattern of water vapor, NO_{x}induced ozone (production), and methane (destruction), as well as of contrails, is shown in Fig. 2a–c for 13 June 2018 at 00:00 UTC over the European region for FL340. Moreover, Fig. 2d shows the pattern of the merged nonCO_{2} aCCFs that combines the individual aCCFs (Fig. 2a–c). Note that to generate the merged aCCFs and to compare the contribution of each species, we adopt typical transatlantic fleet mean values to unify the units of aCCFs in K kg^{−1}(fuel). The approximated conversion factors for NO_{x} emissions and contrails are $\mathrm{13}\times {\mathrm{10}}^{\mathrm{3}}$ kg(NO_{2}) kg^{−1}(fuel) and 0.16 km kg^{−1}(fuel), respectively (Graver and Rutherford, 2018; Penner et al., 1999). It is clear from the merged aCCFs that the contrails have dominant climate effects, which is in line with related studies employing aCCFs (see, e.g., Dietmüller et al., 2022).
To benefit from the spatial and temporal dependency of nonCO_{2} climate effects identified using aCCFs in planning climateoptimal trajectories, we define the following objective expressed in the Lagrangian form for Eq. (4).
for $i\in \mathit{\{}{\text{CH}}_{\mathrm{4}},\text{Cont.},{\text{O}}_{\mathrm{3}},{\text{H}}_{\mathrm{2}}\text{O},{\text{CO}}_{\mathrm{2}}\mathit{\}}$:
Here, $\mathit{\zeta}(\cdot ):=\left[{w}_{x}\right(\cdot \left)\phantom{\rule{0.25em}{0ex}}{w}_{y}\right(\cdot \left)\phantom{\rule{0.25em}{0ex}}T\right(\cdot \left)\phantom{\rule{0.25em}{0ex}}\text{OLR}\right(\cdot \left)\phantom{\rule{0.25em}{0ex}}\text{PV}\right(\cdot \left)\phantom{\rule{0.25em}{0ex}}\mathrm{\Phi}\right(\cdot )\phantom{\rule{0.25em}{0ex}}{\text{r}}_{\mathrm{hum}}$ ${\text{F}}_{\mathrm{in}}(\cdot )\phantom{\rule{0.25em}{0ex}}\text{q}(\cdot )\phantom{\rule{0.25em}{0ex}}{t}_{\mathrm{0}}(\cdot )\phantom{\rule{0.25em}{0ex}}{m}_{\mathrm{0}}(\cdot )]$ is a vector of uncertain variables, i.e., meteorological variables and initial flight conditions. In addition, ${\dot{m}}_{{\mathrm{no}}_{x}}\left(t\right)=\text{FF}(t,\mathit{u})\cdot {\text{EI}}_{{\mathrm{NO}}_{x}}(t,\mathit{x},\mathit{u})$, where ${\text{EI}}_{{\mathrm{NO}}_{x}}(\cdot )$ is the actual NO_{x} emission index in [g(NO_{2}) kg^{−1}(fuel)], and v_{gs} is the ground speed. The NO_{x} emission index varies with many factors such as aircraft type, fuel flow, flight altitude, and synoptical situation. To consider such dependencies, the Boeing fuel flow method 2 (BFFM2), calculating the actual emission index of NO_{x} from the reference conditions, is adopted (DuBois and Paynter, 2006; Jelinek, 2004). The reference conditions are obtained from the International Civil Aviation Organization (ICAO) data bank (https://www.easa.europa.eu/en/domains/environment/icaoaircraftengineemissionsdatabank, last access: 2 July 2023). Notice that the objective regarding climate impact is expressed in the Lagrangian form (i.e., as an integral) since the climate impact needs to be evaluated (and accumulated) along the route, unlike the operating cost, which requires only information on boundary values (e.g., initial and final flight time). Notice that we define the objectives considering the expected performance. One can use other statistics, such as variance, without loss of generality. For instance, to penalize the range of uncertainty in climate impacts, $\mathbb{V}\mathit{\left\{}\text{ATR}\mathit{\right\}}=\mathbb{E}\mathit{\left\{}{\text{ATR}}^{\mathrm{2}}\mathit{\right\}}\mathbb{E}\mathit{\{}\text{ATR}{\mathit{\}}}^{\mathrm{2}}$ should be included in the objective function. For the considered case studies, it will be shown that minimizing the expected objective function will reduce the uncertainty in the most uncertain variable for the climateoptimal routing option, thus providing a robust solution without penalizing the deviation.
2.3 Effects of weatherinduced uncertainty on climate
The dynamical model of aircraft requires weatherrelated variables such as wind and temperature (see, e.g., Sect. 2.2.1). In addition, the nonCO_{2} climate effects of aviation included in the objective function of the optimization problem strongly rely on weather conditions. Since the required meteorological variables are obtained from standard weather forecasts, they are inevitably uncertain. It is worth mentioning that there are other sources of uncertainty affecting the efficiency of the planned climateoptimized trajectories, including uncertainty from climate science (e.g., modeling and estimating aviationinduced climate effects), emission calculation, and also inaccurate modeling of aircraft behavior, which are not within the scope of this paper but have been identified as open problems (see Matthes et al., 2023).
In this paper, the focus is on forecastrelated uncertainties, which will be characterized by employing ensemble prediction system (EPS) weather forecasts, a numerical weather prediction method introduced to deal with uncertainty in weather forecast (Bauer et al., 2015). These are forecasts in which both the initial conditions and the physical parameters of a numerical weather integration model are slightly modified from one member to another and provide N_{EPS} (typically N_{EPS} = 50) different predictions known as ensemble members (Bauer et al., 2015). Each member of an ensemble represents one possible realization of meteorological situations. As the aCCFs take as inputs some meteorological variables, N_{EPS} different aCCFs can be calculated for an EPS weather forecast. For instance, the meteorological variables temperature and relative humidity over ice are required for the aCCF of nighttime contrails. Notice that relative humidity over ice is required for identifying icesupersaturated areas. Feeding N_{EPS} probable realizations of these meteorological variables (i.e., ensemble members), N_{EPS} different aCCFs (i.e., ${\text{aCCF}}_{{\mathrm{Cont}}_{i}}$ for $i=\mathrm{1},\mathrm{\cdots},{N}_{\mathrm{EPS}}$) are calculated. The same applies to system dynamics (considering uncertainty in temperature and wind) and also the NO_{x} emission index (due to the dependency on ambient temperature and specific humidity).
To investigate the degree of uncertainty (or variability) in the meteorological variables provided by an EPS and its effects on the computed aCCFs, we take the standard deviation (SD) from 10 ensemble members of weather data obtained using the ERA5 reanalysis data products (https://cds.climate.copernicus.eu/, last access: 2 July 2023). It should be noted that the reanalysis data products are generated from postprocessing with observations. Thus, the variability among the ensemble members is expected to be lower than the forecast data with more ensemble members yet still valid to illustrate. Figure 3 shows the SD of weather variables required to calculate the aCCFs and aircraft trajectory on 13 June 2018 at 00:00 UTC for FL340. The SD is taken over the normalized variables (with respect to their maximum values) for comparison purposes. The variability of geopotential, temperature, and wind is small compared to potential vorticity, outgoing longwave radiation, and relative humidity. The SDs of the calculated aCCFs based on the ensemble members are illustrated in Fig. 4. Since the aCCF of NO_{x} emissions (i.e., methane and ozone) depends on geopotential and temperature, its SD is small compared to the aCCFs of water vapor and nighttime contrails, which are based on potential vorticity and relative humidity (when applying the icesupersaturated condition), respectively. Notice that the uncertainty in the climate impact of contrails is much higher than water vapor due to the variability of relative humidity in satisfying the persistency condition of contrails (see SD of PCFA in Fig. 4). Due to the considered time (i.e., 00:00 UTC), the aCCF of contrails is based on the formulation of nighttime aCCFs. In spite of negligible uncertainty in the aCCF of NO_{x} and also relatively low uncertainty in the aCCF of water vapor compared to aCCF of contrails, due to the dominant climate impact of contrails, the net nonCO_{2} climate effect is considerably uncertain (see SD of the merged aCCFs in Fig. 4), which must be crucially taken into consideration.
Figure 5 shows how uncertainty in meteorological variables can affect the performance of aircraft trajectories. As can be seen, these uncertainties are accumulated and can considerably degrade the efficiency of the optimized trajectory if not considered in the aircraft trajectory planning a priori. No recent study on the determination of climateoptimized trajectories has considered robustness in the sense of uncertainty in weather forecasts. One of the main reasons for not considering such variations is the computational time that arises within the optimization techniques. In fact, instead of considering one member, the optimization, in this case, is to consider N_{EPS} ensemble members (e.g., 50) and find an optimized trajectory to be optimal with respect to all probable realizations of meteorological variables. Such an increase in dimensions is daunting to cope with by employing classical dynamic optimization approaches such as direct optimal control.
In Sect. 3, we will address the problem of robust climateoptimal trajectory planning with an efficient heuristic method implemented on GPUs.
The aircraft trajectory optimization problem formulated in Sect. 2.2 is solved by employing the method firstly developed by GonzálezArribas et al. (2023), which is a heuristic optimization technique for the structured airspace and is capable of determining the optimized trajectory in four dimensions, i.e., latitude, longitude, altitude, and time.
As mentioned in Sect. 2.1, depending on the problem, some reformulations and approximations are normally made to the dynamic optimization problem to address the required objectives. Here, instead of seeking the optimal control policy (i.e., u^{o}), the goal is to find an optimal flight plan ℱ^{o}, i.e., $({\mathcal{R}}^{\mathrm{o}},{\stackrel{\mathrm{\u203e}}{\text{FL}}}^{\mathrm{o}},{\stackrel{\mathrm{\u203e}}{\text{M}}}^{\mathrm{o}},{\mathcal{C}}^{\mathrm{o}},{\mathcal{D}}^{\mathrm{o}},{d}_{\mathcal{D}}^{\mathrm{o}})$ (see Sect. 2.2.1) that minimizes the objective function given in Eq. (4) and satisfies the aircraft dynamical model (given in Eq. 3), airspace structure, and path and boundary constraints. Such a selection of decision space allows us to directly account for the operational restrictions, such as determining lateral routes that follow the airspace structure.
3.1 Ensemble trajectory integration
To determine the performance of a flight plan and evaluate the cost function Eq. (4), the corresponding trajectories of the aircraft are to be calculated using the aircraft dynamical model (provided in Sect. 2.2.1). As mentioned in Sect. 2.3, aircraft trajectories are affected by uncertainty in meteorological variables, including temperature and wind. Assuming a unique lateral path (constant course) to be tracked in practice with lowlevel controllers in real time, the uncertainty in the wind (both magnitude and direction) will affect ground speed and, consequently, the time the aircraft flies the route (see Fig. 6). In addition, uncertainty in temperature affects fuel consumption because the propulsive and aerodynamic performance of the aircraft and also airspeed have a dependency on temperature (GonzálezArribas et al., 2023). From Eq. (6), one can conclude that the uncertainty in flight time and flight mass can also affect the climate impacts (see also Fig. 5). In addition to the uncertainty in meteorological variables, which is characterized using ensemble weather forecasts, uncertainty associated with initial flight conditions is taken into account within the trajectory optimization problem. In this study, the initial flight time and initial mass of the aircraft are modeled as Gaussian variables, i.e., ${t}_{\mathrm{0}}\sim \mathcal{N}({\overline{t}}_{\mathrm{0}},{\mathit{\sigma}}_{{t}_{\mathrm{0}}})$ and ${m}_{\mathrm{0}}\sim \mathcal{N}({\overline{m}}_{\mathrm{0}},{\mathit{\sigma}}_{{m}_{\mathrm{0}}})$.
To efficiently reflect the effects of wind uncertainty in flight performance variables, instead of time, the distance flown along the route (s) is considered to be the independent variable using $\left(\text{d}t\right)(\text{d}s{)}^{\mathrm{1}}={v}_{\mathrm{gs}}^{\mathrm{1}}$. This is beneficial since the uniqueness of time for all possible realizations of wind means that the position of the aircraft is fixed with respect to time. In this case, the effects of uncertainty cannot be considered efficiently because the range of feasible solutions is limited. The selection of distance flown as the independent variable allows for reflecting wind uncertainty in flight time.
According to the defined objective function (Eqs. 4, 5, 6), the flight performance variables required to evaluate Eq. (4) are flight time, final mass, and climate impact. By using TI(⋅) to denote the integration of the aircraft dynamical model for a given flight plan, weather data, and initial flight conditions, we receive the expected final mass, the expected final time, and the expected ATR required to evaluate the performance of aircraft trajectory (i.e., calculate the objective function given in Eq. 4) as follows:
In the formulated robust optimization problem, the uncertain variables are considered to be continuous random variables with known probability distribution functions. As the weather variables for an ensemble weather forecast are directly represented in a discrete fashion, we assume a discrete probability distribution for the uncertainty. Considering N_{EPS} probable realizations of uncertainty (i.e., $\mathit{\{}{\mathit{\zeta}}_{j}{\mathit{\}}}_{\mathrm{1}}^{{N}_{\mathrm{EPS}}}$), we have
where 𝒲^{j}(⋅) is a set of meteorological variables required for climateoptimal aircraft trajectory planning corresponding to the jth member of an EPS weather forecast, and ${t}_{\mathrm{0}}^{j},{m}_{\mathrm{0}}^{j}\sim {t}_{\mathrm{0}}$, and m_{0} are initial flight conditions that are sampled independently. Generally, different members of the ensemble weather forecasts are considered equally probable. This implies that a specific forecasted weather pattern that has a higher probability will be represented by a larger number of ensemble members. In this study, we use equal weights for each probable realization of weather variables, i.e., with a probability of $\mathbf{P}(\mathcal{W}={\mathcal{W}}^{j})={N}_{\mathrm{EPS}}^{\mathrm{1}}$. Thus, using an unweighted average between all ensemble members, Eq. (8) can be written as
The expected ATR is calculated as
for $i\in \mathit{\{}{\text{CH}}_{\mathrm{4}},\text{Cont.},{\text{O}}_{\mathrm{3}},{\text{H}}_{\mathrm{2}}\text{O},{\text{CO}}_{\mathrm{2}}\mathit{\}}$. For instance, ATR^{j} for ozone and contrails can be calculated as
where x^{j}(t^{j}) and u^{j}(t^{j}) are the state and control variables of the aircraft considering the jth realization of weather variables and the jth sampled initial conditions, $\text{d}s={v}_{\mathrm{gs}}^{j}\cdot \text{d}{t}^{j}$, and
where weather variables such as T^{j} and Φ^{j} are the jth members of the ensemble weather forecast. The actual NO_{x} emission index, i.e., ${\text{EI}}_{{\mathrm{NO}}_{x}}^{j}\left({\mathit{x}}^{j}\right({t}^{j}),{\mathit{u}}^{j}({t}^{j}\left)\right)$, is calculated using BFFM2. As can be seen in Eq. (11), the climate impact due to the NO_{x} emissions depends on the amount of NO_{x} emitted in NO_{x}sensitive regions (i.e., multiplied by the aCCF of NO_{x} emissions), while for contrails, it depends on the distance flown in persistent contrail formation areas.
Heun's method is adopted for integrating the aircraft dynamical model along discretized segments of the route through each phase, i.e., climb, descent, and cruise (GonzálezArribas et al., 2023). Since the calculations are similar for different members, parallelization would be beneficial in reducing computational time. Here, CUDA (Guide, 2013; Klöckner et al., 2012), a tool for generalpurpose computing on the graphics processing unit, is employed to parallelize the computations.
3.2 Performance evaluation of a flight plan
The expected values obtained from Eq. (10) are for a specific flight plan. By these settings, for this flight plan, the cost function given in Eq. (4) can be evaluated with the following equation:
Figure 7 shows how the expected performance is calculated and evaluated for a given flight plan and an ensemble weather forecast.
The objective is to find a flight plan that minimizes Eq. (13). Since the flight plan includes both discrete and continuous decision variables, the optimizer should be capable of solving the optimization within hybrid decision spaces. A classical approach to solving such optimization problems is mixedinteger nonlinear programming, which is mathematically complex and computationally intensive.
3.3 Probabilistic flight plan generation
To solve the optimization with hybrid decision variables in an efficient manner, instead of directly searching for the optimal flight plan, the optimization is conducted in the space of probability distributions defined over flight plans. In other words, instead of directly minimizing min_{ℱ}J(ℱ), the minimization of the following equivalent problem is considered:
which provides the same optimal solution; i.e., if ℱ^{*} is the optimal solution to the original problem, $\mathbf{P}(\mathcal{F}={\mathcal{F}}^{*})=\mathrm{1}$ provides the same results. The optimization is carried out in the space of probability distributions to move to continuous search space. In addition, it can facilitate searching by parameterizing the distribution p with a parameter vector θ∈ℝ^{Θ}, approximating Eq. (14) as
which is generally not identical to Eq. (14), as it relies on whether the parameterization is able to capture a distribution wherein $\mathbf{P}(\mathcal{F}={\mathcal{F}}^{*})=\mathrm{1}$. Thus, the parameterization of the distribution p using the vector θ plays an important role in approximating the original problem. In GonzálezArribas et al. (2023), the probabilisticexecution flight plan (𝒫ℱ) is introduced to parameterize the distribution over the space of possible flight plans. In this approach, the parameter vector θ is defined as follows:
where $\mathrm{{\rm Y}}\in {\mathbb{R}}^{{n}_{\mathrm{sp}}}$ is a vector assigning probability to select each airway, and $\widehat{\text{M}}$, $\widehat{\text{FL}}$, $\widehat{\mathcal{C}}$, $\widehat{\mathcal{D}}$, and d_{𝒟} are the parameterized Mach schedule, flight level, climb profile, descent profile, and distance to go, respectively. For a parameter vector θ of a 𝒫ℱ, flight plans can be randomly generated. For instance, let us explain how the lateral path is sampled from a given vector of parameters Υ. First, assume that all waypoints of the graph are processed and limited to two outgoing airways. This can be done by considering virtual edges of zero length for waypoints with more than two outgoing airways. The vector composed of a set of junctions is defined as $\overline{V}=\mathit{\{}{v}_{k}{\mathit{\}}}_{k=\mathrm{1}}^{{n}_{\mathrm{sp}}}\in {\mathbb{R}}^{{n}_{\mathrm{sp}}}$. At kth junction waypoint, the selection of an airway is done through a random process: for kth entry of a given vector Υ (i.e., υ_{k}) and kth entry of a vector containing uniform random variables ${\mathit{\xi}}_{k}(\in \mathrm{\Xi}\sim U(\mathrm{0},\mathrm{1}\left)\right)$, the airway is selected as
where S(⋅) is a sigmoid function: $S\left(x\right)=\mathrm{0.5}[\mathrm{1}+x(\sqrt{\mathrm{1}+{x}^{\mathrm{2}}}{)}^{\mathrm{1}}]$. Therefore, for a given vector of parameters Υ and a given vector of random variables, a lateral path (restricted to the structure of airspace) can be sampled. The approach to sampling a complete flight plan from p(ℱ;θ), such as sampling Mach schedule and flight level, is presented in GonzálezArribas et al. (2023) (see Algorithm 1). It is worth mentioning that the operational aspects and feasibility of the aircraft trajectory are considered within the generation of flight plans. For instance, continuous Mach adjustment is avoided, and the frequency of flight level changes is limited.
3.4 Optimization: augmented random search
In the probabilisticexecution flight plan approach, to sample the flight plan from a given θ, vectors containing random variables are required. Thus, the 𝒫ℱ associated with θ is stochastic. To evaluate ${\mathbb{E}}_{\mathit{\theta}}\left[J\right(\mathcal{F}\left)\right]:={\mathbb{E}}_{p(\mathcal{F};\mathit{\theta})}\left[J\right(\mathcal{F}\left)\right]$ for a given θ, sampling of multiple flight plans is required. By generating N_{FP} sets of random variables, N_{FP} flight plans (i.e., FP_{j} for $j=\mathrm{1},\mathrm{\cdots},{N}_{\mathrm{FP}}$) can be sampled independently for a given θ. To benefit from the provided N_{EPS} probable realizations of meteorological variables and N_{EPS} samples of exogenous sources of uncertainty (i.e., initial flight time and initial flight mass), we sample N_{EPS} (i.e., N_{FP}=N_{EPS}) potential flight plans for a given θ. In other words, each sampled flight plan is evaluated with one realization of meteorological conditions. In this respect, one can rewrite Eqs. (10) and (13) as
The flight planning problem can now be expressed as the following stochastic optimization problem:
where the objective is to find the optimal value of θ (which parameterizes the probability distribution p(ℱ;θ)) such that a population of randomly sampled flight plans minimizes the expected cost in Eq. (18). Sampling N_{EPS} flight plans and using them in a pairwise manner with ensemble members for trajectory integration to evaluate Eq. (17) is computationally more efficient than sampling a different number of flight plans (N_{FP}≠N_{EPS}) and then evaluating each sampled flight plan with all the ensemble members. This is because, in this case, we only integrate the aircraft trajectory N_{EPS} times instead of performing N_{FP}×N_{EPS} trajectory integrations. In spite of reducing the number of computations, it provides similar results. This is due to the fact that, despite sampling multiple flight plans for the evaluation of the objective function in Eq. (17) for a given θ, as the process goes by, all flight plans converge to a unique flight plan. For instance, let us consider the process of sampling the lateral path presented in Sect. 3.3. The choice of each airway relies on two parameters: υ and ξ. For the first iterations of the optimization algorithm, all outgoing airways from waypoints are almost equally probable. However, with more iterations, the decision variable θ is improved, leading to increasing or decreasing the parameter υ. This parameter converges to a large positive or negative value for the last iterations. For instance, in the case of large υ, ${lim}_{\mathit{\upsilon}\to \mathrm{\infty}}S\left({\mathit{\upsilon}}_{k}\right)=\mathrm{1}$, implying S(υ_{k})≥ξ_{k} for all the sampled flight plans. Thus, in the end, for the optimal value of θ obtained from optimization, we receive a unique flight plan. In this case, for the similar flight plans that are sampled for the last iterations, the expected flight performance variables obtained (from Eq. 17) using N_{FP}×N_{EPS} and N_{EPS} times trajectory integration will be almost the same.
In this work, the V1 version of the augmented random search (ARS) algorithm adopted from Mania et al. (2018) and GonzálezArribas et al. (2023) is employed. This is a gradientlike algorithm, which starts from an initial point θ_{0} and then generates n random search directions ω∈ℝ^{Θ}. In the next step, it evaluates ${\widehat{J}}^{+}:=\widehat{J}(\mathit{\theta}+\mathbf{S}\mathit{\omega})$ and ${\widehat{J}}^{}:=\widehat{J}(\mathit{\theta}\mathbf{S}\mathit{\omega})$, where $\mathbf{S}\in {\mathbb{R}}^{\mathrm{\Theta}\times \mathrm{\Theta}}$ is a diagonal matrix, adjusting the relative variations that are allowed between decision variables. Then, the decision variables θ are improved along all search directions proportional to ${\widehat{J}}^{+}{\widehat{J}}^{}$. As can be concluded, at each iteration, we generate 2n different vectors of decision variables (i.e., n search directions with θ+Sω and θ−Sω), and for each decision variable, we sample N_{EPS} flight plans. Therefore, we need to perform trajectory evaluation (i.e., TI) 2n×N_{EPS} times for each iteration. As these calculations are similar (or very similar) and independent from each other, the parallelization on GPUs is beneficial, enabling very fast function evaluation. The algorithms, details on the optimization approach, and also parallelization on GPUs are provided in GonzálezArribas et al. (2023).
The effectiveness of the proposed optimization algorithm to plan robust climateoptimal aircraft trajectories with respect to uncertain meteorological conditions is analyzed for a flight from Frankfurt to Kyiv for two different days and departure times.

Scenario 1: 13 June 2018 at 00:00 UTC is a scenario in which the aircraft flies through areas favorable for the formation of persistent contrails.

Scenario 2: 10 December 2018 at 12:00 UTC is a scenario with no formation of persistent contrails.
The dominant climate impact of contrails is the main reason for selecting these two scenarios, providing better insight into the mitigation potentials.
For the route graph, the full airspace graph of the considered days is filtered and processed to include all paths from the end of the standard instrument departures of the origin airport to the beginning of the standard instrument arrivals of the destination airport with the maximum length of 104 % of the shortest path length. The considered aircraft is an Airbus model A320214, with the engine CFM565B4/P. Table 2 provides the required parameters to calculate the NO_{x} emission index of the considered aircraft using BFFM2. The initial flight time and initial mass are modeled as Gaussian variables: t_{0}∼𝒩(00:00 UTC, 660) [s] for scenario 1, t_{0}∼𝒩(12:00 UTC, 660) [s] for scenario 2, and m_{0}∼𝒩(61 600, 164) [kg]. The considered standard deviations for initial flight time and initial flight mass are adopted from the studies conducted by Dalmau et al. (2021) and Benjumea (2003), respectively.
As for meteorological input data, due to ease of availability, the ERA5 Reanalysis data products containing 10 ensemble members are adopted for this study. It is worth mentioning that forecast data with more ensemble members can be similarly employed. Simulations are launched on the NVIDIA GeForce RTX 3090 graphics card, providing 10496 CUDA cores at a clock speed of 1.4 to 1.7 GHz.
(DuBois and Paynter, 2006)The weighting parameters of the objective function given in Eq. (4) are selected as ψ_{CST}=α [–] and ${\mathit{\psi}}_{\mathrm{CLM}}=(\mathrm{1}\mathit{\alpha})K$ [USD K^{−1}], where K is a scaling (or conversion) factor determined as
where, for instance, SOC_{climate} is the SOC calculated when the optimization objective is only the climate impact or ATR_{cost} is the ATR when the objective is only SOC. $\mathit{\alpha}\in [\mathrm{0},\mathrm{1}]$ is a weighting parameter that penalizes cost versus climate impact in which α=0 is the pure costoptimal and α=1 is the pure climateoptimal routing strategy. In the simulations, we consider five different values for α in order to explore the tradeoff between operating cost and climate impact represented by SOC and ATR, respectively.
4.1 Scenario 1: formation of persistent contrails
We consider a scenario in which the aircraft flies through warming contrails for the costoptimal routing option. Before presenting the results, the performance of the proposed optimizer in terms of convergence and computational time is analyzed. Since the optimization approach is stochastic, different results may be obtained with different executions. To explore the sensitivity of the optimization method, 50 different runs are performed with similar settings for pure cost (i.e., α=1.0) and pure climateoptimal (i.e., α=0.0) routing options. Then, the objective gap is calculated considering the best performance obtained from different solutions (i.e., the minimum value of J) as the reference. The convergence performances with averaged values as solid lines, and 0, 10, 90, and 100 percentiles are depicted in Fig. 8. For both cases, with around 700 iterations (≈2.8 s) for the costoptimal trajectory and 900 iterations (≈3.6 s) for the climateoptimal one, the estimated objective gaps are quickly reduced up to 1 % of the values of the objective functions J. As the climateoptimal routing option is associated with the inclusion of aCCFs calculated from meteorological variables, the optimization is much more complex, which can be validated in Fig. 8. With around 4000 iterations (16 s), the objective gap is reduced up to 0.5 % of J. Consequently, with a maximum of 4000 iterations, nearly optimal performance can be obtained with +0.5 % maximum deviation around the bestobtained value (i.e., the most optimal case).
Now, we proceed to present the obtained results. The aircraft profiles and climate responses for different routing options are given in Fig. 9. The SOC depends on the flight time and fuel consumption. Therefore, the aircraft for routing strategies with higher values of α, such as $\mathit{\alpha}=\mathrm{1.0},\mathrm{0.8}$, tends to fly at higher altitudes within the vertical constraints because flying at higher altitudes is beneficial for reducing fuel consumption, which contributes a large part of the total operating cost (see Fig. 9a). By analyzing the lateral paths depicted in Fig. 10 with the direction and speed of wind at different flight levels, one can see that the aircraft deviates from the shortest path to benefit from stronger tail winds. For trajectories with lower climate impacts, as can be seen in Fig. 9a, the aircraft flies at relatively lower altitudes compared to the costoptimal routing option, mainly to avoid the formation of persistent contrails (due to warming impacts during nighttime). The climateoptimal routing options reduce the warming effect of contrails. Although the warming climate effects of NO_{x} emissions and water vapor increase with climateoptimal trajectories, the net climate impact decreases. This is because the climate impact of contrails outweighs the impact associated with other species (as discussed in Sect. 2.2.2). The contribution of each species to total climate impact, variability of the obtained climate impacts, and SOC with ranges of uncertainty for different α values and Pareto frontiers are provided in Fig. 11. For a specific case (α=0.2), by accepting an increase of 4.3 % in cost, there is a potential to mitigate the climate impact by 53.0 % considering median performance. In Sect. 2.3, it was shown that the variability of relative humidity among ensemble members is high, leading to high uncertainty in the aCCF of contrails. As expected, the obtained climate impact of contrails is highly uncertain when the aircraft flies through areas sensitive to forming persistent contrails. In contrast, as the aircraft tends to avoid PCFAs, the ranges of uncertainty are reduced, in which, for the complete avoidance that is achieved with α=0.2, the climate impact is almost deterministic. In addition, SOC requires flight time and fuel consumption to represent operating cost in USD, and as it is affected by relatively less uncertain meteorological variables (i.e., wind and temperature compared to relative humidity for the considered case study as analyzed in Sect. 2.3), the uncertainty in its value is small.
By analyzing the contribution of each species to the net ATR for different α, one can conclude that the mitigation potential is achieved mainly by avoiding contrailsensitive areas, which results in slight increases in the climate impact of NO_{x} emissions (see Fig. 11a). However, when the formation of persistent contrails is completely avoided (i.e., with α=0.2), the optimizer tends to reduce NO_{x} emissions mainly by reducing speed to reduce the fuel flow required to calculate the NO_{x} emission index and also total NO_{x} emissions (i.e., NO_{x} emissions = NO_{x} emission index ⋅ fuel consumption). Reducing NO_{x} emissions by flying at lower speeds is achieved at the expense of a considerable increase in flight time and, consequently, SOC. As can be concluded from Pareto frontiers, such a reduction in climate impact for this scenario is not cheap as only 5 % more reduction in climate impact is obtained with almost 3 % more increase in SOC (α=0.0). As the aCCF of contrails is only evaluated in areas favorable for the formation of persistent contrails, typically determined using inequality constraints, it has sharp spatial behaviors (i.e., PCFA (latitude, longitude, altitude, time) ∈ {0, 1}). In addition, contrails have dominant climate effects. Therefore, the optimizer's first choice is to avoid forming persistent contrails, which may be achieved more efficiently than reducing the impacts of other species with relatively lower climate impacts and smooth spatial behaviors. This can be validated in Fig. 11a, as the lowest priority is given to reducing the climate impact associated with NO_{x} emissions (for α=0).
4.2 Scenario 2: no formation of persistent contrails
In the next scenario, we analyze the mitigation potential when no persistent contrails are formed with the costoptimal routing option.
For this case, aircraft profiles and climate responses are depicted in Fig. 12a and b, respectively. As can be seen in Fig. 12a, the optimizer chooses to fly at lower altitudes for routing strategies with higher penalization on climate impact. As no persistent contrails are formed (see Fig. 13a), we depict the lateral paths with the merged aCCFs (calculated using the mean values of the obtained NO_{x} emission index) as the color map at different flight levels in Fig. 13b. As can be seen, flying at lower altitudes is more beneficial in reducing the climate impact of other species (mainly NO_{x}). In addition to lowering cruise altitude, the aircraft flies at lower speeds to reduce the fuel flow and, consequently, fuel consumption, NO_{x} emission index, and NO_{x} emissions. The variability of climate impact and SOC for different α values and Pareto frontiers are given in Fig. 14. By reducing α, the climate impact decreases at the cost of an increase in SOC. For instance, for α=0.2, by accepting a 1.1 % increase in cost, a 16 % reduction in ATR can be achieved. As in the previous case, the relative increase in SOC is considerable for α=0, since the aircraft tends to fly at a relatively lower speed for more reduction in climate impact.
In conclusion, climate impact reduction is achieved at the expense of a higher cost increase than in the previous scenario. Moreover, since no contrails are formed, the uncertainty in climate impact is almost negligible.
This paper presented a methodology to plan a robust climateoptimal aircraft trajectory under uncertain meteorological conditions. The climatesensitive regions were identified using the prototype algorithmic climate change functions (version 1.1). The ensemble prediction system was employed to characterize uncertainty in weather forecasts. A heuristic algorithm was employed and implemented on graphics processing units to solve the proposed robust trajectory optimization in a computationally fast manner. The effectiveness of the proposed approach was explored in two scenarios. Discussion of the obtained results (mainly related to the achieved mitigation potentials, current limitations, and future lines of research) and some general remarks are presented in the following.
The obtained mitigation potentials for the considered scenarios were different due to the variability of meteorological conditions. In both cases, the climateoptimal routing options could reduce the climate impacts. The costoptimal trajectories flew at higher altitudes compared to climateoptimal ones, as flying at higher altitudes is beneficial for reducing fuel consumption. This is also in line with related studies in the literature (e.g., Yamashita et al., 2020). Due to the dominant climate impact and nonsmooth spatial behavior of contrails, the mitigation potential obtained for the scenario with contrail effects (i.e., scenario 1) was higher than the scenario with no formation of persistent contrails. The nonsmooth spatial behavior of the contrail climate impact is related to the conditions of PCFA. Due to the high variability among the ensemble members of relative humidity over ice needed to determine PCFA, the climate impact of contrails was highly uncertain. However, for the cases without contrail formation, the total climate impact was almost deterministic. This is because the variability in the other weather variables was almost negligible. For both scenarios, α=0.2 seems to be a reasonable choice, since the climate impacts were reduced at the expense of acceptable increases in the operating cost, and the results were almost deterministic.
In spite of considering the ensemble members in trajectory planning, a unique (or deterministic) flight plan is determined. This reflects the operational feasibility and applicability of this method since, in the flight planning context, the requirement is to determine a unique lateral route in latitude and longitude that starts and ends at predefined points in space and follows the real structure of airspace as well as having a fixed altitude profile and a fixed airspeed schedule. In this case, the effects of uncertainty are reflected in the aircraft performance variables. For instance, let us consider the first scenario. For α=1.0, the lateral path, speed schedule, and flight level were determined in a costoptimal manner (see Figs. 9a, 10). In our approach, the optimized trajectory is assumed to be tracked as close as possible in practice with the system's lowlevel controllers in real time. Aiming to optimize a unique flight plan, the proposed method provides some ranges for the aircraft performance variables, such as fuel consumption, flight time, and climate impact, due to different probable realizations of weather conditions. For instance, the climate impact is expected to lie within the determined ranges (see Fig. 9b) if the considered ensemble members could acceptably predict future weather conditions.
One of the next steps should be analyzing the feasibility of such a routing strategy for real traffic scenarios. In fact, air traffic management (ATM) is a complex multiagent system that cannot be represented by individual elements but by their collective behavior at the network scale. It was shown in the paper that for the climateoptimal routing options, the aircraft tends to fly at relatively lower altitudes compared to the costoptimal one. Such behavior to avoid climatesensitive areas may result in more congested areas, raising some challenges, particularly increased workload, complexity, and conflicts. Thus, the mitigation potentials reported at the microlevel may not be achievable considering real traffic scenarios. Therefore, after generating climatically optimal flight plans, one needs to assess the resulting effects at the network scale and perform a resolution (typically modeled as an optimization problem) to restabilize the ATM system by compensating for the negative impacts while keeping the modified trajectories as close as possible to inputted climateoptimized ones. The assessment of manageability of climateoptimal trajectory planning in an ATM system is called macroscale analysis and lies outside the scope of this paper (see Baneshi et al., 2023, for a study in this area).
In this study, we only considered the minimization of the expected performance, e.g., expected climate impact. However, the concept of robustness is mainly related to having less uncertain results (i.e., also minimizing the uncertainty range). In the case of robustness to meteorological uncertainty, we need to find a flight plan that avoids areas of airspace with high variability among the ensemble members. For instance, in GonzálezArribas et al. (2018), the dispersion of the arrival time is minimized by avoiding regions with high variability among the ensemble members of wind (characterized by SD). In this study, we observed that the most uncertain variable is the climate impact of contrails. Since, for both scenarios, only warming contrails were formed, minimization of the expected values directly led to avoiding uncertain PCFAs, and as can be seen from the results, for the climateoptimized trajectories, we obtained robust solutions (almost deterministic results). However, during the daytime, different behavior is expected for the cases with the cooling contrails. This is because, to reduce the expected climate impact, aircraft will tend to fly in uncertain PCFAs to benefit from cooling effects. In the next versions of ROOST, such scenarios will be taken into account, and controlling the dispersion of all flight performance variables will be addressed by including their variance in addition to the averaged values as objectives in the objective function.
Regarding the computational time and convergence performance, it was shown that they are scenariodependent. For more complex problems, such as the case including climate impacts quantified by using aCCFs, the optimizer required more iterations to enhance the convergence compared to the costoptimal routing option. It is worth mentioning that the distance between the origin and destination, available route graphs, and also parameters within the optimization algorithm can change the convergence performance and computational time. The number of iterations is a userdefined parameter that needs to be specified based on the required performance and availability of computational resources. In the performed simulations, we considered 4000 iterations. By looking at the Pareto frontiers, it is clear that the optimizer was able to find nearly optimal solutions. Thanks to the parallelization on GPUs, the computational time for achieving a nearly optimal solution is promising. There are several controlling parameters within the optimization algorithm of ROOST, including the number of search directions, the augmented random search (ARS) step size, and the Nesterov velocity factor (see GonzálezArribas et al., 2023, for a description of these controlling parameters) that can affect the convergence performance. In our future work, we aim to examine the effects of all these parameters and propose an optimal selection of them. Moreover, adaptive (scenariodependent) stopping criteria will be proposed, helping to optimize aircraft trajectories more efficiently in the sense of computational time.
To explore the tradeoff between climate impact and the operating cost, Pareto frontiers were generated. By changing the weighting parameter α, different Paretooptimal solutions were obtained. However, with this approach, a specific value for α does not necessarily result in a similar cost increase and climate impact mitigation potential for different scenarios (e.g., for α=0.2, the climate impact is reduced by 53 % and 16 % by accepting a 4.3 % and 1.1 % increase in the operating cost for scenarios 1 and 2, respectively). This approach is suitable for analyzing the mitigation potential for a single flight. However, it is not the most efficient way to study the Paretooptimal solutions of the aggregated results of a large number of flights. For such cases, having the flexibility to directly optimize flights, requesting a certain range reduction in climate impact or allowing a specific range for the increased operating costs would be beneficial. In this respect, one can aggregate climateoptimal trajectories having, for instance, a 0.5 % to 1.0 % increase in the operating cost. In future versions of ROOST, we will add this feature to the optimization tool by defining some path constraints. Such an aggregation of results is doable with α; however, one needs to generate more points in the Pareto frontier in order to classify the results based on a certain percentage increase in cost or a percentage decrease in climate impact. Scaling up this trajectorylevel analysis to the network scale will increase the computational time by a factor of the considered α values.
As was explored in the paper, the mitigation of climate impact within the flight planning context is achieved only by accepting some extra costs due to avoidance of highly climatesensitive regions, which is in line with related studies in the literature (e.g., Yamashita et al., 2020, 2021; Lührs et al., 2021; Niklaß et al., 2019). One potential approach to offsetting these additional costs is implementing fees and taxes that account for the climate impact caused by aviation in order to motivate airliners to adopt climateoptimal trajectories. Currently, there is no climate policy for the aviationinduced nonCO_{2} climate effects in the planned marketbased instruments such as the emission trading scheme (ETS). However, several studies have proposed frameworks to include costs of aviationinduced nonCO_{2} climate effects in the operating cost, such as the concept of equivalent CO_{2} emissions and the concept of climatecharged airspace (Niklaß et al., 2021). If taxes are applied to the nonCO_{2} climate effects, climateoptimal trajectories can also be economically beneficial.
The aCCFs used in this study represent a prototype formulation. The aCCF algorithms were developed for meteorological summer and winter conditions with a focus on the North Atlantic flight corridor. Thus, the usage of the aCCFs for different seasons and regions needs special caution. However, further development of the aCCFs and an expansion of their geographic scope and seasonal representation represent ongoing research.
The robust aircraft trajectory optimization technique presented in the paper is released as an opensource Python library called ROOST V1.0 (Robust Optimization of Structured Trajectories). It is developed at https://github.com/AircraftOperationsLab/roost (last access: 2 July 2023) and is available via the DOI https://doi.org/10.5281/zenodo.7121862 (Simorgh, 2022). It is distributed under the GNU Lesser General Public License (Version 3.0). All the results presented in the paper were obtained using ROOST. It should be noted that the optimizer ROOST uses BADA4.2 (license granted for the activities developed within the FlyATM4E project) to represent the aerodynamic and propulsive performance of the aircraft. Due to restrictions imposed by the BADA license, the current version (in the GitHub repository) is incomplete, as three Python scripts related to the aircraft performance model used have been excluded (i.e., bada4.py, apm.py, and badalib.cu). We are currently assessing the existing opensource aircraft performance models in order to make the complete library available to the public. In principle, it is possible to use other performance models as long as the functional specification (not necessarily the internal implementation) is the same as the BADA pointmass model, i.e., the drag polar is a function of the same variables, and so on. Otherwise, small modifications to the code would probably have to be applied. Potential alternatives include the OpenAP model (Sun et al., 2020) and the model proposed by Poll and Schumann (2021).
The ERA5 datasets used in this study can be freely accessed from the respective repositories after registration. ERA5 data were retrieved from the Copernicus Climate Data Store (https://doi.org/10.24381/cds.adbb2d47, Hersbach et al., 2023, 2020).
Conceptualization, AS and MS; development of the library (ROOST V1.0), DGA; incorporation of climate impacts to the library, AS; algorithmic climate change functions, SD, SM, and FY; writing – original draft, AS and MS; writing – review and editing, AS, MS, DGA, SM, VG, SD, SB, HY, FL, BL, MMM, FY, and FC. All authors have read and agreed to the published version of the paper.
At least one of the (co)authors is a member of the editorial board of Geoscientific Model Development. The peerreview process was guided by an independent editor, and the authors also have no other competing interests to declare.
Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
This research was carried out as a part of the EU project FlyATM4E.
FlyATM4E has received funding from the SESAR Joint Undertaking under the European Union's Horizon 2020 research and innovation program (grant no. 891317). The JU receives support from the European Union's Horizon 2020 research and innovation program and the SESAR JU members other than the union.
This paper was edited by Andrea Stenke and reviewed by two anonymous referees.
AMSCouncil: Enhancing weather information with probability forecasts, B. Am. Meteorol. Soc., 89, 1049–1053, 2008. a
Appleman, H.: The formation of exhaust condensation trails by jet aircraft, B. Am. Meteorol. Soc., 34, 14–20, 1953. a
Baneshi, F., Soler, M., and Simorgh, A.: Conflict assessment and resolution of climateoptimal aircraft trajectories at network scale, Transport. Res. DTr. E., 115, 103592, https://doi.org/10.1016/j.trd.2022.103592, 2023. a
Bauer, P., Thorpe, A., and Brunet, G.: The quiet revolution of numerical weather prediction, Nature, 525, 47–55, https://doi.org/10.1038/nature14956, 2015. a, b, c
Benjumea, A. C.: Aspectos antropométricos de la población laboral española aplicados al diseño industrial, Instituto Nacional de Seguridad e Higiene en el Trabajo, Ministerio de Trabajo e inmigración, ISBN 8474256550, 2003 (in Spanish). a
Bickel, M., Ponater, M., Bock, L., Burkhardt, U., and Reineke, S.: Estimating the effective radiative forcing of contrail cirrus, J. Climate, 33, 1991–2005, 2020. a
Bonami, P., Olivares, A., Soler, M., and Staffetti, E.: Multiphase mixedinteger optimal control approach to aircraft trajectory optimization, J. Guid. Control, 36, 1267–1277, 2013. a
Campbell, S., Neogi, N., and Bragg, M.: An optimal strategy for persistent contrail avoidance, in: AIAA Guidance, Navigation and Control Conference and Exhibit, 18–21 August 2008, Honolulu, Hawaii, 6515, 2008. a, b
Dahlmann, K., Grewe, V., Frömming, C., and Burkhardt, U.: Can we reliably assess climate mitigation options for air traffic scenarios despite large uncertainties in atmospheric processes?, Transport. Res. DTr. E., 46, 40–55, 2016. a, b
Dalmau, R., Ballerini, F., Naessens, H., Belkoura, S., and Wangnick, S.: An explainable machine learning approach to improve takeoff time predictions, J. Air Transp. Manag., 95, 102090, https://doi.org/10.1016/j.jairtraman.2021.102090, 2021. a
Dietmüller, S., Matthes, S., Dahlmann, K., Yamashita, H., Simorgh, A., Soler, M., Linke, F., Lührs, B., Meuser, M. M., Weder, C., Grewe, V., Yin, F., and Castino, F.: A python library for computing individual and merged nonCO_{2} algorithmic climate change functions: CLIMaCCF V1.0, Geosci. Model Dev. Discuss. [preprint], https://doi.org/10.5194/gmd2022203, in review, 2022. a, b, c, d, e, f
DuBois, D. and Paynter, G. C.: “Fuel Flow Method2” for Estimating Aircraft Emissions, SAE Transactions, 1–14, 2006. a, b
Franco Espín, A., Rivas Rivas, D., and Valenzuela Romero, A.: Optimal Aircraft Path Planning in a Structured Airspace Using Ensemble Weather Forecasts, 8th SESAR Innovation Days, 3–7 December 2018, Code 157125, 2018. a
Frömming, C., Grewe, V., Jöckel, P., Brinkop, S., Dietmüller, S., Garny, H., Ponater, M., Tsati, E., and Matthes, S.: Climate cost functions as a basis for climate optimized flight trajectories, Air Traffic Semin., 239, 1–9, 2013. a
Gallo, E., Navarro, F., Nuic, A., and Iagaru, M.: Advanced Aircraft Performance Modeling for ATM: Bada 4.0 Results, in: 2006 ieee/aiaa 25TH Digital Avionics Systems Conference, IEEE, 1–12, https://doi.org/10.1109/dasc.2006.313660, 2006. a
Gierens, K. M., Lim, L., and Eleftheratos, K.: A review of various strategies for contrail avoidance, Open Atmospheric Science Journal, 2, 1–7, 2008. a
GonzálezArribas, D., Soler, M., and SanjurjoRivo, M.: Robust aircraft trajectory planning under wind uncertainty using optimal control, J. Guid. Control, 41, 673–688, 2018. a, b
GonzálezArribas, D., Baneshi, F., Andrés, E., Soler, M., Jardines, A., and GarcíaHeras, J.: Fast 4D flight planning under uncertainty through parallel stochastic path simulation, Transport. Res. CEmer., 148, 104018, https://doi.org/10.1016/j.trc.2023.104018, 2023. a, b, c, d, e, f, g, h, i, j
Graver, B. and Rutherford, D.: Transatlantic airline fuel efficiency ranking, 2017, International Council on Clean Transportation, report, http://www.theicct.org (last access: 2 July 2023), 2018. a
Grewe, V. and Dahlmann, K.: How ambiguous are climate metrics? And are we prepared to assess and compare the climate impact of new air traffic technologies?, Atmos. Environ., 106, 373–374, 2015. a
Grewe, V., Champougny, T., Matthes, S., Frömming, C., Brinkop, S., Søvde, O. A., Irvine, E. A., and Halscheidt, L.: Reduction of the air traffic's contribution to climate change: A REACT4C case study, Atmos. Environ., 94, 616–625, 2014a. a
Grewe, V., Frömming, C., Matthes, S., Brinkop, S., Ponater, M., Dietmüller, S., Jöckel, P., Garny, H., Tsati, E., Dahlmann, K., Søvde, O. A., Fuglestvedt, J., Berntsen, T. K., Shine, K. P., Irvine, E. A., Champougny, T., and Hullah, P.: Aircraft routing with minimal climate impact: the REACT4C climate cost function modelling approach (V1.0), Geosci. Model Dev., 7, 175–201, https://doi.org/10.5194/gmd71752014, 2014.b. a, b
Guide, D.: CUDA C programming guide, NVIDIA, July, 29, 31, https://docs.nvidia.com/cuda/cudacprogrammingguide/ (last access: 2 July 2023), 2013. a
Hansen, J., Sato, M., Ruedy, R., Nazarenko, L., Lacis, A., Schmidt, G. A., Russell, G., Aleinov, I., Bauer, M., Bauer, S., Bell, N., Cairns, B., Canuto, V., Chandler, M., Cheng, Y., Del Genio, A., Faluvegi, G., Fleming, E., Friend, A., Hall, T., Jackman, C., Kelley, M., Kiang, N., Koch, D., Lean, J., Lerner, J., Lo, K., Menon, S., Miller, R., Minnis, P., Novakov, T., Oinas, V., Perlwitz, Ja., Perlwitz, Ju., Rind, D., Romanou, A., Shindell, D., Stone, P., Sun, S., Tausnev, N., Thresher, D., Wielicki, B., Wong, T., Yao, M., and Zhang, S.: Efficacy of climate forcings, J. Geophys. Res.Atmos., 110, https://doi.org/10.1029/2005JD005776, 2005. a
Hartjes, S., Hendriks, T., and Visser, D.: Contrail mitigation through 3D aircraft trajectory optimization, in: 16th AIAA Aviation Technology, Integration, and Operations Conference, 13–17 June 2016, Washington, D.C., USA, 3908, 2016. a, b
Hersbach, H., Bell, B., Berrisford, P., Hirahara, S., Horányi, A., MuñozSabater, J., Nicolas, J., Peubey, C., Radu, R., Schepers, D., Simmons, A., Soci, C., Abdalla, S., Abellan, X., Balsamo, G., Bechtold, P.,Biavati, G., Bidlot, J., Bonavita, M., De Chiara, G., Dahlgren, P., Dee, D., Diamantakis, M., Dragani, R., Flemming, J., Forbes, R., Fuentes, M.,Geer, A., Haimberger, L., Healy, S., Hogan, R. J., Hólm, E., Janisková, M., Keeley, S., Laloyaux, P., Lopez, P., Lupu, C., Radnoti, G., de Rosnay, P., Rozum, I., Vamborg, F., Villaume, S., and Thépaut, J.N.: The ERA5 global reanalysis, Q. J. Roy. Meteor. Soc., 146, 1999–2049, https://doi.org/10.1002/qj.3803, 2020. a
Hersbach, H., Bell, B., Berrisford, P., Biavati, G., Horányi, A., Muñoz Sabater, J., Nicolas, J., Peubey, C., Radu, R., Rozum, I., Schepers, D., Simmons, A., Soci, C., Dee, D., and Thépaut, J.N.: ERA5 hourly data on single levels from 1940 to present, Copernicus Climate Change Service (C3S) Climate Data Store (CDS) [data set], https://doi.org/10.24381/cds.adbb2d47 (last access: May 2022), 2023. a
Jelinek, F.: The Advanced Emission Model (AEM3)Validation Report, Ratio, 306, 1–13, 2004. a
Klöckner, A., Pinto, N., Lee, Y., Catanzaro, B., Ivanov, P., and Fasih, A.: PyCUDA and PyOpenCL: A scriptingbased approach to GPU runtime code generation, Parallel Comput., 38, 157–174, 2012. a
Lee, D. S., Fahey, D. W., Skowron, A., Allen, M. R., Burkhardt, U., Chen, Q., Doherty, S. J., Freeman, S., Forster, P. M., Fuglestvedt, J., Gettelman, A., De Leon, R. R., Lim, L. L., Lund, M. T., Millar, R. J., Owen, B., Penner, J. E., Pitari, G., Prather, M. J., Sausen, R., and Wilcox, L. J.: The contribution of global aviation to anthropogenic climate forcing for 2000 to 2018, Atmos. Environ., 244, 117834, https://doi.org/10.1016/j.atmosenv.2020.117834, 2021. a, b, c, d
Lim, Y., Gardi, A., and Sabatini, R.: Optimal aircraft trajectories to minimize the radiative impact of contrails and CO_{2}, Energy Proced., 110, 446–452, 2017. a
Lührs, B., Niklass, M., Froemming, C., Grewe, V., and Gollnick, V.: Costbenefit assessment of 2D and 3D climate and weather optimized trajectories, in: 16th AIAA Aviation Technology, Integration, and Operations Conference, 13–17 June 2016, Washington, D.C., USA, 3758, 2016. a, b
Lührs, B., Linke, F., Matthes, S., Grewe, V., and Yin, F.: Climate impact mitigation potential of European air traffic in a weather situation with strong contrail formation, Aerospace, 8, 50, https://doi.org/10.3390/aerospace8020050, 2021. a, b, c, d, e
Mania, H., Guy, A., and Recht, B.: Simple random search of static linear policies is competitive for reinforcement learning, in: Advances in Neural Information Processing Systems, edited by: Bengio, S., Wallach, H., Larochelle, H., Grauman, K., CesaBianchi, N., and Garnett, R., Neural Information Processing Systems Foundation, Inc. (NeurIPS), ISBN 9781510884472, 1800–1809, 2018. a
Matthes, S., Schumann, U., Grewe, V., Frömming, C., Dahlmann, K., Koch, A., and Mannstein, H.: Climate optimized air transport, in: Atmospheric physics, edited by: Schumann, U., Springer, 727–746, https://doi.org/10.1007/9783642301834_44, 2012. a
Matthes, S., Grewe, V., Dahlmann, K., Frömming, C., Irvine, E., Lim, L., Linke, F., Lührs, B., Owen, B., Shine, K., Stromatas, S., Yamashita, H., and Yin, F.: A concept for multicriteria environmental assessment of aircraft trajectories, Aerospace, 4, 42, https://doi.org/10.3390/aerospace4030042, 2017. a, b
Matthes, S., Lührs, B., Dahlmann, K., Grewe, V., Linke, F., Yin, F., Klingaman, E., and Shine, K. P.: Climateoptimized trajectories and robust mitigation potential: Flying ATM4E, Aerospace, 7, 156, https://doi.org/10.3390/aerospace7110156, 2020. a, b
Matthes, S., Dietmüller, S., Dahlmann, K., Frömming, C., Yamashita, H., Grewe, V., Yin, F., and Castino, F.: Algorithmic climate change functions (aCCFs) V1.0A: Consolidation of the approach and note for usage, Geosci. Model Dev. Discuss., submitted, 2023. a, b, c
Niklaß, M., Gollnick, V., Lührs, B., Dahlmann, K., Froemming, C., Grewe, V., and van Manen, J.: Costbenefit assessment of climaterestricted airspaces as an interim climate mitigation option, J. Air Transport., 25, 27–38, 2017. a, b
Niklaß, M., Lührs, B., Grewe, V., Dahlmann, K., Luchkova, T., Linke, F., and Gollnick, V.: Potential to reduce the climate impact of aviation by climate restricted airspaces, Transport Policy, 83, 102–110, 2019. a, b
Niklaß, M., Grewe, V., Gollnick, V., and Dahlmann, K.: Concept of climatecharged airspaces: a potential policy instrument for internalizing aviation's climate impact of nonCO_{2} effects, Climate Policy, 21, 1066–1085, 2021. a, b
Penner, J. E., Lister, D., Griggs, D. J., Dokken, D. J., and McFarland, M.: Aviation and the global atmosphere: a special report of the Intergovernmental Panel on Climate Change, ISBN10 0521663008, ISBN13 9780521663007, 1999. a
Poll, D. and Schumann, U.: An estimation method for the fuel burn and other performance characteristics of civil transport aircraft in the cruise. Part 1 fundamental quantities and governing relations for a general atmosphere, Aeronaut. J., 125, 257–295, 2021. a
Ponater, M., Grewe, V., Sausen, R., Schumann, U., Pechtl, S., Highwood, E., and Stuber, N.: Climate sensitivity of radiative impacts from transport systems, in: Proceedings of an International Conference on Transport, Atmosphere and Climate (TAC), Proceedings of the TACConference, 26 to 29 June 2006, Oxford, UK, 190–196, 2007. a, b
Rap, A., Forster, P. M., Haywood, J. M., Jones, A., and Boucher, O.: Estimating the climate impact of linear contrails using the UK Met Office climate model, Geophys. Res. Lett., 37, L20703, https://doi.org/10.1029/2010GL045161, 2010. a
Scherer, C.: Global Market Forecast Cities, Airports & Aircraft 20192038, AIRBUS SAS, 31707, ISBN 9782955438246, 2019. a
Schmidt, E.: Die entstehung von eisnebel aus den auspuffgasen von flugmotoren, Schriften der Deutschen Akademie der Luftfahrtforschung, Verlag R. Oldenbourg, München, Heft 44, 5, 1–15, 1941. a
Simorgh, A.: AbolfazlSimorgh/roost: Initial Release (v0.1.0), Zenodo [code], https://doi.org/10.5281/zenodo.7121862, 2022. a, b
Simorgh, A., Soler, M., GonzálezArribas, D., Matthes, S., Grewe, V., Dietmüller, S., Baumann, S., Yamashita, H., Yin, F., Castino, F., Linke, F., Lührs, B., and Meuser, M. M.: A Comprehensive Survey on Climate Optimal Aircraft Trajectory Planning, Aerospace, 9, 146, https://doi.org/10.3390/aerospace9030146, 2022. a, b, c, d, e, f, g, h, i
Soler, M., Zou, B., and Hansen, M.: Flight trajectory design in the presence of contrails: Application of a multiphase mixedinteger optimal control approach, Transport. Res. CEmer., 48, 172–194, 2014. a
Sridhar, B., Ng, H. K., and Chen, N. Y.: Aircraft trajectory optimization and contrails avoidance in the presence of winds, J. Guid. Control, 34, 1577–1584, 2011. a, b
Sun, J., Hoekstra, J. M., and Ellerbroek, J.: OpenAP: An opensource aircraft performance model for air transportation studies and simulations, Aerospace, 7, 104, https://doi.org/10.3390/aerospace7080104, 2020. a
van Manen, J. and Grewe, V.: Algorithmic climate change functions for the use in ecoefficient flight planning, Transport. Res. DTr. E., 67, 388–405, 2019. a
Vitali, A., Battipede, M., and Lerro, A.: MultiObjective and MultiPhase 4D Trajectory Optimization for Climate MitigationOriented Flight Planning, Aerospace, 8, 395, https://doi.org/10.3390/aerospace8120395, 2021. a
WMO: Guidelines on ensemble prediction systems and forecasting, edited by: Svoboda, M., Hayes, M., and Wood, D. A., World Meteorological Organization Weather Climate and Water, 1091, ISBN 9789263110916, 2012. a
Yamashita, H., Yin, F., Grewe, V., Jöckel, P., Matthes, S., Kern, B., Dahlmann, K., and Frömming, C.: Newly developed aircraft routing options for air traffic simulation in the chemistry–climate model EMAC 2.53: AirTraf 2.0, Geosci. Model Dev., 13, 4869–4890, https://doi.org/10.5194/gmd1348692020, 2020. a, b, c, d, e, f, g
Yamashita, H., Yin, F., Grewe, V., Jöckel, P., Matthes, S., Kern, B., Dahlmann, K., and Frömming, C.: Analysis of aircraft routing strategies for north atlantic flights by using AirTraf 2.0, Aerospace, 8, 33, https://doi.org/10.3390/aerospace8020033, 2021. a, b, c
Yin, F., Grewe, V., Frömming, C., and Yamashita, H.: Impact on flight trajectory characteristics when avoiding the formation of persistent contrails for transatlantic flights, Transport. Res. DTr. E., 65, 466–484, 2018a. a, b
Yin, F., Grewe, V., van Manen, J., Matthes, S., Yamashita, H., Linke, F., and Lührs, B.: Verification of the ozone algorithmic climate change functions for predicting the shortterm NO_{x} effects from aviation enroute, in: International Conference on Research in Air Transportation (ICRAT), 26–29 June 2018, Castelldefels, Barcelona, Spain, https://repository.tudelft.nl/islandora/object/uuid:3eab3d00272842208d3e0eef4cae8a11?collection=research (last access: 2 July 2023), 2018b. a
Yin, F., Grewe, V., Castino, F., Rao, P., Matthes, S., Dahlmann, K., Dietmüller, S., Frömming, C., Yamashita, H., Peter, P., Klingaman, E., Shine, K., Lührs, B., and Linke, F.: Predicting the climate impact of aviation for enroute emissions: The algorithmic climate change function submodel ACCF 1.0 of EMAC 2.53, Geosci. Model Dev. Discuss. [preprint], https://doi.org/10.5194/gmd2022220, in review, 2022. a, b, c