，，，，

1.College of Civil Aviation，Nanjing University of Aeronautics and Astronautics，Nanjing 211106，P.R.China；2.Beijing Capital International Airport Co.，Ltd.，Beijing 100621，P.R.China

Abstract:This paper proposes an optimization model for the airport ground movement problem（GMP）based on bilevel programming to address taxi conflicts on the airport ground and to improve the operating safety and efficiency.To solve GMP，an iterative heuristic algorithm is designed.Instead of separately investigating each problem，this model simultaneously coordinates and optimizes the aircraft routing and scheduling.A simulation test is conducted on Nanjing Lukou International Airport（NKG）and the results show that the bilevel programming model can clearly outperform the widely used first-come-first-service（FCFS）scheduling scheme in terms of aircraft operational time under the precondition of none conflict.The research effort demonstrates that with the reduced operating cost and the improved overall efficiency，the proposed model can assist operations of the airports that are facing increasing traffic demand and working at almost maximum capacity.

Key words：airport ground movement（GMP）；aircraft routing and scheduling；bilevel programming；iterative heuristic；air transportation

0 Introduction

Airports are experiencing ever more congestions on ground movements，due to the long-term growth of air transportation.It has increased air traffic controllers’workload，safety risks，operating costs and decreased passenger comfort.More noise pollution and exhaust emissions have also impeded airports’sustainable growth.Currently，increasing the throughput capacity and efficiently using existed infrastructure are two major methods in airport surface operations to relieve congestions［1］.As the former is always restricted by stringent land use and financial investment［2］，optimizing the usage of airport terminal areas has become the leading approach.Airport surface resources consist of gates，taxiways and runways.As taxiways are the physical links between gates and runways，optimization of ground movement on taxiways directly relates to other problems，like gate assignment and runway sequencing［3-5］.Therefore，the taxiway operation is crucial to the holistic airside efficiency.

Atkin et al.［6］attributed the airport ground movement problem（GMP）to a fundamental routing and scheduling problem.An aircraft is instructed to taxi from the source to its destination following a specific path（taxi routes）in a timely manner，e.g.pushback delays，holding patterns and speed profiles.In the case where only a few aircraft need to be dealt with at one time，a shortest path algorithm，like Dijkstra’s algorithm or A*，can be used to generate taxi routes.For larger hub airports with more intensive traffic flows on the ground，an alternative routing and scheduling strategy is required to avoid conflicts between moving aircraft as well as to increase overall efficiency.

Some studies have dealt with GMPs using the mixed integer linear programming（MILP），where objective functions and constraints are formulated as linear expressions with the presence of both integers and continuous variables.Time variables are either discretized through an exact position approach，usually under a time-space network model［7］，or described as continuous ones when binary variables are used to represent the ordering of passing taxiway intersections or sectors［8］.In Ref.［3］，a comprehensive MILP model were incorporated with the receding horizon scheme in conformance with the dynamic nature of airport GMPs，and an iterative algorithm was designed to fasten the computation.As a GMP is essentially an NP-hard problem，the enlargement of the problem scale will dramatically increase the complexity of the algorithm.Other practices in MILP for GMPs can be found in Refs.［9-11］.

Considering the difficulty of developing an exact algorithm in terms of MILP formulation，some researchers adopted either heuristics or metaheuristics in searching for an acceptable approximate optimization rather than a perfect theoretical one.These approaches usually compensate for the loss of optimality by relatively short time span.Nogueira et al.［12］used the ant colony algorithm to optimize taxi routes of aircraft.They searched optimal taxi routes and minimized taxi time in accordance with safety constraints.Weiszer et al.［13］presented an integrated optimization approach for airport ground operation.They adopted an active routing（AR）strategy in the model to obtain the shortest taxi time.Other studies involving the application of heuristics in GMPs can be found in Refs.［14-18］.

In recent studies on airport GMPs，multiple aspects of interest have been considered，including the incorporation with runway sequencing［3-4，19-20］，pollution emission［1，21］and the issue of uncertainty［22-25］.These efforts essentially focused on the routing and scheduling problem，despite their different research scopes.

