A Cooperative Co-Evolutionary Optimisation Model for Best-Fit Aircraft Sequence and Feasible Runway Configuration in a Multi-Runway Airport

A careful arrival and departure sequencing of aircraft can reduce the inter-arrival/departure time, thereby opening up opportunities for new landing and/or take-off slots, which may increase the runway throughput. This sequence when serviced with a suitable runway configuration may result in an optimal aircraft sequence with a runway configuration that can process the maximum number of aircraft within a given time interval. In this paper, we propose a Cooperative Co-evolutionary Genetic Algorithm (CCoGA) to find the combined solution of a best-fit sequence with a feasible runway configuration for a given traffic demand at an airport. The aircraft sequence and the runway configuration are modelled as individual species, which can cooperatively interact with each other. Therefore, we computationally evolve the best possible combination of aircraft sequence (arrival and departure) and the feasible runway configuration. The proposed CCoGA algorithm is evaluated for Chicago O’Hare International Airport runway layout and resulting configurations. Arrival and departure traffic demand is modelled through a Poisson distribution. Two different arrival/departure sequencing methods, i.e., constraint position shifting with one, two and N-position shifting and first come first serve, are modelled. Runway configuration and traffic sequence (arrivals and departure) are modelled as two species, which are evolved co-operatively, through the CCoGA algorithm, to achieve the optimal traffic sequencing with a feasible runway configuration. Time-space diagrams are presented for the best-evolved population of arrival-departure sequence and runway configuration to illustrate the possibility of using available departure slots between arrivals to maximize capacity. Arrival-departure capacity envelopes are then presented to illustrate the trade-off between the arrivals and departures, given a runway configuration for each sequencing method. Results demonstrate the high mutual dependence between arrival-departure sequence and the runway configuration, as well as its effect on overall runway capacity. The results also demonstrate the viability of using evolutionary computation-based methods for modelling and evaluating complex problems in the air transport domain.


Introduction
The steady and continued growth of air traffic brings significant challenges to airport capacity and safety. Air traffic demand is expected to more than double in Europe and the U.S. [1,2]. According to the International Air Transport Association (IATA) 20-year Air Passenger Forecast (IATA 2016), by 2036, 7.2 billion passengers will travel by air annually. Over the years, airports have emerged as a bottleneck in the global air transportation network [3]. They are the source and sink and perhaps the most constrained resource in the air transport network. The reasons are partly due to the difficulty of expansion of the airport infrastructure due to cost, land use and environmental concerns. This fixed resource requires a careful consideration of both the arrival and departure sequence and the most suitable runway configuration (combination of runways in use) at a given airport to optimise given objectives, subject to a variety of operational constraints [4].
The FAA divides aircraft into three classes (i.e., heavy, medium (large) and small), based on the maximum take-off weight capacity and specifies the separation requirements during Instruments Flight Rules (IFR) approaches. For an airport with several runways, a large number of combinations of simultaneously active runways, assignments of aircraft types, and movements are possible; further potentially complicated by the existence of noise restrictions or certain weather conditions. Each of these combinations is called a runway configuration, and this may change several times a day given the complex interplay of factors. For example, Boston/Logan Airport has five runways, which can operate in about 20 different configurations [5].
Ad-hoc handling of arrival and departure in such a highly constrained environment results in inefficient usage of runway and loss of capacity at an airport; which may further results in delays; consequently, causing a huge economic loss [6]. The Air Traffic Controller (ATC) adopts various criteria for sequencing aircraft in the terminal airspace for increasing the runway capacity at busy times. The most common approaches are First Come First Serve (FCFS), Constraint Position Shift (CPS), a stream of arrivals and a stream of departures. A shift (forward or backward) of an aircraft from the FCFS has additional advantages such as similar wake category aircraft are positioned side-by-side in the sequence that reduces the total headway time, consequently increasing the total movement on the runway. Constraint shifting maintains some sense of equity among aircraft by not deviating too much from the original sequence.
Although runway capacity may be affected by a very large number of diverse factors, four main factors play a key role such as separation standard (determined by aircraft wake characteristics, communication/navigation/surveillance infrastructure and weather conditions), runway configuration, aircraft mix (heavy, medium and small) and sequencing strategies applied by ATCs [7]. For example, the impact of weather (such as low visibility conditions) may result in increased separation standards, few runway configurations operable (only those runways with suitable landing/take-off aids available) and limited ATC sequencing strategies (separate runways for arrivals and departures) being applicable.
Therefore, the question of the optimal arrival/departure sequence for maximizing runway capacity cannot be addressed in isolation without the careful consideration of runway configuration. These are important because the runway throughput is estimated based on the total span of the landing or the taking-off sequence, which depends on the consecutive wake vortex separation and the runway configuration used. Furthermore, the geometric layout and orientation of active runways may vary the runway capacity. For instance, the mutual dependency of two active, but converging runways may decrease the capacity (departures are held until the arriving aircraft vacates the runway to avoid a missed approach procedure leading to a loss of separation with departing aircraft as in Las Vegas Airport). In addition, runway assignment to an aircraft may vary the throughput because different runway assignment changes the position of the aircraft in the stream. Thus, choosing the most suitable runway based on the traffic pattern is quite challenging even for a highly experienced ATC. An optimal aircraft sequence coupled with the right runway configuration may result in increased runway throughput.
Identifying the optimal combination of an aircraft sequence and a runway configuration that can maximize the runway throughput is a nonlinear optimisation problem, which is NP-hard [8,9]. The nonlinearity and interdependency among components in this problem make the classical mathematical assumptions of linearity, homogeneity, normally independent and identically distributed observations obsolete. To address such a problem, a heuristic approach that can balance exploring the search space and exploit the knowledge gained during the search process may be well suited. However, an exhaustive search or single solution-based heuristic cannot solve this nonlinear NP-hard problem in less than exponential time with Pareto optimality [10,11]. To deal with this complexity, soft computing methods are generally used. Meta-heuristics might be an efficient way to produce acceptable solutions by trial and error to a complex nonlinear problem. This search technique has the potential to intensify and diversify the subcomponents those are related to each other [12,13]. In addition, meta-heuristic approaches have the capability to iteratively produce high-quality solutions by allowing the local search to escape from local optima.
For the aircraft sequencing problem, some meta-heuristic techniques have been proposed in the literature [14][15][16][17]. These approaches consider maximizing runway capacity by finding optimal sequences for arrival and/or departure. However, an integrated approach that considers both arrival and departure sequencing coupled with runway configuration is absent from the literature. In this paper, we propose an evolutionary computation-based approach [18] to search for an optimal arrival and departure sequence, as well as to identify the best-fit runway configuration that can maximize the runway throughput capacity.
A novel meta-heuristic approach is proposed for solving two sub-problems that can obtain a coordinated solution, while each sub-problem co-evolves cooperatively. A concept diagram of the proposed approach is presented in Figure 1. The problem is decomposed into subproblems, i.e., optimal aircraft sequence, best-fit runway configuration and the suitable runway assignment. The optimal sequence is selected based on the inter-arrival/departure time. Runway configuration is selected based on the throughput of the configuration based on the aircraft sequence. The runways are allocated to the aircraft based on the runway length and the aircraft wake category. A simple coordination mechanism (round-robin optimisation of subcomponents) is used to evolve the solution. An aircraft sequence of mixed traffic (i.e., heavy, medium and small) where the arrivals/departures follow the Poisson arrival process and the generated runway configurations with an assigned runway of each flight are fed into the model. Two independent optimisation algorithms co-operate with each other to find a feasible configuration and a best-fit sequence in each iteration. The server estimates the capacity of the configuration based on the obtained sequence. The safe separation standard is maintained during the capacity estimation. The estimated capacity is fed back to the model as a reward to the feasible configuration and the best-fit sequence for the next generation. This is modelled as a cooperative co-evolutionary algorithm based on a genetic algorithm.  The remainder of the paper is organized as follows: Section 2 presents the background study. The problem formulation is described in Section 3. Section 4 outlines the methodology and proposed algorithm. Section 5 presents the experimental design. Results and the analysis are discussed in Section 6. In Section 7, we present our conclusion and future work.