To our best knowledge，a critical issue that has not been adequately addressed in this domain is the integrating and coordinating both routing and scheduling problems simultaneously，instead of previous attempts on only one problem.Biased focus on either side may lose the global optimality of the overall GMP and the airport airside capacity is potentially decreased.Therefore，we propose a bilevel programming model dealing with routing and scheduling simultaneously from the collaborative perspective.To the best of our knowledge，no other research has introduced the use of bilevel programming on GMP though this method has been applied to other transport optimization problem，as illustrated in Refs.［26-29］.An iterative heuristic is then designed to solve the problem in an efficient way.

1 Problem Description and Bilevel Modeling

1.1 Ground movement problem

As aforementioned，airport GMPs can be defined as the work of assigning taxi routes and the timings of arriving at taxiway sectors or intersections for both arrival and departure flights.The running procedure for an arrival flight is taxiing from runway exit to a predetermined gate following a specific path，while a departure flight taxies from the gate to an assigned runway threshold.In the taxiing process，a prescribed safety separation should be strictly maintained between any two aircraft in order to avoid conflicts.Fig.1 demonstrates the taxischeduling flow in a large airport.

Fig.1 Operation flow on airport ground

In order to clarify the research category and standardize the research content，certain assumptions and simplifications are defined as follows，though it may compromise the fidelity of the model.

（1）Airport ground is abstracted to a directed graphG=(V，E)，whereVandErepresent nodes and edges，respectively.Any edge(v1，v2)∈Econnecting nodesv1andv2has a direction，which means aircraft are only allowed to taxifromv1tov2.

（2）In order to constrain the problem scale within a manageable range，it is assumed that the taxi route for each aircraft is selected from a set of predetermined feasible routes.

（3）Aircraft ground taxi is a continuous operation process，which involves massive operating data.In this paper，the ground taxiresearch is conducted within a busy period during peak hours.

（4）The moving speed of an aircraft is assumed to be constant during taxiing.The process is continuous and unimpeded，which means any stopping and waiting is not allowed.Taxi conflicts are avoided by controlling the time of entering taxiway and selecting the taxi routes.

（5）Some adjacent stands are integrated as an apron.This model only requires aircraft to be assigned to certain aprons.The number of stands in each apron is sufficient for all the departure and arrival aircraft within the planned period，and the stand for each aircraft is assigned in advance.

（6）The departure aircraft taxies from the stand to the runway threshold，and the arrival aircraft taxies from runway exit to the stand.Both the runway threshold and the exit for each aircraft are fixed.

1.2 Bilevel pr ogr amming model

In the ground movement bilevel programming model，the decision variable at the upper level is the start time of taxiing.By controlling the start time of taxiing，it aims at reducing conflicts during taxiway operation and aircraft’s waiting time before taxiing.The decision variable at the lower level is the taxi routes assigned for aircraft.The taxi conflicts and travel distance can be reduced after selecting proper routes.Notations used in this model are presented in Table 1.

Table 1 Notations

The objective of the upper-level programming is to minimize the total waiting time and the number of conflicts for all arrival and departure aircraft.The decision variable is the start time of taxiingTsi.During the scheduling process，the targets are minimizing waiting time of arrival and departure aircraft and lowering the ground operation costs，as well as resolving conflicts and improving efficiency.Departure aircraft are supposed to start taxiing as early as possible to reduce the delay at stands.Meanwhile，the arrival aircraft should land as early as possible to minimize the costs of waiting and impacts on followup flights.Thus，the objective function and constraints at the upper level are formulated as follows

Eq.（1）is the objective function at the upper level.Since the cost of waiting is higher in the air than in the stands，the waiting cost of arrival aircraftcais set to be greater than that of departure aircraftcd.δstands for the number of taxiing conflicts happened，andMis a big positive number used as penalty for conflicts.In Eq.（2），the aircraft taxi routes calculated at the lower level are transferred to be a constraint at the upper level.Eq.（3）ensures the departure aircraft must not start taxiing until the estimated pushback time.Eq.（4）ensures the departure aircraft must enter the runway prior to the estimated departure time.Eq.（5）ensures the arrival aircraft must start taxiing after the estimated arrival time.

The objective of the lower-level programming is to minimize the total taxi time of all aircraft by properly assigning taxi routes for the aircraft.The decision variable is the feasible taxi route for aircraftf i，N i=.The aircraft is supposed to start taxiing according to the time passed from the upper level.The objective function and constraints at the lower level are formulated as follows

Eq.（6）is the objective function at the lower level，whereis the time of aircraftfiarriving at the first node，andthe time of aircraftfiarriving at the last node.Eq.（7）transfers the start time of taxiing at the upper level as a constraint at the lower level.For each aircraft，Eq.（8）ensures the taxi time from one node to another is shorter than or equal to the time difference between their arrivings，which reflects the constraint of the continuous operation.Eq.（9）is used to detect conflicts between two aircraft in any node or taxiway chain，in whichτijpis the time-equivalent safety interval andτijp=(T iq-T ip)/L pq·dsep.Eqs.（10—13）are used to detect rearend conflicts and head-on conflicts.Any dissatisfaction for the inequalities in Eqs.（9—13）will lead to an increasing number of conflict pointsδ.

2 Iterative Solution Method

To solve the bilevel programming model，the first step is to initialize the decision variable at the upper level.Each aircraft is randomly given an initial solution within reasonable ranges，that is，the start time of taxiing for both departure and arrival aircraft is provided.Then at the lower level，within the constraints of this parameter，proper taxi routes are scheduled to optimize the objective function.The solution is fed to the upper level.Afterwards，the upper level makes the better decision once again based on the optimized solution from the lower level and returns the result to the lower level.This process will be repeated till the stopping criterion is fulfilled.Fig.2 illustrates it in detail.

Fig.2 Flow chart of solving process

As the genetic algorithm（GA）has the advantages of less computation，rapid convergence and few parameters，this paper imitates the iterative and evolutionary process of it.A solution to ground scheduling problem based on bilevel programming is shown as follows.

Step 1Initialization.Suppose generationg=0，and the waiting time at the upper level is randomly assigned within reasonable ranges.Then the total waiting timeTDis calculated and the waiting time solution{T1，…，T i，…，T n}andTDare passed to the lower-level programming.To initiate taxi routes at the lower level，routesR iare randomly selected for each aircraft under the constraints of waiting time{T1，…，T i，…，T n}，and then the taxi routes solution{R1，…，R i，…，Rn}are formed.Finally，the total travel distanceLis calculated.

Step 2The time arriving at each node in the selected route is calculated and the number of conflictsδis detected through constraints（9—13）.The routes solution，total travel distance and the number of conflicts are then handed over to the upper level.

Step 3Suppose generationg=g+1.The upper level duplicates the waiting time solution inwgroups and adjust the time in each group by introducing a random numberτ∈(-10，10).Eq.（14）limits the variation range of the waiting time.The waiting timeT iof each aircraftfiwill mutate with a specified probabilityσ

Step 4The total waiting timeT wDis calculated for each group.Then，based on the time solution，as well as combined with latest route solution，the upper level calculates the time arriving at each node in the route and detects the number of conflictsδw.

Step 5Among the duplicated groups in which≤TDandδw≤δ，the optimal group is selected as the new waiting time solution at the upper level.The total waiting time and the number of conflicts in this group are passed to the lower level.

Step 6The lower level duplicates the taxi routes solution inwgroups and adjusts each route with a specified probabilityσ.In order to reduce the computing complexity，the alternative taxi routes are selected from a predetermined feasible route set.

Step 7The lower level calculates the total travel distanceL wfor each group.Then，the time arriving at each node in the selected route are obtained based on each route solution and relations with the waiting time from the upper level.Through this approach，the number of conflictsδwis detected.

Step 8Among the duplicated groups in whichL w≤Landδw≤δ，the lower level selects the optimal group as its latest taxi route solution.Similarly，the routes，the total travel distance and the number of conflicts are returned to the upper level.

Step 9Ifg≥gen，where “gen” means the times of generation，the loop ends；otherwise goes back to Step 3.

Step 10If the number of conflictsδ=0 at this time，the solutions to taxi routes and waiting time are taken as the final result of the bilevel programming.Otherwise the solving process fails.

3 Computational Test

In this paper，the ground movement of Nanjing Lukou Airport（NKG）is chosen as a studying objective.We select part of the airport and abstract several adjacent stands as three integrated aprons G1，G2 and G3.After omitting some irrelevant taxiways and nodes，the taxiway layout is drawn with 37 nodes，48 edges，two runways and three integrated aprons，as shown in Fig.3，where TML means terminal.

Fig.3 Layout of NKG taxiway network layout