Background
In the aircraft sequencing, scheduling or runway configuration literature, it is often found that there is a research gap between the result from proposed algorithms and the method used in practice. The reasons are overly simple modelling of the problem and plain assumptions. As a result, an ideal solution might be achieved from the model, however fail to obtain expected performance in reality. To date, very few methods have been proposed in the reviewed literature for the multi-runway airport capacity. Approaches are lacking for solving large problem instances incorporating aircraft mixed sequence, runway assignment and the best-fit runway configuration.
In order to maintain smooth airport operation, runway utilization is essential. With the continued growth in traffic and constraints in building new runways (Kansai International Airport, Japan) or expanding existing runaways due to environmental constraints (Melbourne International Airport, Australia), runway capacity has emerged as the principal bottleneck of an Air Traffic Management System (ATMS) [19,20]. Maximizing the capacity of existing runways is one of the preferred approaches to handle this scenario.
Researchers have been inspired to apply several other approaches to enhance the runway capacity. Dynamic programming algorithms for runway scheduling under Constrained Position Shifting (CPS) and other system constraints have been proposed [21]. The CPS strategy shows that an aircraft's position in the optimised sequence does not deviate significantly from its position in the FCFS sequence. CPS ensures a degree of fairness in the aircraft operation. An Ant Colony System (ACS) algorithm incorporating Receding Horizon Control (RHC) was proposed in [22]. RHC reduces the computational time and helps to make the ACS algorithm tolerant to an uncertain dynamic environment. However, this approach considers a single objective that might be insufficient to capture the dynamic air and ground resources. A mixed integer programming approach was proposed for selecting a best-fit runway configuration [1,23,24]. This approach can determine the balance of arrival/departure to be served at any moment. Dear [25] proposed a mixed integer linear programming approach to investigate different performance indicators of the Aircraft Scheduling Problem (ASP). The purpose was to examine the formulation of the ASP. A dynamic programming approach [26,27] attempted to solve the aircraft landing problem. The objective was to minimize the cost in terms of landing delay. To incorporate the uncertainty in the airport capacity, an aggregate mathematical model for air traffic flow management has been proposed [28]. This approach considers the trade-off between airport arrivals and departures, as well as the coordination among the hub airports. A generalized dynamic programming is proposed for optimal scheduling of two runways [29]. In this approach, the same type of aircraft is grouped together as a batch, and precedence constraints are applied to sequence the aircraft.
Many exact and heuristic algorithms have been proposed to maximize runway capacity by optimising the aircraft scheduling, sequencing and landing [26,30,31]. A single runway-based integrated arrival/departure management approach has been proposed by [32]. This approach utilizes the compact time-indexed formulation defined through clique inequality that is solved by mixed-integer programming. The pruning rules have been introduced to reduce the complexity of the runway sequencing problem [33]. The pruning principles are integrated into a dynamic program that generates optimal sequences with low computational time.
However, the exact or heuristic algorithms perform well only for the problems that have small instances or a simpler problem. To address the complex problem, i.e., a problem dependent on several sub-problems, numerous metaheuristics have been proposed to search for good quality solutions. Meta-heuristics have two classes, single solution-based methods [34][35][36][37] and population-based methods. Single solution-based methods include, but are not restricted to tabu search [38][39][40] and simulated annealing [41][42][43]. A single solution-based method iteratively applies the generation of the neighbourhood solutions from the current single solution. This process iterates until a given stopping criterion, e.g., number of iterations. Therefore, single solution meta-heuristics might not achieve a good solution for complex nonlinear problems.
Population-based methods include, but are not restricted to ant colony optimisation [22,44], genetic algorithms [16,[45][46][47] and particle swarm optimisation [48]. Unlike single solution meta-heuristics, population-based meta-heuristics generate the initial population, then the replacement phase is started by selecting a new population from the previous population. This operation iterates until a given stopping criterion.
A robust optimization based on an artificial bee colony has been proposed to identify an optimal schedule by evaluating the robustness of feasible solutions [49]. The approach reduces the computational effort in the iterative relaxation procedure for min-max regret optimisation. To obtain a robust aircraft sequencing and scheduling model, a metaheuristic approach based on artificial bee colony has been proposed [36,42]. The approach considers improving the resilience of runway operation and the robustness of the makespan optimisation of aircraft sequencing and scheduling. Few hybrid approaches attempt to utilize the problem's interesting properties (e.g., receding horizon control) to solve the problem. Few approaches use the rolling horizon technique to solve large instances in a short computation time compatible with real-time applications [22,50]. All of these approaches claimed a significant improvement compared to the commonly-used air traffic control technique, e.g., FCFS.
Population-based meta-heuristics, e.g., the Evolutionary Algorithm (EA), have been used successfully on many numerical and combinatorial optimisation problems. Plain EA (e.g., GA) cannot perform as described when the problem space is large (i.e., high-dimensional) and no key objective measure exists to determine the fitness of the individuals [54]. However, the Co-Evolutionary Algorithm (CEA)-based meta-heuristic that follows the decomposition (i.e., divide and conquer) strategy performs well in high-dimensional problems [55,56]. In the decomposition method, a high-dimensional objective vector splits into smaller subcomponents that can be handled by conventional EAs. Since interdependencies may exist between subcomponents, co-adaptation is imperative in modelling such interdependencies during optimisation. Co-evolution can be regarded as a dynamic approach to implement the divide-and-conquer strategy [57]. Among the available EA's, CEA, the Cultural Algorithm (CA) [58] and Particle Swarm Optimisation (PSO) [59][60][61] are all promising methods.
CEAs have unique characteristics to handle the large nonlinear complex NP-hard problem where a trade-off solution is essential. For example, co-adaptation of subproblems, a divide and conquer strategy is very effective to obtain a near optimal solution that an exact method cannot achieve. The advantages of CEA-based meta-heuristics include solving large-scale nonlinear NP-hard problems and achieving the trade-off solution within a short computation time. This effective property is attributed to the decomposition of a large-scale problem into a set of smaller subproblems. This property of population-based (co-evolution) meta-heuristics motivates us to apply them to solving the complex, nonlinear NP-hard problem of mixed aircraft sequencing scheduling and best-fit runway configuration. Airport runway operations are performed in a highly interactive environment where the traffic scenario may change very quickly. The model proposes CCoGA, a meta-heuristic based on the cooperation of two co-evolution algorithms. The model optimises the interdependency of arrival and departure mixed traffic and the assignment of the suitable runway, as well as the runway configuration.

Research Gap and Contributions of the Research
Most of the literature reviewed in Section 2 is based on the simple model that considers an exact solution to solve the Aircraft Sequencing and Scheduling (ASS) or the runway assignment problem. The proposed model incorporates multiple objectives such as on-air tactical strategies (i.e., mixed aircraft sequencing) and the best-fit ground resources (i.e., runway configuration) to obtain a trade-off solution from many alternatives. The proposed model is able to take into account the terminal airspace constraints (i.e., safe separation, safe airspeed and suitable runway assignment) and the performance indicators (runway capacity).
The proposed approach uses the detail modelling of inter-arrival-departure time and runway occupancy time to estimate the headway between two consecutive runway operations. It also uses the best-fit runway configuration according to the mixed sequence in terms of capacity. However, our model does not deal with ground movements and ground delay. The contribution of the work is not the classical aircraft sequencing and scheduling problem. Rather, it deals with the coordinated management of aircraft mix sequence (arrival and departure), suitable runway assignment and the best-fit runway configuration that can achieve the trade-off solutions together. The trade-off solution means the concerted sequence and runway configuration active at that time horizon can provide the maximum throughput capacity.

Problem Formulation
We first formulate the runway capacity estimation model by adopting it from the literature [62,63], and then, we develop a runway configuration model.

Runway Capacity Estimation
Runway capacity is measured in several alternatives ways. All of them are intended to estimate the number of aircraft movements that can be performed on the runway system during the specific unit of time. A simple technique that can be found in the literature is to estimate the runway capacity by estimating the average inter-arrival-departure time, i.e., headway between consecutive aircraft [62]. We assume that aircraft arrival and departure follow the Poisson process, meaning that aircraft arrival and departure occur in a random fashion in the initial sequence. We modelled the Poisson arrival process for generating the initial aircraft sequence. The aircraft enters the FAF or RT with a time-varying rate of λ within the specific time horizon t as justified by [64,65]. Therefore, two possible scenarios are perceived, i.e., the opening and closing scenario. In an opening scenario, the approach speed of the lead aircraft is higher than the trailing aircraft; vice versa for the closing case. Table 1 lists all the variables, notations and the units used in this research. The headway (t a ij ) of two consecutive approaching aircraft i and j can be estimated mathematically as follows given that there is no position error.
where V a i , V a j are the lead and trail aircraft speed, respectively, γ is the common approach length, ROT is the runway occupancy time and δ ij is the separation between consecutive aircraft. In the literature, different separation matrices can be found in terms of distance [66,67] and time [21]. Most of the matrices are considered according to different aircraft classes. It is common that distance is used as a separation matrix between arriving aircraft and time is used as a separation matrix for departing aircraft. The proposed model uses both the distance and time as separation matrices. Table 2 lists the separation distances and runway occupancy time between two consecutive arrivals and departures of aircraft, respectively.

Notations Explanation and Units
A Set of all aircraft, i, j ∈ A R Set of all runways, r 1 , r 2 , r n ∈ R R o Set of all runway orientations, l a1a2 , l d1d2 , l a1d1 , l a1d2 , l a2d1 , l a2d2 ∈ R o N c Set of all available runway configurations in an airport at time t t a ij , t d ij The headway of two arriving and departing aircraft i and j (s) Approach speed of the leading and trailing aircraft, respectively (nmi/h) Similarly, headway (t d ij ) for two consecutive departing aircraft i and j can be estimated using the following equation.
where ρ is the length of the runway and ROT i is the runway occupancy time of aircraft i. The relative distance among aircraft changes frequently, especially in the terminal airspace. This position error may arise due to changing wind speed, The limitation of the radar system and navigation aids at the airport. In practice, a safety buffer is added to the minimum separations between aircraft to make up for imperfections in the ATC system. Usually, ATC allows a buffer of an additional 10 s between each aircraft for safety (10 s implies about 1/3 nautical mile. longitudinal separation). Air traffic controllers experience less intense traffic situations to cope with when separations are maintained at greater than the safety distance (i.e., by including a buffer safety distance) [68]. From a modelling perspective, the size of the buffer depends on the probability of the violation of the minimum separation and the aircraft speed [62]. The buffer time of two approaching or departing aircraft can be estimated as follows, Here, σ 0 is the standard deviation of the in-trail delivery error, q v is the value of the cumulative standard normal at the probability of violation P v and P v is the probability of violation of the minimum separation criteria between two aircraft. Therefore, the total headway for arrival becomes, Therefore, at any given time t, the inter-arrival-departure time can be expressed as a function of ξ, where A is the set of aircraft that is indicated as a sequence.

Runway Configuration
The definition of the runway configuration greatly impacts the capacity of the runway. Many airports in the world have multiple runways. Therefore, a large number of configurations with a different geometric layout, orientation and movements is possible.
For a multi-runway airport, the number of possible configurations is large. The availability of a certain configuration at most major airports might be restricted by factors such as wind, noise and installed instruments (e.g., Instrument Landing System (ILS)) prevailing at any particular time. For instance, a runway cannot be operated in the presence of strong crosswinds. Similarly, a runway might not be usable due to poor visibility conditions if it is not sufficiently instrumented for the operations.
In practice, the runway configuration is chosen based on experience. A configuration represents the combination of the runways that are active at a given time. Each configuration is associated with the arrival and the departure runways and their orientation. For instance, A 1 and A 2 are the active arrival runways and D 1 and D 2 are the active departure runways in a configuration, so the configuration becomes A 1 , A 2 |D 1 , D 2 . The active runways are represented on the basis of their bearing in degrees from Magnetic North divided by 10. In addition, parallel runway pairs are indicated by an L (Left), R (Right) and C (Centre) along with the numbers. For example, an operating runway configuration 31L|31R indicates that aircraft arrivals are handled on left runway 31L, which aligns 310 • from Magnetic North, and departures are handled on a parallel runway (right) 31R, which aligns 310 • from Magnetic North.
An ATC takes into account the capacity of each available configuration in making decisions concerning the best runway configuration to use at any given time. Manual approaches are simple procedures, often guided by common sense, that are meant to provide good, but cannot necessarily guarantee maximum capacity. Arriving-departing aircraft sequencing and best-fit runways are critical aspects of ensuring smooth inflow and outflow of traffic into and out of airports [69]. Therefore, this research has been motivated to address the nonlinear problem considering arrival, departure flight sequencing, runway allocation for the flights and the best-fit runway configuration for capacity maximization.
In the proposed model, a set of assumptions has been made for selecting a configuration from the available runway orientations (i.e., independent parallel, dependent parallel, cross, converging and diverging runways) [70]. The following assumptions are embraced in developing the model.

•
The crosswind/tailwind component for each runway in a configuration does not exceed the maximum threshold.

•
Arriving aircraft use the ILS. The arrival and departure direction of a configuration with a runway must be the same, i.e., θ a = θ d .

•
When multiple arrival/departure runways are active in a parallel orientation, one runway is dedicated for arrival, while another is dedicated for departure. However, in mixed mode operation, the runway direction should be the same.

•
The converging or diverging angle between two active runways should be greater than or equal to fifteen degrees, i.e., φ c ≥ 15 Runway orientation (l) is defined for the arrival-arrival, departure-departure and arrival-departure runway pairs of an active configuration. In a particular runway configuration, the maximum number of active runways can be four; with a maximum of two runways for arrival and two runways for departure used for an airport [71].
Before finding a feasible configuration, an aircraft sequence A, the set of runways R and the set of the runway orientation (R o ) are defined. Let an airport has N c runway configurations. At any time t, the configuration c k be chosen such that, where r a 1 k , r a 2 k and r d 1 k ,r d 2 k are the arrival and departure runway, respectively. The orientations between active runways for the k-th configuration are l k a 1 a 2 , l k where a 1 , a 2 , d 1 , d 2 denote the arrival and departure, respectively. Each configuration is evaluated by estimating the throughput for a given aircraft sequence as follows.
where ψ(c k ) is the throughput for the aircraft sequence a ∈ A and c k is the configuration.

Objective Function Formulation
The purpose of this study is to develop a model for obtaining an optimal aircraft sequence incorporated with a suitable runway configuration, which can maximise the runway throughput capacity. To meet this goal, an objective function (12) is formulated for obtaining an optimal sequence by minimizing the headway between consecutive aircraft in the sequence while optimally selecting the best-fit runway configuration.
A runway configuration identifies the set of runways that can be used under certain operating conditions for a given time constraint and has a major influence on a runway's operating efficiency and the overall capacity of an airport's runway system. For example, runway length is an influencing factor for allocating the runway to the flight. Heavy aircraft requires a longer run to take off, while medium or small aircraft needs the minimum required length. Therefore, the longer runway is assigned to the heavy aircraft when multiple runways are available in an active configuration. An objective function is designed (13) for obtaining the suitable runway configuration taking into account that it can achieve maximum capacity with the optimal aircraft sequence. The objective functions are formulated with suitable constraints.
At any time t, if the optimal configuration is c k for the aircraft sequence A, the objective function can be expressed as, where ζ(c k ) is the sum of inter-arrival-departure times. For a given time interval t, if the optimal aircraft sequence is A, a constrained optimisation function for the runway configuration can be expressed as: subject to : t a ij ≥ δ ij and t d ij ≥ ROT i where ψ(t) is the throughput at t with configuration c k while the aircraft sequence is A. ξ ij is the required headway for aircraft i and j. δ ij is the separation between aircraft i and j. H R l i and M,S R l j are the length of runway assigned for heavy and medium/small aircraft, respectively. ROT i is the runway occupancy time of leading aircraft i. The number of aircraft in the sequence is N, while τ is the unit of time for simulation. θ a and θ d represent the arrival and departure direction, respectively. φ c and φ d indicate the converging and diverging angle of the active runway when two runways are operated simultaneously. c k (n) presents the number of runways in configuration k.

Methodology
In this research, a population-based meta-heuristic approach based on a genetic algorithm is used to address the complex nonlinear NP-hard problem. The methodology co-adapts the aircraft sequence and the runway configuration simultaneously. To the authors' best knowledge, this is the first attempt to co-adapt the aircraft sequence and the available runway configuration within a time horizon for maximum runway capacity. The flow diagram in Figure 2 illustrates the proposed model. The optimisation process starts with the initialization of the sub-components as follows: • a population of aircraft sequences • a population of runway configurations The process for arrivals or departures is characterized as a Poisson distribution with a specified mean of arrival or departure rate [62]. The initial population of aircraft sequence is generated based on the Poisson arrival process at the Final Approach Fix (FAF) for arriving aircraft and at the runway threshold for departing aircraft. The initial population of runway configurations is initialized assuming the runway set of O'Hare International Airport, while feasible configurations are chosen based on the rule set discussed in Section 3.2. A set of individuals is selected from the population based on the fitness of each individual through a tournament selection. The selected individual (aircraft sequence) is decoded as a genetic chromosome and sent to CCoGA. Similarly, a set of individuals is selected from the population of runway configurations and decoded into a chromosome and sent to CCoGA in which genetic operation, i.e., crossover and mutation, is applied to evolve the population. In a crossover operation, the genetic material (i.e., gene) of two chromosomes is swapped to produce a new individual.
The crossover process establishes exchanging ATC actions between two chromosomes (parents) and generating new offspring (children). The upper and lower limit of possible values for each chromosome are given, which reduce the size of the search space. Based on initial experimentations, Simulated Binary Crossover (SBX) is used with high cross-over probability (e.g., 90%) [72]. The SBX crossover helps in generating offspring near the parents. Therefore, the crossover guarantees that the extent of the children is proportional to the extent of the parents and also favours that near parent individuals are monotonically more likely to be chosen as children than individuals distant from the parents.
The mutation process includes random modifications of the execution time associated with each ATC action in a chromosome. Mutation maintains the genetic diversity of one generation of a population to the next generation. The mutation here allows jumps in the search space of a decision variable, thus giving better chances of escaping from local optima. A tournament selection process involves running several tournaments among a few individuals chosen randomly from the population. The winner of each tournament (the one with the best fitness) is selected for crossover. If the tournament size is larger, weak individuals have a smaller chance to be selected. In this study, the tournament size is selected based on 10% of the total population size [73].  Figure 2. The process flow of the Cooperative Co-evolutionary Genetic Algorithm (CCoGA): two separate populations co-evolve mutually to achieve the best-fit sequence and runway configuration.

Selection of species
Aircraft sequence and runway configuration evolve individually through crossover and mutation as explained in Sections 4.1 and 4.2. The fitness of the evolved aircraft sequences and runway configurations is evaluated in terms of throughput capacity and yields a new population. After each evolution, the fitness of each individual is estimated and the fitness recorded for the next generation.
Again, a new population is generated from the elite individuals, and the fitness of each individual is determined. The new population is sent to the next generation of evolution. The process is repeated until the solution has converged.
The solution of the model is obtained by combining the cooperative solution of two subcomponents. This means that the fitness of a candidate solution is evaluated on the basis of the combined fitness of each subpopulation. A candidate solution from each population is selected and fed into the CCoGA algorithm, which simulates and evaluates the traffic scenario and configurations. The objective value is then used to update the new generation. At the end of each generation of constrained maximization and minimization, the fitness of the variables is updated, and the process is repeated until convergence.

Genetic Representation of Aircraft Sequence
A suitable chromosome representation must be defined for an evolutionary algorithm to work. A chromosome is a set of genes that defines a proposed solution to the problem. The set of all solutions is known as the population. Each subcomponent of the problem is represented as a genetic chromosome. The chromosome has a set of parameters (genes) characterizing the temporal aircraft type (heavy, medium and small), the speed of each wake category and the assigned runway of an active configuration. Figure 3 illustrates the encoding for aircraft sequences. The values of the genes within the chromosome are as follows: • (F n ) Three different gene values that represent the aircraft wake category (heavy, medium and small). Each wake category is randomly selected based on the traffic distribution as shown in Table 2. • (V) Each wake category of aircraft maintains a standard speed during landing or taking off. Speed determines the separation maintained with the aircraft, which is following or leading. • (R a ) Each aircraft is assigned a runway before landing or taking off. Assigning different active runways changes the position of the aircraft in the sequence.  Evolutionary techniques change the genes and their position in the chromosome through various operators, e.g., in this study, selection and mutation. This genetic operation propagates the characteristics of parents to offspring. In this research, genetic algorithm properties are adapted to obtain the most feasible solution among the available alternatives. For example, changes of an assigned runway might change the flight position in the sequence. Consequently, the separation minima and the headway time might change from the FCFS sequence. Figure 4 shows the genetic representation of a runway configuration. The genes of the chromosome are as follows:  Similarly, the runway configuration is estimated during the co-evolution. For each configuration, the throughput is projected for the selected aircraft sequence. When the genes of the chromosome are updated due to evolution, the runway allocation of the flight also changes accordingly. The evolved chromosome affects the throughput of the configuration and the selected aircraft sequence. After each evolution, the best sets of chromosomes are allowed to participate in the next generation. During the evolution of each generation, the constraints are continued in each change.

Experimental Design
The terminal airspace operation is illustrated in Figure 5. The Initial Approach Fix (IAF) is the point where the initial approach segment of an instrument approach begins. The last segment when coming into the airport is the final approach segment, which begins at the Final Approach Fix (FAF). The FAF is usually the outer marker on localizer approaches. For the arriving aircraft, it is assumed that all the landing aircraft sequence begins from the FAF; while the departure sequence commences from the runway threshold. It is also assumed that all the landing aircraft are sequenced at the FAF and all the departing aircraft sequenced at the runway threshold of the departure runway. The runway configuration might change any time based on the influencing factors or the ATC's strategic plan, e.g., noise abatement. However, in this paper, a fixed time interval (i.e., 30 min) is considered for measuring the efficiency of each configuration in terms of throughput. To well capture the real-world scenario of the busy airport, we have considered the traffic distribution and the consecutive arrival and departure at the aircraft. In this study, it is assumed that all the configurations and selected runways satisfy the FAA safety rules [5]. In the experiments, arriving and departing aircraft are stochastically generated based on the Poisson arrival process. The experiments are designed for a busy period of the airport when arrival and departure demand is high. In total, fifty aircraft (i.e., arriving and departing) are generated in the initial sequence. Fifty different runway configurations are generated randomly from the five possible orientations. The time horizon we consider for the simulation is 30 min. Figure 6 illustrates the runway layout of O'Hare International Airport where many runway configurations and five orientations are present. The O'Hare International Airport runway layout is used as a baseline input for runway configuration because all possible orientations are available in the runway layout. The configurations are generated assuming that cross-wind and tail-wind components do not exceed the maximum limit and satisfy Equations (16)- (19). For instance, 9L|9R represents an independent parallel runway configuration where runway centre lines are spaced by a minimum of 760 m. The experiments are calibrated by estimating the capacity of each available runway configuration and the inter-arrival-departure time of the combined aircraft sequence. The simulation parameters are presented in Table 3. The population size was chosen to be 50 from [64]. In the evolution of runway configuration, three operators are utilized (e.g., selection, crossover and mutation). The mutation rate is 0.035 for the runway configuration due to the small chromosome size.  The parameters of the evolutionary algorithms are tuned depending on the problem characteristics and encoded chromosomes. An elitism strategy is incorporated in each subproblem of the algorithm for constructing a new population. This allows the best organism from the current generation to carry over to the next generation [74,75]. It guarantees the solution quality by applying the recombination of elite individuals of the previous generation and the created new individuals of the current generation. For tuning the parameter and sensitivity analysis, multiple runs of the algorithms with different probabilities are carried out. A too small mutation rate causes unnecessary delay for the solution to converge, while a too large mutation rate causes a much faster convergence, which might find only the local optima. The genetic operator parameters are summarized in Table 3. The proposed methodology is developed and implemented using the Microsoft C# programming language, and the simulation is run on OSX 10.7.3, processor core i5, 1867 MHz, memory 8 GB.

Result Analysis and Discussion
The performance of the model was measured for the traffic scenario at a busy time. A stream of consecutive arriving and departing aircraft and all the possible runway orientations were applied in the model. The model integrated the real-world constraints and assumptions to capture the complex traffic scenario.

Optimal Aircraft Sequence
First, Figure 7 shows the evolution versus required process time per sequence. Three constraint techniques were applied to CCoGA for the combined solution. Multiple (20-times each) runs were conducted to observe the consistency of the model. The average simulation time up to 300 generations for each constraint position technique was 16.29 min. The three constraint position techniques were integrated into the CCoGA to understand the impact of position shift, i.e., 1-CPS, 2-CPS and N-CPS of the aircraft sequence. FCFS is the most commonly-used technique air traffic controllers use for sequencing the aircraft. In an FCFS, the aircraft lands or takes-off in order of the scheduled arrival or departure times at the runway. The ATC only enforces maintaining the minimum separation. Since the aircraft cannot change their position in the initial sequence of CCoGA-FCFS, the total inter-arrival-departure time of the sequence remained constant, as shown in Figure 7. A major disadvantage of the FCFS technique is that it may lead to reducing the overall capacity due to the large spacing requirements. At the busy airport, ATC would like to process the maximum number of aircraft, since long waiting at the arriving aircraft contributes to congestion and increases the huge workload of the ATC. Since ATC does have a certain degree of flexibility to shift the aircraft quickly, therefore another beneficial technique widely practised is 1-CPS, where an aircraft can shift one position left or right from the original sequence. The proposed model integrates the 1-CPS to observe the impact in combination with the best-fit runway configuration. Note that, by shifting only one position from the original (i.e., FCFS) sequence, the total process time decreased significantly. The solution converges very quickly in 1-CPS because aircraft can shift their position only one position.

Number of Co-Evolution
Similarly, we have been encouraged to make more room for shifting the aircraft. In 2-CPS, ATCs get more spaces to accommodate aircraft from the FCFS position. However, aircraft still have a restriction to move more than two positions. As a result, 2-CPS sequence obtained a slightly lower process time than 1-CPS. The solution converged slower than 1-CPS for the 2-CPS sequence.
N-CPS is a technique where aircraft can shift their position anywhere in the sequence. Apparently, it might be challenging for ATC; however, in the dynamic environment of terminal airspace, it is quite possible to merge the aircraft in any position, especially the arriving aircraft given that minimum safety is maintained. Furthermore, due to the invention of the point-merge technique with the Precision-Area Navigation (P-RNAV) and Continuous Descent Approach (CDA) technique, aircraft can merge any point of the final approach fix. Note in Figure 7 that, by integrating the N-CPS technique with CCoGA, a significant improvement was achieved.
Second, a time-space diagram is commonly used to show the trajectories of individual aircraft in time and in space. It describes the relationship between the position of aircraft in an arrival or departure stream and the time as the flight progresses along the runway or away from the runway. The obtained aircraft sequences from the different constraint position techniques are presented as the form of a time-space diagram. Figure 8 shows the FCFS sequence that follows the Poisson arrival process. Three different colours indicate the three different wake category of aircraft (H, M and S).
Since the FCFS sequence follows the random (i.e., Poisson arrival) process, the total process time for the sequence is 5696.72 s, which is the highest compared to the other sequences. The scattered point represents the arrival and departure throughput of the evolving aircraft sequence and the configuration. The solid lines in Figure 12 show the capacity envelope for each constraint case. Since the CCoGA 2-CPS technique has higher room to adjust aircraft position, the inter-arrival-departure time has reduced more compared to 1-CPS. However, the runway configuration plays a vital role in throughput, and consequently, many 1-CPS sequences obtained the same throughput as 2-CPS in Figure 12. Similarly, CCoGA N-CPS obtained a significantly higher throughput than 1-CPS and 2-CPS.
Another interesting point to note here is that there was a trade-off between arrival and departure throughput. When a configuration gets a large number of arrivals, at that point, departure capacity reduces. On the other hand, when the departure capacity increases, the corresponding arrival capacity decreases.

Best-Fit Runway Configuration
The choice of the runway configuration and the assignment of the runway for arriving and departuring flights are crucial decisions for ATC. The proposed model provides multiple solutions of runway configurations available for the mix sequence. Table 4 summaries the feasible runway configuration, corresponding orientation and the throughput that is achieved when that configuration is simulated. The throughput and the feasible configurations presented were obtained from the CCoGA N-CPS study. There was a total of 50 aircraft in the scenario, and the throughput was estimated at 30-min intervals. All the configurations generated from the model were feasible with maximum throughput capacity. Table 4. Feasible runway configuration and throughput: 30 best runway configurations that accommodate the maximum numbers of arriving and departing flights.

Configuration
Orientation Throughput CV  DP  IR  IR  IR  IR  30  55  85  9R, 9C | 4R, 10C  DP  CV  IR  IR  IR  IR  30  55  85  28R, 4R | 27C, 9R  CV  DP  IR  IR  IR  IR  30  55  85  28R, 10C | 27C, 9R  DP  DP  IR  IR  IR  IR  30  55  85  28R, 10C | 28R, 9C  DP  IR  IR  IR  IR  IR  30  44  Note the impact that the runway orientation has on the capacity of runway configuration. Independent parallel runways obtained the highest throughput because runways were mutually independent for operation. When two dependent parallel runways operated simultaneously, the arrival must be at least 1.5 nmi apart during the final approach. One runway had priority (e.g., only arrival) and was called the primary runway, and another one had lower priority (e.g., only departure) and was called the secondary runway. In such a primary/secondary combination runway, capacity reduces compared to independent parallel runways. Likewise, the configuration with cross orientation obtained an equal number of arrivals and departures as the throughput. Cross runways were inter-dependent; therefore, when a runway was used, the corresponding cross-runway could not be used until the first cross runway was free. Therefore, ATC had a huge scope to choose the best combination according to the requirements as the proposed model provided many alternative options. Here, the thirty best feasible configurations were presented from the final evolution that outlined a trade-off among the configurations in terms of runway throughput.
The proposed model provided a trade-off solution considering aircraft sequence and runway configuration, meaning the optimised sequence and runway configuration active at a given time horizon that can provide the maximum throughput capacity. Figure 13 illustrates the obtained configurations versus runway throughput. We notice that the proposed model achieved a better  It also shows the increase in runway throughput by inserting departures between consecutive arrivals. Note that different arrival runway configurations accommodate a significant number of departing aircraft. However, few configurations cannot accommodate any departing aircraft because the required runway occupancy time for departure is higher than the inter-arrival gap. Therefore, a decision can be made to utilize the best-fit configuration based on the demand and the flight assignment to the designated runway.

Conclusions
In this paper, we proposed a cooperative co-evolution genetic algorithm that can maximize the runway throughput capacity by optimising arrival and departure sequences with best-fit and feasible runway configurations under constraints. The proposed method incorporates multiple objectives such as on-air tactical strategies (i.e., mixed aircraft sequencing) and the best-fit ground resources (i.e., runway configuration) to obtain a trade-off solution from alternative sequencing strategies. The proposed model is able to take into account the terminal airspace constraints (i.e., safe separation, safe airspeed and suitable runway assignment) and the performance indicators (runway capacity).
The proposed approach uses the detailed modelling of inter-arrival-departure time and runway occupancy time to estimate the headway between two consecutive runway operations. It also uses the best-fit runway configuration according to the mixed sequence in terms of capacity. It deals with the coordinated management of aircraft mix sequence (arrival and departure), suitable runway assignment and the best-fit runway configuration that can achieve the trade-off solutions together. The trade-off solution means the optimised sequence and runway configuration active at a given time horizon can provide the maximum throughput capacity.
The simulation results show that the optimal sequence for a configuration significantly reduces the inter-arrival-departure time, i.e., an increased number of slots to allocate more aircraft, thereby increasing the runway capacity. Time-space diagrams for arrival and departure sequencing under different sequencing strategies given a runway configuration might be an effective visual decision support tool for ATC. The arrival-departure capacity envelope is then presented to illustrate the trade-off between the arrivals and departures, given a runway configuration. Results demonstrate the mutual dependence between arrival-departure sequence and runway configuration, as well as its effect on runway capacity. The results also demonstrate the viability of using the Co-operative co-evolution approach for maximizing the runway throughput capacity with consideration of runway configuration. N-CPS sequencing (which allows ATCs to shift an aircraft n-positions in a sequence) with an independent parallel runway may give high throughput in the above scenarios, but is highly unlikely in the real world given the workload it entails on an ATC. However, with the advent of point-merge approaches (where aircraft are not required to line up and lock the sequence 10 nm from the runway), this is no longer a constraint.
As future work, we will be extending the cooperative co-evolutionary method by considering wider traffic and runway constraints such as taxiway-exits, noise restriction and a heterogeneous traffic environment where diverse air vehicles, such as UAVs, may operate.

Conflicts of Interest:
The authors declare no conflict of interest.