The flight schedule data are chosen between 8：00 and 8：15 during rushing hours，as listed in Table 2.The estimated pushback and landing time of departure and arrival flights are converted into the start time of taxiing.The estimated departure time is then converted into 5 min after starting taxiing.σis set to be 0.4 andwis set to be 100.

Table 2 Aircraft timetable

4 Result Analysis

4.1 Simulation results

Simulation is conducted on MATLAB software.We obtain the optimal solution=186 s at the upper level，and=43 250 m at the lower level.The evolutionary processes of optimal solutions at both upper and lower levels are illustrated in Fig.4.For the curve of the upper level programming（ULP），we can find that with the number of iterations increasing，the value of total waiting time declines from the large initial value.The first five generations have a rapid descent as the initial solution has been randomly given a large value，resulting in a higher number of conflicts and a function value.Then，the total waiting time keeps going down as iteration moves forward.The total waiting time starts leveling off since the 90th generation，reaching the optimal solution within the maximum number of iterations.

Fig.4 Optimal solution evolution process at the upper and the lower levels

From the evolutionary curve of the lower level programming（LLP），the total travel distance shows a stepped downward trend as iteration goes further.The value decreases generally faster during the first 100 generations.This phase gains better evolutionary efficiency and a faster decrease of the total length.Whereas only minor reductions appear during the remaining，which is probably related to the increasing difficulty of searching for a better solution as iteration moves forward.

Fig.5 illustrates the evolutionary process of the number of conflicts and the relevant LLP object function value during the 200 generations.As shown in Fig.5，there are six conflict points initially.The number of conflicts decreases to zero until the 113th generation.The result corresponds with the general evolutionary trend of the object function value.Since conflicts are strictly prohibited in the designed ground movement model，only conflicts-free solutions can be considered as feasible，i.e.，solutions after the 113th generation.We can conclude that it is adequate to set 200 as the maximal number of iterations as the number of conflicts stalls for over 80 iterations.The setting ensures computational efficiency and the acquirement of stable optimal solutions.

Fig.5 Evolution process of conflict points and total travel distance

Fig.6 illustrates the taxiing and the waiting time of each aircraft in the optimal result.It is indicated from the travel distance that all the aircraft successfully taxi to their assigned destinations and have balanced distances.Every waiting aircraft calls for rather short waiting time，which is within the acceptable limit.Only No.4 and No.6 need to wait for 106 s and 80 s，respectively，accounting for 12.5%of the total.The flight schedule only suffers a minor disturbance.It is proved by simulation results that the total time can be reduced substantially at the cost of a little extra waiting time.It helps to improve taxiing efficiency and lower operation cost.

Fig.6 Total time of the aircraft

We perform simulation by the first-come-firstservice（FCFS）algorithm on the same flight date.Table 3 compares the simulation results.

Table 3 Result comparison between the FCFS scheme and the bilevel programming

In Table 3，for the FCFS scheme，the total travel distance of the 16 aircraft is 35 330 m，and 2 208.1 m on average.As aircraft are operated according to the flight plan，there is no waiting time.In contrast，for the bilevel programming，the total travel distance increases by 2.0%to 36 050 m，and 2 253.1 m on average.The total waiting time is 186 s，and 11.6 s per aircraft，which are acceptable.The intention of avoiding conflicts inevitably leads to the presence of waiting time in the bilevel programming.As a consequence，the taxiing efficiency and the safety performance are improved，and there is no conflict.Although the FCFS scheduling can avoid waiting when the aircraft pushes back，the total wasted time is much higher than that of the bilevel programming，as more time is spent on avoiding conflicts on the taxiway.It is anticipated conservatively that each conflict takes 30 s to be resolved.Given the assumption that the aircraft moves on at a constant speed of 10 m/s，the operational time is the sum of the travel time，the waiting time and the time for solving the conflicts.It turns out the bilevel programming surpasses the FCFS scheme in terms of the average operational time，as illustrated in Table 3.

Fig.7 further compares the operational time of each aircraft in these two methods.Scatters above the dashed line（slope equals to 1）represents the aircraft that have a longer operational time in FCFS than in the bilevel programming，or vice versa.The operational time in the bilevel programming is generally shown to be shorter due to its significant improvement in avoiding conflicts.Only two aircraft have an increased operational time after optimization.Alhough it may need further investigation afterwards，a possible explanation would be that their taxi routes are relatively distant from the main traffic flow on the ground，by which they are less likely to conflict with other moving aircraft.

Fig.7 Aircraft operational time comparison between FCFS and the bilevel programming

From the above comparisons，we can conclude that the bilevel programming outperforms the FCFS scheduling in terms of the total travel distance，the taxi time and the number of conflicts.This is the because the bilevel programming regularly feedbacks the routes and the waiting time，optimizing the taxi routes as well as the reducing taxi conflicts in both the spatial and the temporal ways.For example，spatially，the aircraft selects a relatively short route and taxies to the destination without passing by congested areas；or temporally，the aircraft chooses the shorter waiting time and then starts to taxi after conflicts.In contrast，the FCFS scheduling is unable to select shorter routes，or to avoid conflicts except waiting.In conclusion，the bilevel programming is superior to the FCFS scheduling in terms of the total travel distance，the taxi time and the number of conflicts.

4.2 Analysis of heuristic par ameter s

Looking at the solving method in the bilevel programming，the value of mutation probabilityσand the number of duplicationswhave a significant impact on the program execution and results.The following part analyses their impacts on the simulation results.

These two parameters are studied through controlling variables.First，we hold the number of duplicationsw=100 and change the value of mutation probability from 0.2，0.4 to 0.6.Ensuring all other conditions being equal，we repeat the simulation test and average the test values.Some representative simulation results are chosen for the analysis，as shown in Table 4.

Table 4 Impacts of mutation pr obability

In Table 4，the results are very close when the mutation probability equals to 0.2 and 0.4，while the average travel distance and the average waiting time increase to some extent as the mutation probability rises to 0.6.The program executing time is very similar among these three test groups，which is acceptable.Whereas the most remarkable impact of the mutation probability on this method is convergence.When the mutation probability is 0.2，the algorithm presents a slow convergence.It is not until the 180th generation that the optimal solution is found，or even no conflicts-free solutions can be acquired at the end.When the mutation probability is 0.6，the optimal solution appears at the 10th generation，while too fast convergence tends to be trapped into local optimum and lacks stability.Finally，the algorithm convergence lingers when the mutation probability is 0.4.It has a moderate convergence while the speed is not too slow to find a satisfactory optimal solution within the number of iterations set before.From the above analysis，the result is better when the mutation probability equals to 0.2 or 0.4，rather than 0.6.

Multiple test results illustrate that when the mutation probability is too small，the program will have a slow convergence，or even no conflicts-free solutions can be acquired within the maximum generation set in advance.This is caused by a lower possibility of generating new excellent individuals through mutation.Alternatively，the crossover operation within the population can hardly produce more good individuals，leading to the algorithm’s slowness of convergence and weakness in solving ability.In addition，a too large mutation probability may enhance the randomness of individuals within the whole population，resulting in early convergence and local optimum.The movement plan is not fully optimized as well.

Similarly，we hold the mutation probabilityσ=0.4 and assignwwith 50，100，and 200.As all other conditions being equal，we repeat the simulation tests and average the test values.The results are shown in Table 5.

Table 5 Impacts of number of duplications

In Table 5，the aircraft’s average travel distance and average waiting time show a downward trend with the number of duplications increasing.Additionally，the program executing time has a positive correlation with the number of duplications and increase linearly.

The increase of the number of duplications enlarges the population scale，so better solutions are more likely to be acquired.The aircraft’s average travel distance tends to decrease with the rise of the number of duplications.However，it does not mean “bigger is better” due to the demerits of larger data size and more time spent on computation.So often it takes a fair amount of time to make a little progress.In practical use，it is recommended to select a value in the middle after several tests and comparisons，and to find the optimal results with acceptable executing time.

5 Conclusions

This paper studies optimization of GMPs and discusses the integration of both routing and scheduling.By introducing the bilevel programming，an optimization model is established.An iterative heuristic algorithm is designed for the model.Simulation results indicate that the solving process has good convergence and executing efficiency，and is able to optimize taxi routes and waiting time simultaneously.As the application of the bilevel programming to GMPs is rarely studied at present，research results in this paper can be theory guidance to its future development.The solving process of the bilevel programming is a complete NP-hard problem.When the scale of the problem is enlarged，computational time will significantly increases as well.So the algorithm used in this paper may not be capable of solving the problem anymore.In conclusion，research direction in the future should focus on finding an efficient globally optimal algorithm for the large-scale bilevel programming.