Acessibilidade / Reportar erro

A NEW MODELLING FOR THE AIRCRAFT LANDING PROBLEM AND MATHEURISTIC APPROACH TO SOLVE THE PROBLEM WITH A LARGE NUMBER OF AIRCRAFT

ABSTRACT

Air traffic management has become increasingly complex due to the increasing use of air transport. One of the main management bottlenecks is planning the efficient use of runways for takeoff and landing. This paper aims to investigate the Aircraft Landing Problem, which seeks to minimize earliness and tardiness in aircraft landing time, assigning them to a runway to land and sequencing them. A new mathematical formulation based on Job Shop was proposed for the problem, comparing it with four mathematical formulations in the literature; three directly comparable and another containing a particularity that does not allow a direct comparison with the other formulations. Computational tests were performed on 49 instances of the literature using the Gurobi Optimizer optimization package. These mathematical formulations commonly used for the ALP present difficulties in finding the optimal solution when the number of aircraft to land is large, i.e., more than 50 aircraft. Therefore, we proposed a matheuristic to solve instances with a greater number of aircraft than the Gurobi Optimizer cannot solve optimally. This matheuristic first finds an initial solution using relax-and-fix (RF) and then fix-and-optimize (FO) improves the found solution. Comparisons were also made using the first feasible solution obtained by Gurobi and then was improved with FO. Among the matheuristic variations, the one that obtained the best result was the combination of RF with FO and this also showed efficiency in relation to the work in the literature that uses FO.

Keywords:
aircraft landing problem; mixed integer programming; matheuristic; relax-and-fix; fix-and-optmizer

1 INTRODUCTION

Air transport is considered one of the safest forms of transportation to carry cargo and people. The first commercial flight took place on the first day of 1914, approximately 7 months before the beginning of the first world war. This flight, which took place in the United States, left the city of St. Petersburg, Florida, towards Tampa, in the same state, covering 34 kilometers in 23 minutes. At that time, people perhaps did not even imagine that one day there could be congestion of commercial flights that needed the optimization of a sequence of landings and takeoffs.

In 2019, Beijing International Airport in China received more than 100 million passengers, making it the second busiest airport in the world1 1 http://english.www.gov.cn/archive/statistics/202003/11/content WS5e6861e1c6d0c201c2cbe089.html . In 2019, São Paulo International Airport, in Guarulhos, the largest Brazilian airport, received just over 43 million passengers 2 2 https://www.gru.com.br/pt/institucional/informacoes-operacionais/movimentacao-aeroportuaria .

As a consequence of the high and increasing volume of passengers and flights, airport runways became one of the main bottlenecks for air traffic management and one of the main factors determining airport capacity. The construction of new runways to decongest airports is not always available due to high investment costs. Therefore, the alternative is to optimize existing infrastructures.

The Aircraft Landing Problem (ALP) can be roughly defined as the assignment of the runway for an aircraft to land optimizing landing time and minimizing costs due to earliness and tardiness. Solutions to the ALP must find sequences of aircraft’s landing respecting a time window, i.e., a time interval to perform their landing. Furthermore, the solution must take into account the separation time between the landing of one aircraft and the landing of the next. In the ALP, one or more runways (single or multiple runways) can be considered. According to Salehipour & Ahmadian (2017SALEHIPOUR A & AHMADIAN M. 2017. A heuristic algorithm for the Aircraft Landing Problem. In: 22nd International Congress on Modelling and Simulation. pp. 1344-1349. Hobart - Austrália.), the computational effort to find the optimal solution becomes smaller when using multiple runways.

The ALP is considered an NP-Hard (Girish, 2016GIRISH BS. 2016. An efficient hybrid particle swarm optimization algorithm in a rolling horizon framework for the aircraft landing problem. Applied Soft Computing, 44: 200 - 221.) problem, i.e., the computational effort to find an optimal solution grows exponentially with the number of aircraft. To deal with this limitation, several methods to find efficient solutions have been proposed in the literature, ranging from exact to meta-heuristic methods.

The classic formulation for the ALP was proposed in Beasley et al. (2000BEASLEY JE, KRISHNAMOORTHY M, SHARAIHA YM & ABRAMSON D. 2000. Scheduling Aircraft Landings-The Static Case. Transportation Science, 34(2): 180-197.), considering single and multiple runways and optimally solving instances with up to 50 aircraft. Since then, variations have emerged in both the constraints and the objective functions. In Beasley et al. (2000BEASLEY JE, KRISHNAMOORTHY M, SHARAIHA YM & ABRAMSON D. 2000. Scheduling Aircraft Landings-The Static Case. Transportation Science, 34(2): 180-197.), Salehipour et al. (2013SALEHIPOUR A, MODARRES M & NAENI LM. 2013. An efficient hybrid meta-heuristic for aircraft landing problem. Computers & Operations Research , 40(1): 207 - 213.) and Faye (2015FAYE A. 2015. Solving the Aircraft Landing Problem with time discretization approach. European Journal of Operational Research, 242(3): 1028 - 1038.), the authors attempted to minimize earliness and tardiness for multiple runways, with 10-50, 10-50 and 10-44 aircraft, respectively. More recently, in Ikli et al. (2020IKLI S, MANCEL C, MONGEAU M, OLIVE X & RACHELSON E. 2020. Coupling Mathematical Optimization and Machine Learning for the Aircraft Landing Problem. In: ICRAT 2020, 9th International Conference for Research in Air Transportation. ICRAT 2020. Tampa, United States.), the authors carried out experiments for single runways with 10-44 aircraft. In other articles, for single runways, Beasley et al. (2001BEASLEY JE, SONANDER J & HAVELOCK P. 2001. Scheduling aircraft landings at London Heathrow using a population heuristic. Journal of the Operational Research Society, 52(5): 483- 493.) and Ji et al. (2016JI XP, CAO XB & TANG K. 2016. Sequence searching and evaluation: a unified approach for aircraft arrival sequencing and scheduling problems. Memetic Computing, 8: 109-123.), additional expenses are considered, such as indirect costs caused by landing before or after the target landing time. In the first, the authors used 10-44 aircraft and in the second they used only one instance, with real data, for 20 aircraft, from Heathrow airport.

Exact methods, such as mixed-integer programming (MIP) ((Beasley et al., 2000BEASLEY JE, KRISHNAMOORTHY M, SHARAIHA YM & ABRAMSON D. 2000. Scheduling Aircraft Landings-The Static Case. Transportation Science, 34(2): 180-197.), (Faye, 2015FAYE A. 2015. Solving the Aircraft Landing Problem with time discretization approach. European Journal of Operational Research, 242(3): 1028 - 1038.), (Prakash et al., 2018PRAKASH R, PIPLANI R & DESAI J. 2018. An optimal data-splitting algorithm for aircraft scheduling on a single runway to maximize throughput. Transportation Research Part C: Emerging Technologies, 95: 570-581.)) and dynamic programming ((Balakrishnan & Chandran, 2006BALAKRISHNAN H & CHANDRAN B. 2006. Scheduling Aircraft Landings Under Constrained Position Shifting, .), (Briskorn & Stolletz, 2014BRISKORN D & STOLLETZ R. 2014. Aircraft landing problems with aircraft classes. Journal of Scheduling, 17(1): 31-45.), (Lieder et al., 2015LIEDER A, BRISKORN D & STOLLETZ R. 2015. A dynamic programming approach for the aircraft landing problem with aircraft classes. European Journal of Operational Research , 243(1): 61-69.)) have been used to solve ALP. Other optimization methods employed for such are meta-heuristics, mainly Genetic Algorithms (GAs), Ant Colony Optimization (ACO) and Local Search (LS). GAs-based solutions are more common, and can be found in Stevens (1995STEVENS G. 1995. An Approach To Scheduling Aircraft Landing Times Using Genetic Algorithms. Tech. rep.. Department of Computer Science. RMIT University.), Beasley et al. (2001BEASLEY JE, SONANDER J & HAVELOCK P. 2001. Scheduling aircraft landings at London Heathrow using a population heuristic. Journal of the Operational Research Society, 52(5): 483- 493.) and Cheng et al. (1999CHENG VHL, CRAWFORD LS & MENON PK. 1999. Air traffic control using genetic search techniques. In: Proceedings of the 1999 IEEE International Conference on Control Applications (Cat. No.99CH36328), vol. 1. pp. 249-254.). Regarding the other two metaheuristics, examples of their use can be found in Bencheikh et al. (2009BENCHEIKH G, BOUKACHOUR J, KHOUKHI F & ALAOUI A. 2009. Hybrid method for aircraft landing scheduling based on a Job Shop formulation. International Journal of Computer Science and Network Security, 9: 78 - 88.) (ACO) and in Soomer & Franx (2008SOOMER M & FRANX G. 2008. Scheduling aircraft landings using airlines preferences. European Journal of Operational Research , 190(1): 277-291.) (LS). Other solutions based on matheuristics were published in Vadlamani & Hosseini (2014VADLAMANI S & HOSSEINI S. 2014. A novel heuristic approach for solving aircraft landing problem with single runway. Journal of Air Transport Management, 40: 144-148.), Salehipour & Ahmadian (2017SALEHIPOUR A & AHMADIAN M. 2017. A heuristic algorithm for the Aircraft Landing Problem. In: 22nd International Congress on Modelling and Simulation. pp. 1344-1349. Hobart - Austrália.) and Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.). Finally, Veresnikov et al. (2019aVERESNIKOV G, EGOROV N, KULIDA E & LEBEDEV V. 2019a. Methods for Solving of the Aircraft Landing Problem. I. Exact Solution Methods. Automation and Remote Control, 80: 1317-1334.), Veresnikov et al. (2019bVERESNIKOV G, EGOROV N, KULIDA E & LEBEDEV V. 2019b. Methods for Solving of the Aircraft Landing Problem. II. Approximate Solution Methods. Automation and Remote Control , 80: 1502-1518.) and Ikli et al. (2021IKLI S, MANCEL C, MONGEAU M, OLIVE X & RACHELSON E. 2021. The aircraft runway scheduling problem: A survey. Computers & Operations Research, 132: 105336.)) provides an extensive review of the main approaches and methods to solve instances of the ALP.

Vadlamani & Hosseini (2014VADLAMANI S & HOSSEINI S. 2014. A novel heuristic approach for solving aircraft landing problem with single runway. Journal of Air Transport Management, 40: 144-148.) presented a matheuristic based on the Adapted Large Neighborhood Search (ALNS) to solve the single-lane ALP. The authors divided the ALP into two subproblems, a scheduling subproblem and a feasibility subproblem. The initial solution was generated randomly and in the first subproblem, the ALNS was applied finding an almost optimal landing solution. Afterwards, the solution found was used by the second subproblem and CPLEX was used to verify the feasibility. Computational tests were performed on instances of Beasley (1990BEASLEY JE. 1990. OR-Library: Distributing Test Problems by Electronic Mail. The Journal of the Operational Research Society, 41(11): 1069-1072.) from 10 to 150 aircraft and compared with the results of Salehipour et al. (2013SALEHIPOUR A, MODARRES M & NAENI LM. 2013. An efficient hybrid meta-heuristic for aircraft landing problem. Computers & Operations Research , 40(1): 207 - 213.). The approach found an optimal solution for all instances, and the comparative work did not find an optimal solution for instances with 50 and 100 aircraft. Using the hybrid meta-heuristic of Salehipour et al. (2013SALEHIPOUR A, MODARRES M & NAENI LM. 2013. An efficient hybrid meta-heuristic for aircraft landing problem. Computers & Operations Research , 40(1): 207 - 213.) to obtain an initial landing solution, Salehipour & Ahmadian (2017SALEHIPOUR A & AHMADIAN M. 2017. A heuristic algorithm for the Aircraft Landing Problem. In: 22nd International Congress on Modelling and Simulation. pp. 1344-1349. Hobart - Austrália.) solved the single-lane ALP. After this initial solution, the authors applied an LS that allowed a percentage of the aircraft to change their initial position. After having fixed this landing sequence, the CPLEX was used. The computational experiments were performed on the instances of Beasley (1990BEASLEY JE. 1990. OR-Library: Distributing Test Problems by Electronic Mail. The Journal of the Operational Research Society, 41(11): 1069-1072.) and compared with the heuristics of Pinol & Beasley (2006PINOL H & BEASLEY J. 2006. Scatter Search and Bionomic Algorithms for the aircraft landing problem. European Journal of Operational Research , 171(2): 439 - 462.) and Salehipour et al. (2013SALEHIPOUR A, MODARRES M & NAENI LM. 2013. An efficient hybrid meta-heuristic for aircraft landing problem. Computers & Operations Research , 40(1): 207 - 213.). The results showed that the approach in question was better for 10 of the 13 analyzed instances. Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.) solved the ALP with multiple runways using matheuristic proposing a relax-and-solve (R&S) algorithm. The authors used the algorithm of Salehipour et al. (2013SALEHIPOUR A, MODARRES M & NAENI LM. 2013. An efficient hybrid meta-heuristic for aircraft landing problem. Computers & Operations Research , 40(1): 207 - 213.) to generate the initial solution. Given this solution, in the “relax” part, the aircraft could change its position in the landing sequence and in the “solve” part, a sequence with a viable solution was obtained. Computational tests were performed on instances of Beasley (1990BEASLEY JE. 1990. OR-Library: Distributing Test Problems by Electronic Mail. The Journal of the Operational Research Society, 41(11): 1069-1072.) with up to 500 aircraft and comparisons were performed with thirteen approximate solution methods from the literature and CPLEX. Only the R&S did not find the best known solution for an instance in a time of 60.00 seconds.

State-of-the-art mathematical formulations solve instances with up to 50 aircraft without much computational effort. However, they have difficulty in providing good solutions when a large number of aircraft, more than 50, is considered. This study provides an alternative to deal with these situations. It proposes a new mathematical formulation, based on Job Shop, for ALPs with single and multiple runways. To validate its performance, experiments were carried out with 49 benchmark instances and the results were compared with four formulations in the literature. The main objective is to better sequence the set of aircraft ready to land, seeking the least possible deviation from the target landing time and respecting the separation times between aircraft.

We proposed a matheuristic that divides the problem into smaller subproblems, requiring less computational effort. An initial solution is obtained using the Relax-and-fix (RF) heuristic, and then optimized by the Gurobi solver, producing the first feasible solution. Next, the Fix-and- Optimize (FO) heuristic is applied to this solution to improve its current values. The approach is simple but provides a high-quality solution to the problem. Afterwards, the parameter calibration mechanism from RF is applied to this solution, looking for the best solution value for the upper bounds of the ALP to be solved.

The article is organized as follows. In Section 2, we provide a description of the problem, related works and introduce the proposed mathematical formulation. In Section 3, we detail the proposed matheuristic, able to solve instances with a large number of aircraft. In Section 4, we present the computational experiments and compare the results obtained with the proposed formulation with those obtained by state-of-the-art alternatives, we discuss the computational results obtained with the matheuristic and we compared the best results with the First-Come First-Served (FCFS) rule. Finally, in Section 5, we report our main conclusions and point out future work directions.

2 AIRCRAFT LANDING PROBLEM

The aircraft landing problem (ALP) considers the number of aircraft A, in which each aircraft aA has a time interval for landing, called a time window. Let E a be the earliest landing time of the aircraft a. This corresponds to the landing time if the aircraft flies at its maximum speed. The target landing time T a of the aircraft a is the expected landing time when a flies at its ideal speed, known as the cruise speed. The latest landing time L a of the aircraft a is the maximum time allowed to keep the aircraft flying before landing. Deviation from the landing of T a incurs extra costs, g a and h a which are penalties (≥ 0) for the situation when the aircraft lands before or after T a , respectively. Figure 1 represents this time window for each aircraft a .

Figure 1
Time window of an aircraft a.

Additionally, separation times between aircraft are considered for a safe landing. S aa ′ is the minimum separation time between aircraft a and a , with aa and s aa ′ is the minimum separation time between aircraft a and a in different runways, with a ≠ a . Assume R to be the set of runways, with each runway rR. Table 1 summarizes the notations used in this paper.

Table 1
The parameters used in the paper.

Flight controllers commonly use the First-Come First-Served (FCFS) procedure to specify the landing order (Ahmadian & Salehipour, 2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.). In this procedure, each aircraft lands according to the order of its arrival times on the runway, respecting minimum separation times between successive landings. According to Balakrishnan & Chandran (2006BALAKRISHNAN H & CHANDRAN B. 2006. Scheduling Aircraft Landings Under Constrained Position Shifting, .) FCFS is not the ideal method to optimally assign an aircraft to landing on a runway. Figure 2 illustrates an example of what a landing sequence of five aircraft would look like, considering that the target landing time of each of them is 140, 150, 135, 180, and 155 time units, respectively. In addition, a unit of time is considered as the minimum time interval between the landing time of each aircraft and the landing time of the following aircraft. The earliest and latest landing times of the landing for each aircraft were considered as ten units of time before and after the target landing time, respectively (Figure 2a.). In this example (Figure 2b.), FCFS and a simple sequencing based on the landing time in increasing order were used. Using FCFS, two aircraft (actual landing time in red) land late and three at the target landing time (actual landing time in blue). Using the increasing order, all aircraft landed at their target time.

Figure 2
Comparison of landing sequencing using FCFS and using an increasing order.

Four formulations found in the literature (F1, F2, F3, and F4) will be used to compare the formulation proposed in this work (F5). These formulations were chosen because they present competitive results within the literature, in addition to F1 being the classic formulation of the problem. Formulations F1, F2, F3, and F5 are directly comparable, as they consider the same objective. They use the same objective function to minimize earliness and tardiness in aircraft landing, and the constraint sets are equivalent. F4 does not have the same set of solutions as the other four formulations, as it considers the Constrained-Position Shifting (CPS) set of constraints, i.e., the number of displacements from the initial landing position is limited to m times. Therefore, two variations of it (F4a and F4b) will be used in this study for appropriate comparisons. The formulations have time windows and separation restrictions. They mainly differ in their set of separation constraints.

Table 2 presents a summary of the main differences of the mathematical formulations, considering the number of indices of the decision variables, a difference in the approach of separation constraints and general comments.

Table 2
Main differences of mathematical formulations.

The decision variables used in F5 are presented. These variables are of the binary and continuous types:

Binary variables:

  • δaar: takes value 1 if the aircraft a lands before the aircraft a on the runway r with a≠ a′ and rR , and 0 otherwise;

  • yar: takes the value 1 if the aircraft a lands on the runway rR, and 0 otherwise.

Continuous Variables:

  • xa: the landing time of the aircraft a;

  • αa: the amount of time the aircraft a lands before the T a ;

  • βa: the amount of time the aircraft a lands after T a .

The mathematical formulation proposed in this work is based on the formulation for Job Shop by Manne (1960MANNE AS. 1960. On the Job-Shop Scheduling Problem. Operations Research, 8(2): 219-223.). The Job Shop generally consists of ordering a set of jobs to be performed by a set of machines. Each machine can process only one job at a time, the jobs have different routes in the machines and all the machines can start and finish the jobs, in the same way as in the ALP when the aircraft land exclusively on a runway and any one can start or finish the job route. Thus, the jobs are the aircraft and the machines the runways available for landing in order to assign and sequence the routes. Works such as Ernst et al. (1999ERNST A, KRISHNAMOORTHY M & STORER R. 1999. Heuristic and exact algorithms for scheduling aircraft landings. Networks, 34: 229-241.), Bianco et al. (2006BIANCO L, DELL’OLMO P & GIORDANI S. 2006. Scheduling models for air traffic control in terminal areas. Journal of Scheduling, 9: 223-253.) and Bencheikh et al. (2009BENCHEIKH G, BOUKACHOUR J, KHOUKHI F & ALAOUI A. 2009. Hybrid method for aircraft landing scheduling based on a Job Shop formulation. International Journal of Computer Science and Network Security, 9: 78 - 88.) also formulated the ALP using Job Shop, differentiating themselves between them and F5, mainly, in the approach of separation time constraints between the aircraft.

F5 is given by the objective function (1), subject to constraints in typical (2) to (13).

Min a A ( g a α a + h a β a ) (1)

Subject to

E a x a L a , a A (2)

r R y a r = 1 , a A (3)

x a ' x a + S a a ' - M ( 1 - y a r ) - M ( 1 - y a ' r ) - M ( 1 - δ a a ' r ) , a A , a ' A , a ' > a , r R (4)

x a x a ' + S a ' a - M ( 1 - y a r ) - M ( 1 - y a ' r ) - M δ a a ' r , a A , a ' A , a ' > a , r R (5)

x a = T a - α a + β a , a A (6)

α a T a - x a , a A (7)

β a x a - T a , a A (8)

y a r { 0 , 1 } , a A , r R (9)

x a 0 , a A (10)

0 α a T a - E a , a A (11)

0 β a L a - T a , a A (12)

δ a a ' r { 0 , 1 } , a A , a ' A , a ' > a , r R (13)

The objective function (1) aims to minimize the cost of deviation from the target times (T a ). The constraints in (2) indicate the time window of each aircraft a to be respected. In (3), the constraints ensure that an aircraft lands on only one runway. The constraints on (4) and (5) are those of separation in relation to the aircraft a and a′. In (4), the aircraft a′ can only land after the aircraft a, considering the separation time between them, if they land on the same runway r (y ar = y ar = 1) and a′ land after a on runway r (δ aar = 1). In constraint (5), it represents the same considering that the aircraft a lands after the aircraft a on the same runway and a does not land after a on the runway r (δ aar = 0). The constraints in (6), (7) and (8) relate the decision variables x a with T a , α a and β a for linearization of the objective function. The domain of variables refers to constraints on (9) to (13).

If |A| denotes the number of aircraft and |R| denotes the number of runways, in the table 3 presents the total number of continuous and binary decision variables and constraints for each formulation. The total number of constraints excludes the limits of the decision variables.

Table 3
Relation of the total number of decision variables and constraints of each formulation.

3 APPLICATION OF MATHEURISTIC

This section presents the proposed matheuristic to solve the ALP with instances with a larger number of aircraft. Matheuristics are effective methods for dealing with generalized assignment problems and NP-hard problems (Maniezzo et al., 2021MANIEZZO V, BOSCHETTI MA & STU¨ TZLE T. 2021. Matheuristics. Springer International Publishing.). With the increase in the number of aircraft available for landing, it will have a high computational cost. Thus, applying matheuristics generates a better solution with less computational time.

The proposed matheuristic works as follows. First, the initial solution is obtained using the relax- and-fix (RF) heuristic and the first feasible solution obtained by Gurobi (G). Afterwards, an improvement procedure based on the decomposition of binary variables, fix-and-optimize (FO) is applied (RF + FO or G + FO), and Gurobi solved minor subproblems.

3.1 Generating initial solution using Relax-and-Fix

According to Farhadi (2016FARHADI F. 2016. Heuristics and Meta-heuristics for Runway Scheduling Problems. pp. 141- 163.), exploring the linear relaxation of NP-hards problems has a low computational cost to obtain viable solutions for these problems. It is a common strategy to find upper bounds, relaxing the binary domain or integer variables.

The RF heuristic is used to reduce the complexity of the problem when significantly increasing the number of aircraft and generating an initial solution for the ALP. In summary, RF optimizes the problem with different fixed, binary, and relaxed variables in each iteration until a final solution is found.

In this work, when RF is used, the sets of binary variables Q of a MIP are partitioned into D distinct sets Q 1, . . . , Q D . Thus, the variables related to Q are relaxed. The number of iterations is defined by the value of D. In the d = 1 iteration, the subproblem MIP1 has all variables related to Q 1 kept in their original domain and Gurobi solves this subproblem. In the next iterations, d = 2, . . . , D, the subproblem MIPd related to variables Q d are kept in their original domain and variables related to Q d−1 are fixed together with the solution obtained from MIPd−1 if the solution is feasible; if not, the algorithm stops. In the end, in the D iteration, the solution of MIPD is the heuristic solution to the problem; if the solution is infeasible, the heuristic failed. The timeout for solving each MIPd subproblem is equal to Timeout \ D. If MIPd does not use up all of its allocated time, it is used in the next iteration.

To partition the set Q, we consider the overlap (O) Od = Qd Qd1 0, with d = 1, . . . , D − 1. In the d iteration, with 1 < d < D, the variables related to Q d−1 \ O d−1 are fixed. The variables related to Q d remain in the original domain and the rest are relaxed. The window (W) indicates the number of elements in Q d . The Q D may have a smaller or equal number of elements because it is the last set. Higher values of W lead to larger subproblems, but a smaller amount of them leads to higher computational effort. The values of O d are related to the number of subproblems, larger numbers will lead to more subproblems, increasing computational time. Algorithm 1 presents the pseudocode for initial solution generation.

Algorithm 1
Pseudocode for initial solution generation.

As an example of how RF partitions work, consider the instance with one runway and 100 aircraft (A = 1, 2, 3, . . . , 100). Choosing W = 20 and O = 5, we have D = 7 iterations (1, . . . , 7) and consequently, seven partitions in Q = Q 1, . . . , Q 7 and six overlap O = O 1, . . . , O 6 (in bold):

  • Q1 = [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20];

  • Q2 = [16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35];

  • Q3 = [31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50];

  • Q4 = [46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65];

  • Q5 = [61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80];

  • Q6 = [76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95];

  • Q7 = [91 92 93 94 95 96 97 98 99 100].

The domains of the variables would look like this, considering the example above, assuming that all iterations get a feasible solution: d = 1: variables in Q 1 kept in the original domain and the other [21, 22,..., 100] are relaxed; d = 2: variables in Q 2 kept in the original domain, Q 1 \ O 1 = [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15] are set with the solutions found and the other [36, 37,..., 100] are relaxed; and so on until D = 7.

3.2 Improving the solution with Fix-and-Optimize

The FO heuristic was introduced by Pochet & Wolsey (2006POCHET Y & WOLSEY L. 2006. Production Planning by Mixed Integer Programming, .) under the name exchange and was renamed by Helber & Sahling (2010HELBER S & SAHLING F. 2010. A fix-and-optimize approach for the multi-level capacitated lot sizing problem. International Journal of Production Economics, 123(2): 247-256.) to solve the dynamic multi-level capacitated lot sizing problem (MLCLSP). The FO algorithm improves a sequence by relaxing a subsequence of the current sequence, which has a subset of consecutive aircraft, and reconstructs a viable sequence using (Ahmadian & Salehipour, 2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.) optimization techniques From an initial solution, a subset of variables is fixed and the resulting subproblems are optimized. As in RF, the set of binary variables Q is partitioned into D subsets Q 1, . . . , Q D .

At each iteration d, 1 ≤ dD, a subproblem is defined, in which the variables Q d remain in their original domain and the others are fixed to the incumbent solution, i.e., the best feasible solution found until the moment. The subproblem obtained is solved and a new best solution is found, this becomes the new incumbent solution. The algorithm stops when a best solution is not found after a stipulated maximum number of iterations or after reaching the execution time limit. The variables were partitioned by the number of aircraft, as this is the largest set of ALP parameters for the analyzed cases. In Algorithm 2, we present the pseudocode for FO.

Algorithm 2
Pseudocode for the FO heuristic.

4 COMPUTATIONAL EXPERIMENTS

The mathematical formulations and matheuristics were implemented in Python 3 using the optimization package Gurobi Optimizer version 9.1.2. The computational tests were performed on a computer with a 2.70 GHz Intel Core i7-7500U processor, 16 GB of RAM and Linux Ubuntu 16.04 LTS operating system.

4.1 Comparison of Mathematical Formulations

This section presents the computational results obtained with the proposed formulation based on the Job Shop (F5). Furthermore, it was compared with the formulations of Beasley et al. (2000BEASLEY JE, KRISHNAMOORTHY M, SHARAIHA YM & ABRAMSON D. 2000. Scheduling Aircraft Landings-The Static Case. Transportation Science, 34(2): 180-197.) (F1), Salehipour et al. (2013SALEHIPOUR A, MODARRES M & NAENI LM. 2013. An efficient hybrid meta-heuristic for aircraft landing problem. Computers & Operations Research , 40(1): 207 - 213.) (F2), Faye (2015FAYE A. 2015. Solving the Aircraft Landing Problem with time discretization approach. European Journal of Operational Research, 242(3): 1028 - 1038.) (F3) and the two F4 variations (F4a and F4b) from Ikli et al. (2020IKLI S, MANCEL C, MONGEAU M, OLIVE X & RACHELSON E. 2020. Coupling Mathematical Optimization and Machine Learning for the Aircraft Landing Problem. In: ICRAT 2020, 9th International Conference for Research in Air Transportation. ICRAT 2020. Tampa, United States.), to assess which one has the best computational performance. As the stopping criteria, 3,600 seconds were assigned to each instance as a total timeout for the execution of the B&C.

The tests considered 13 instances3 3 http://people.brunel.ac.uk/∼mastjjb/jeb/orlib/airlandinfo.html of Beasley (1990BEASLEY JE. 1990. OR-Library: Distributing Test Problems by Electronic Mail. The Journal of the Operational Research Society, 41(11): 1069-1072.) with the number of aircraft (A ) ranging from 10 to 500. Furthermore, each instance has the earliest landing time (E a ), the target landing time (T a ), the latest landing time (L a ) and the cost penalties (g a and h a ) for each aircraft a and the minimum separation time (S aa ′ ) between the aircraft a and a .

The separation time s aa ′ between the aircraft on different runways was considered as zero, as found in the literature. The number of runways R was assigned until the value of the objective function is equal to zero, i.e., all aircraft land at their target time. The value for M is the smallest possible value equal to L a + S aa ′ − E a ′ . F4a has the relaxation of CPS constraints, using m = A − 1. Thus, each aircraft can land in all possible positions. In F4b, the CPS constraints are removed. This allows a fairer comparison of these formulations (F4a and F4b) with the others, which do not use CPS constraints and consider separation constraints in a different way.

Table 4 shows the computational results obtained by F5. It shows that 40 optimal solutions (GAP equal to 0.00%) are obtained from the 49 analyzed instances, reaching the maximum execution time in 9 instances (time equal to 3,600.00 seconds). The computational time is small for instances with up to 50 aircraft when compared to other instances with a larger number of aircraft. The GAP found in these instances ranges from 19.17% to 73.32%. F5, on the other hand, presented an average GAP of 8.02% and average execution time equal to 693.95 seconds.

Table 4
Computational results of the Job Shop based formulation (F5).

Table 5 shows the percentage of the difference in the value of the solution found by the compared formulations. For each formulation, considering f as the solution value and Best Solution as the lowest solution value among all the formulations, we represented the quality of its solution by its average percentage, which is calculated as follows:

Table 5
Percentage difference between the comparable mathematical formulations.

Difference ( % ) = f ( i ) - f ( BestSolution ) f ( i ) × 100 , i = F 1 , F 2 , F 3 , F 4 a , F 4 b and F 5 (14)

Given by the difference between the solution value the average percentage of F3 (0.04%) and F5 (0.08%) were the smallest of all formulations. Additionally, these F3 and F5 formulations found the best solution in 46 and 45 of the 49 analyzed instances.

According to the obtained results, the six formulations (F1, F2, F3, F4a, F4b and F5) found the optimal solution (GAP equal to 0%) for the first eight instances (Airland1 to Airland8), which total 25 instances with the variation of runways R. They correspond to a maximum of 50 aircraft, making formulations with fewer variables and constraints. The other instances (Airland9 to Airland13) have a larger number of aircraft (100 to 500), making their solution by exact methods difficult. It is important to emphasize that the computational cost increases with the growth in the number of aircraft. However, there was a considerable decrease in computational effort when two or more runways were used in the instances.

Table 6 summarizes the test results obtained for a time budget of 3,600 seconds for all evaluated formulations. OptSol and ViSol are the amount of optimal solution and the amount of viable solution found by the formulations, respectively.

Table 6
General comparison between mathematical formulations.

This table shows that the formulations F1 and F3 found the largest number of optimal solutions (41 instances), i.e., GAP equal to 0.00%. F5 finds just one instance less (40 instances) in the number of optimal solutions. The average percentage GAP, the variation of the GAPs and the average time in seconds of execution of the F3 were the smaller compared to the others.

In view of the remarks made earlier and the results presented in Tables 5 and 6, F1, F3 and F5 showed superiority to other comparable formulations in this study. In order to perform an improved performance analysis, authors such as Dolan & Moré (2002DOLAN E & MORÉ J. 2002. Benchmarking Optimization Software with Performance Profiles. Mathematical Programming, 91(2): 201 - 213.) proposed a tool known as Performance Profiles to facilitate the visualization and interpretation of results obtained in computational experiments by different approaches. Figure 3 presents the performance of formulations F1, F3 and F5 considering the value of the objective function. For the proper analysis, P is the set of 49 instances, S are the 3 formulations analyzed and t is the evaluation metric. Thus, the quality of the r ps solution is calculated in (15),

Figure 3
Performance of mathematical formulations considering the value of the objective function.

r p s = t p s m i n { t p s : s S } (15)

i.e., the solution obtained by each formulation divided by the minimum value of the solution among the three formulations. In this equation, t ps is the value of the solution by the formulation in S to solve the instance p.

Where Np is the total number of instances, the probability P s (τ) that each formulation is equal to or better than τ is given by:

P s ( τ ) = | p P : r p s τ | N p (16)

The τ = 1.0000 represents the percentage of instances that the formulation in S gives the best/smallest solution. Thus, we observe that, in relation to the objective function, the F3 and F5 formulations obtained approximately 94% of these solutions in τ = 1.0000. This means that these formulations presented a better solution in 46 of the 49 instances compared to the other formulations. Formulation of F1 have approximately 88% of the solutions. Furthermore, F3 was able to converge faster to the probability (P s (τ)) of 100%. Then, we had F5 converging on τ = 1.0160 and then F1 only in the last τ.

Finally, we can conclude that the F3 and F5 formulations are competitive and have shown good computational results. They stood out as they presented the highest number of the smallest solutions found (94%) compared to the F1 formulation as the aircraft landing problem aim to minimize the earliness and tardiness of landings. F3 stood out as it had the highest number of optimal solutions, in addition to lower average time and average percentage GAP compared to the other five formulations. In addition F3 and F5 presented only three and four instances, respectively, with values higher than the best value of the solution compared to the other formulations analyzed.

4.2 Computational results from matheuristics

This section presents the results obtained by the matheuristics proposed in this work. The F5 was used together with the RF and FO to solve the ALP in instances with a number between 100 and 500 aircraft. Computational tests were performed in 9 of the 49 instances in which the F5 formulation did not find the optimal solution in a runtime of 3,600.00 seconds. These instances have a greater number of aircraft, making it difficult to resolve them by exact methods.

In order to define good values of W and O in the RF, we perform the calibration of these through extensive tests to obtain a quality solution in a viable runtime. In Appendix 5 APPENDIX A. RESULTS OBTAINED BY VARYING THE VALUES OF W AND O In this section we present the detailed computational results obtained for the calibration of the parameters W and O for the ALP. The RF does not always result in the global optimum for the complete mathematical formulation of the ALP, due to the value of the solution being given through the subproblems obtained with the linear relaxations. These solutions are considered upper bounds for the problem, raising the importance of defining good parameters for the subproblems. The subproblem sizes (W ) and overlap (O) are the essential parameter decisions to obtain good bounds. First, the value for W was considered less than or equal to 50 because it is the number of aircraft of the instances that were solved with small computational times and returning the optimal value of the solution. Thus, consider W containing 20, 25, 30, 35, 40, 45 and 50 aircraft and varying the O by 5, 10 and 15. In the Tables A1, A2, A3, A4, A5, A6 and A7 we present the results obtained in a time limit of 600.00 seconds. Table A1 Parameter calibration: W = 20. Instances A R W 20 0 5 W 20 0 10 W 20 0 15 Solution Value Time (s) Solution Value Time (s) Solution Value Time (s) Airland9 100 1 5,611.70 64.32 5,611.70 56.68 5,611.70 44.93 Airland10 150 1 12,658.96 472.41 12,477.24 442.77 12,765.15 444.17 2 1,143.70 3.00 1,143.70 3.86 1,178.21 7.71 Airland11 200 1 12,440.12 68.93 12,418.32 55.13 12,418.32 58.10 Airland12 250 1 16,596.57 565.48 16,249.68 559.33 16,161.42 565.47 2 1,696.59 11.46 1,695.62 16.21 1,695.62 45.41 Airland13 500 1 37,743.68 586.97 37,405.35 587.08 37,969.88 595.35 2 4,023.05 86.57 3,942.80 126.78 3,942.80 258.27 3 719.89 118.27 690.43 186.07 749.09 370.50 Average: 219.71 225.99 265.55 Table A2 Parameter calibration: W = 25. Instances A R W 25 0 5 W 25 0 10 W 25 0 15 Solution Value Time (s) Solution Value Time (s) Solution Value Time (s) Airland9 100 1 5,611.70 67.15 5,611.70 63.13 5,611.70 66.39 Airland10 150 1 12,726.00 450.53 12,450.94 442.48 12,385.59 443.62 2 1,155.90 2.89 1,143.70 3.34 1,143.70 4.84 Airland11 200 1 12,418.32 183.31 12,418.32 109.86 12,418.32 159.31 Airland12 250 1 16,527.09 554.46 16,222.25 565.28 16,161.42 556.97 2 1,708.81 9.79 1,708.81 19.16 1,695.62 26.96 Airland13 500 1 37,714.25 560.89 37,574.43 541.15 37,674.88 593.11 2 3,920.48 75.36 3,920.39 93.35 3,920.39 169.12 3 698.62 83.00 767.69 120.75 673.85 173.47 Average: 220.82 217.61 243.76 Table A3 Parameter calibration: W = 30. Instances A R W 30 0 5 W 30 0 10 W 30 0 15 Solution Value Time (s) Solution Value Time (s) Solution Value Time (s) Airland9 100 1 5,625.92 210.90 5,635.23 183.45 5,635.23 95.39 Airland10 150 1 12,563.92 500.41 12,763.55 450.53 12,646.96 480.58 2 1,143.70 3.02 1,144.98 6.64 1,143.70 4.87 Airland11 200 1 12,418.32 327.88 12,418.32 436.78 12,418.32 415.05 Airland12 250 1 16,241.29 600.11 16,247.44 554.55 16,132.58 565.38 2 1,695.62 20.31 1,695.62 10.28 1,695.62 29.37 Airland13 500 1 37,754.64 588.32 37,700.84 601.72 37,902.51 584.98 2 3,921.60 82.65 3,920.39 101.74 3,920.39 559.06 3 673.85 67.43 673.85 97.75 673.85 125.38 Average: 266.78 271.49 317.78 Table A4 Parameter calibration: W = 35. Instances A R W 35 0 5 W 35 0 10 W 35 0 15 Solution Value Time(s) Solution Value Time(s) Solution Value Time(s) Airland9 100 1 5,611.99 153.23 5,625.92 213.93 5,625.92 250.24 Airland10 150 1 12,466.87 480.52 12,429.73 500.31 12,616.38 450.44 2 1,182.22 4.47 1,143.70 4.31 1,143.70 13.11 Airland11 200 1 12,440.12 436.26 12,440.12 490.06 12,440.12 460.79 Airland12 250 1 16,491.69 534.05 16,241.29 600.12 16,358.23 554.47 2 1,695.62 14.86 1,695.62 24.73 1,695.62 30.63 Airland13 500 1 37,625.87 591.79 37,827.55 599.14 37,735.22 600.21 2 3,920.39 122.09 3,920.39 136.13 3,920.39 580.33 3 673.85 75.33 673.85 71.07 673.85 104.50 Average: 268.07 293.31 338.30 Table A5 Parameter calibration: W = 40. Instances A R W 40 0 5 W 40 0 10 W 40 0 15 Solution Value Time (s) Solution Value Time (s) Solution Value Time (s) Airland9 100 1 5,611.99 242.94 5,635.23 202.61 5,635.23 199.66 Airland10 150 1 12,619.01 480.30 12,511.03 480.41 12,638.95 500.36 2 1,149.44 9.06 1,143.70 24.57 1,143.70 54.61 Airland11 200 1 12,488.84 451.96 12,494.10 440.71 12,440.12 349.83 Airland12 250 1 16,221.72 525.60 16,253.74 534.00 16,203.33 600.21 2 1,695.62 21.57 1,695.62 20.71 1,695.62 73.17 Airland13 500 1 38,161.09 567.49 37,858.25 585.33 37,966.50 601.78 2 3,960.67 261.08 3,920.39 132.72 3,920.39 206.82 3 706.33 53.40 706.33 57.46 673.85 71.10 Average: 290.38 275.39 295.28 Table A6 Parameter calibration: W = 45. Instances A R W 45 0 5 W 45 0 10 W 45 0 15 Solution Value Time(s) Solution Value Time(s) Solution Value Time(s) Airland9 100 1 5,611.70 203.91 5,611.70 219.11 5,635.23 210.30 Airland10 150 1 12,686.96 450.27 12,471.05 480.16 12,519.15 480.27 2 1,143.70 31.31 1,143.70 14.56 1,143.70 32.82 Airland11 200 1 12,471.81 539.82 12,440.12 523.25 12,440.12 445.13 Airland12 250 1 16,550.16 514.98 16,148.82 525.43 16,345.00 533.95 2 1,695.62 38.26 1,695.62 41.36 1,695.62 33.78 Airland13 500 1 38,165.77 561.69 38,099.89 564.96 37,932.25 594.11 2 3,920.39 428.19 3,920.39 513.49 3,920.39 315.83 3 676.40 55.48 676.40 53.65 673.85 175.83 Average: 313.77 326.22 313.56 Table A7 Parameter calibration: W = 50. Instances A R W 50 0 5 W 50 0 10 W 50 0 15 Solution Value Time (s) Solution Value Time (s) Solution Value Time (s) Airland9 100 1 5,635.23 200.63 5,635.23 208.48 5,635.23 233.34 Airland10 150 1 12,581.48 450.27 12,356.13 450.28 12,407.26 480.16 2 1,143.70 127.56 1,143.70 27.84 1,143.70 39.72 Airland11 200 1 12,418.32 486.67 12,418.32 600.10 12,418.32 519.17 Airland12 250 1 16,531.52 600.20 16,434.86 514.98 16,261.72 525.43 2 1,709.87 52.65 1,695.62 75.69 1,695.62 84.49 Airland13 500 1 38,666.29 552.73 38,065.71 570.10 38,465.66 563.80 2 3,989.43 270.97 3,920.39 455.30 3,920.39 563.00 3 673.85 42.15 673.85 46.38 673.85 53.08 Average: 309.32 327.68 340.24 For the proper analysis, in Table A8 we compare the values of the solution obtained for each variation of W and O through the difference being the Solution Value of the variation by the Best Solution among all the variations divided by the Solution Value of the variance multiplied by 100. For the proper analysis, in Table A8, we compare the values of the solution obtained for each variation of W and O, considering f as the solution value, which is calculated as follows: Table A8 Percentage difference of the values of the solutions obtained between the variations of W and O. Instances A R Best solution W 20 O 5 W 20 O 10 W 20 O 15 W 25 O 5 W 25 O 10 W 25 O 15 Airland9 100 1 5,611.70 0.00 0.00 0.00 0.00 0.00 0.00 Airland10 150 1 12,356.13 2.39 0.97 3.20 2.91 0.76 0.24 2 1,143.70 0.00 0.00 2.93 1.06 0.00 0.00 Airland11 200 1 12,418.32 0.18 0.00 0.00 0.00 0.00 0.00 Airland12 250 1 16,132.58 2.80 0.72 0.18 2.39 0.55 0.18 2 1,695.62 0.06 0.00 0.00 0.77 0.77 0.00 Airland13 500 1 37,405.35 0.90 0.00 1.49 0.82 0.45 0.72 2 3,920.39 2.55 0.57 0.57 0.00 0.00 0.00 3 673.85 6.40 2.40 10.04 3.55 12.22 0.00 Average: 1.70 0.52 2.05 1.28 1.64 0.13 Instances A R Best solution W 30 O 5 W 30 O 10 W 30 O 15 W 35 O 5 W 35 O 10 W 35 O 15 Airland9 100 1 5,611.70 0.25 0.42 0.42 0.01 0.25 0.25 Airland10 150 1 12,356.13 1.65 3.19 2.30 0.89 0.59 2.06 2 1,143.70 0.00 0.11 0.00 3.26 0.00 0.00 Airland11 200 1 12,418.32 0.00 0.00 0.00 0.18 0.18 0.18 Airland12 250 1 16,132.58 0.67 0.71 0.00 2.18 0.67 1.38 2 1,695.62 0.00 0.00 0.00 0.00 0.00 0.00 Airland13 500 1 37,405.35 0.93 0.78 1.31 0.59 1.12 0.87 2 3,920.39 0.03 0.00 0.00 0.00 0.00 0.00 3 673.85 0.00 0.00 0.00 0.00 0.00 0.00 Average: 0.39 0.58 0.45 0.79 0.31 0.53 Instances A R Best solution W 40 O 5 W 40 O 10 W 40 O 15 W 45 O 5 W 45 O 10 W 45 O 15 Airland9 100 1 5,611.70 0.01 0.42 0.42 0.00 0.00 0.42 Airland10 150 1 12,356.13 2.08 1.24 2.24 2.61 0.92 1.30 2 1,143.70 0.50 0.00 0.00 0.00 0.00 0.00 Airland11 200 1 12,418.32 0.56 0.61 0.18 0.43 0.18 0.18 Airland12 250 1 16,132.58 0.55 0.75 0.44 2.52 0.10 1.30 2 1,695.62 0.00 0.00 0.00 0.00 0.00 0.00 Airland13 500 1 37,405.35 1.98 1.20 1.48 1.99 1.82 1.39 2 3,920.39 1.02 0.00 0.00 0.00 0.00 0.00 3 673.85 4.60 4.60 0.00 0.38 0.38 0.00 Average: 1.26 0.98 0.53 0.88 0.38 0.51 Instances A R Best solution W 50 O 5 W 50 O 10 W 50 O 15 Airland9 100 1 5,611.70 0.42 0.42 0.42 Airland10 150 1 12,356.13 1.79 0.00 0.41 2 1,143.70 0.00 0.00 0.00 Airland11 200 1 12,418.32 0.00 0.00 0.00 Airland12 250 1 16,132.58 2.41 1.84 0.79 2 1,695.62 0.83 0.00 0.00 Airland13 500 1 37,405.35 3.26 1.73 2.76 2 3,920.39 1.73 0.00 0.00 3 673.85 0.00 0.00 0.00 Average: 1.16 0.44 0.49 Difference(%) = f ( W i O j ) - f ( Best Solution ) f ( i ) × 100 , i = 20, 25, 30, 35, 40, 45 and 50, j = 5, 10 and 15 (20) We observed that the best variance was with W = 25 and O = 15 among the values tested with an average execution time of 243.76 seconds. In this variation of parameter values, all solutions found are optimal (no instance with time ≥ 600.00 seconds) and only 3 solutions are greater than the best solution among all variations. In Figure A1, the values of these average differences are better observed. Figure A1 Average percentage difference of the solutions for the variations in the values of W and O. In Figure A2, we observe the behavior of the variations of O in relation to the average percentage difference with the values of W . In general, higher values for O lead to smaller differences, i.e., better upper bounds. Furthermore, in Figure A3, the average time spent, in seconds, of the values of O according to the variations of W , was analyzed. Noting that higher values for O generate a higher computational cost, on average. Figure A2 Behavior of the O parameter in relation to the average percentage difference. Figure A3 Behavior of the O parameter in relation to the average time percentage. Finally, we add that larger values for W generate fewer but larger subproblems, which become more difficult to optimally solve in a short time. Moreover, larger values of O increase the number of subproblems, requiring more computational effort. , details of this calibration are presented. We use W worth 25 and O worth 15 to generate the initial solution.

In Table 7, we find the computational results obtained through the resolution of F5 with Gurobi in 3,600.00 seconds, of RF in 600.00 seconds, the combination of RF with FO in 1,200.00 seconds, and the combination of the first feasible solution obtained by Gurobi with FO (G + FO) in 1,200.00 seconds. It is worth noting that the runtime spent combining RF with FO was the time taken to generate the initial solution using RF plus the time taken by the FO to improve the solution. The maximum number of iterations for the FO used was 1,000. The percentage difference of the solution, considering f as the value of the solution and the Best Solution as the lowest value between the methods, which is calculated as follows:

Table 7
Computational results obtained by matheuristics.

Difference(%) = f ( i ) - f ( Best Solution ) f ( i ) × 100 , i = F5, RF, RF + FO and G + FO (17)

The RF heuristic finds the best solution value in 7 out of 9 instances, finding the other two smallest solutions combined with FO. The percentage difference is, on average, 0.11% for RF and 0.00% for RF + FO highlighting the importance of RF to obtain a good upper bound and the use of adequate values for W and O. The use of the first feasible solution of Gurobi combined with FO (G + FO) was the heuristic that returned the worst solution values, with average difference of 31.34% of the best solution found.

Furthermore, the computational effort comparing the heuristic methods with the exact method was significantly lower. The best method, RF + FO, reduces approximately 90% of the average execution time compared to the computational times spent to solve using F5 and presenting better or equal solution values.

Table 8 shows the percentage difference of solution values compared to Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.). The authors found the best or equal solutions compared to the literature with execution times around 60.00 seconds. Considering f as the solution value, the percentage difference is calculated as follows:

Table 8
Percentage difference of objective solution value of solving methods with Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.).

Difference(%) = f ( i ) - f ( Ahmadian & Salehipour (2022) ) f ( i ) × 100 , i = F5, RF, RF + FO and G + FO (18)

The heuristic that combines RF with FO only found the highest solution value for the Airland13 instance with a runway for reaching the stopping criterion of 1,200.00 execution seconds. This instance was the most computationally complex as it had 500 aircraft to be sequenced on just one runway, Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.) also does not find the optimal solution for this instance showing an average difference of 0.18% on average for all 9 instances. Using only the RF, this average percentage difference is also relatively small, which is 0.29%.

We used the Performance Profiles of Dolan & Moré (2002DOLAN E & MORÉ J. 2002. Benchmarking Optimization Software with Performance Profiles. Mathematical Programming, 91(2): 201 - 213.) for the performance analysis of matheuristics with Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.). Considering P = 9 and S = 4, G + FO was not considered because it presented the highest average difference. In Figure 4, we present the performance of the methods.

Figure 4
Performance of matheuristics considering the value of the objective function.

We observed that, in relation to the objective function, RF + FO obtained approximately 89% of these solutions in τ = 1.00, i.e., it presented a better solution in 8 of the 9 instances compared to those of Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.). In RF and F5, they presented approximately 67% (6 instances) and 44% (4 instances) of the solutions, respectively. It is worth mentioning that F5 was solved by the exact method.

4.3 Comparison with FCFS rule

Air traffic controllers sequence the aircraft according to the FCFS rule, in Table 9, we made a comparison with the best values of the solution found between the exact method and with matheuristics with the value of the solution using FCFS to verify how optimization methods can improve costs at airports. The percentage improvement is calculated:

Table 9
Percentage improvement between the best solutions found in this work and the FCFS rule.

Improvement(%) = f ( FCFS ) - f ( Best Solution ) f ( FCFS ) × 100 , (19)

We observe that FCFS finds a solution equal to the optimization methods in four instances of the 49 (Airland6 with one and three lanes, Airland7 with one and two lanes). In the other instances, an improvement of up to 100.00% in the values of the solutions (11 instances) is achieved. This improvement ranged from 17.32% to 100.00% with an average percentage of 72.28%. The numbers obtained show even more the efficiency of the proposed solution methods.

5 CONCLUSIONS

The Aircraft Landing Problem (ALP) was addressed in this work, which aims to minimize earliness and tardiness of aircraft ready for landing by assigning which runway to land and sequencing them. A new mathematical formulation was proposed to solve the ALP, based on Job Shop (F5). Comparisons were performed with four formulations found in the literature, one of which was adapted, generating two formulations for direct comparison with the proposal. F5 proved to be competitive with the others for finding a greater number of smaller solutions to the problem, 94% of the 49 instances of the literature analyzed.

Furthermore, for the instances that F5 did not find an optimal solution in the stipulated execution time (9 instances), we proposed an efficient matheuristic to solve the ALP with the largest number of aircraft available for landing. First finding an initial solution using relax-and-fix (RF) and the first solution obtained by Gurobi. Then, applying fix-and-optimize (FO) to improve the solution. RF extensive tests were performed to calibrate the parameters W and O. Computational experiments were performed and compared with the results obtained with the mathematical formulation (F5), with the RF, with the combination of the initial solution found by the RF and then the FO (RF + FO) and with the combination of the first solution found with the Gurobi and then FO (G + FO). RF finds the best solution in 7 out of 9 instances and combined with FO, finds the other two best solutions in an average runtime of 335.47 seconds. Comparisons were also performed with the results of Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.), not finding the best solution in just one instance (Airland13 with one runway) with RF + FO. Ahmadian & Salehipour (2022AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.) also does not find the optimal solution for this instance and the execution times were around 60.00 seconds. In addition, we performed comparisons with the best solution values found by the exact method and matheuristics with the FCFS rule. The results obtained evidence of the efficiency of the solution methods, improving an average of 72.28% compared to FCFS.

Future work may be carried out for dynamic ALP, when some landing times are unknown at the beginning of the optimization, and new objective functions may also be considered and analyzed according to the formulations. In addition, the study and implementation of other matheuristics that proved to be efficient, presenting a high quality solution for instances with a larger number of aircraft.

Acknowledgements

The authors would like to thank the financial support of the São Paulo Research Foundation (FAPESP) grant 2017/21488-2 and grant 2013/07375-0.

References

  • AHMADIAN MM & SALEHIPOUR A. 2022. Heuristics for flights arrival scheduling at airports. International Transactions in Operational Research, 29(4): 2316-2345.
  • BALAKRISHNAN H & CHANDRAN B. 2006. Scheduling Aircraft Landings Under Constrained Position Shifting, .
  • BEASLEY JE. 1990. OR-Library: Distributing Test Problems by Electronic Mail. The Journal of the Operational Research Society, 41(11): 1069-1072.
  • BEASLEY JE, KRISHNAMOORTHY M, SHARAIHA YM & ABRAMSON D. 2000. Scheduling Aircraft Landings-The Static Case. Transportation Science, 34(2): 180-197.
  • BEASLEY JE, SONANDER J & HAVELOCK P. 2001. Scheduling aircraft landings at London Heathrow using a population heuristic. Journal of the Operational Research Society, 52(5): 483- 493.
  • BENCHEIKH G, BOUKACHOUR J, KHOUKHI F & ALAOUI A. 2009. Hybrid method for aircraft landing scheduling based on a Job Shop formulation. International Journal of Computer Science and Network Security, 9: 78 - 88.
  • BIANCO L, DELL’OLMO P & GIORDANI S. 2006. Scheduling models for air traffic control in terminal areas. Journal of Scheduling, 9: 223-253.
  • BRISKORN D & STOLLETZ R. 2014. Aircraft landing problems with aircraft classes. Journal of Scheduling, 17(1): 31-45.
  • CHENG VHL, CRAWFORD LS & MENON PK. 1999. Air traffic control using genetic search techniques. In: Proceedings of the 1999 IEEE International Conference on Control Applications (Cat. No.99CH36328), vol. 1. pp. 249-254.
  • DOLAN E & MORÉ J. 2002. Benchmarking Optimization Software with Performance Profiles. Mathematical Programming, 91(2): 201 - 213.
  • ERNST A, KRISHNAMOORTHY M & STORER R. 1999. Heuristic and exact algorithms for scheduling aircraft landings. Networks, 34: 229-241.
  • FARHADI F. 2016. Heuristics and Meta-heuristics for Runway Scheduling Problems. pp. 141- 163.
  • FAYE A. 2015. Solving the Aircraft Landing Problem with time discretization approach. European Journal of Operational Research, 242(3): 1028 - 1038.
  • GIRISH BS. 2016. An efficient hybrid particle swarm optimization algorithm in a rolling horizon framework for the aircraft landing problem. Applied Soft Computing, 44: 200 - 221.
  • HELBER S & SAHLING F. 2010. A fix-and-optimize approach for the multi-level capacitated lot sizing problem. International Journal of Production Economics, 123(2): 247-256.
  • IKLI S, MANCEL C, MONGEAU M, OLIVE X & RACHELSON E. 2020. Coupling Mathematical Optimization and Machine Learning for the Aircraft Landing Problem. In: ICRAT 2020, 9th International Conference for Research in Air Transportation. ICRAT 2020. Tampa, United States.
  • IKLI S, MANCEL C, MONGEAU M, OLIVE X & RACHELSON E. 2021. The aircraft runway scheduling problem: A survey. Computers & Operations Research, 132: 105336.
  • JI XP, CAO XB & TANG K. 2016. Sequence searching and evaluation: a unified approach for aircraft arrival sequencing and scheduling problems. Memetic Computing, 8: 109-123.
  • LIEDER A, BRISKORN D & STOLLETZ R. 2015. A dynamic programming approach for the aircraft landing problem with aircraft classes. European Journal of Operational Research , 243(1): 61-69.
  • MANIEZZO V, BOSCHETTI MA & STU¨ TZLE T. 2021. Matheuristics. Springer International Publishing.
  • MANNE AS. 1960. On the Job-Shop Scheduling Problem. Operations Research, 8(2): 219-223.
  • PINOL H & BEASLEY J. 2006. Scatter Search and Bionomic Algorithms for the aircraft landing problem. European Journal of Operational Research , 171(2): 439 - 462.
  • POCHET Y & WOLSEY L. 2006. Production Planning by Mixed Integer Programming, .
  • PRAKASH R, PIPLANI R & DESAI J. 2018. An optimal data-splitting algorithm for aircraft scheduling on a single runway to maximize throughput. Transportation Research Part C: Emerging Technologies, 95: 570-581.
  • SALEHIPOUR A & AHMADIAN M. 2017. A heuristic algorithm for the Aircraft Landing Problem. In: 22nd International Congress on Modelling and Simulation. pp. 1344-1349. Hobart - Austrália.
  • SALEHIPOUR A, MODARRES M & NAENI LM. 2013. An efficient hybrid meta-heuristic for aircraft landing problem. Computers & Operations Research , 40(1): 207 - 213.
  • SOOMER M & FRANX G. 2008. Scheduling aircraft landings using airlines preferences. European Journal of Operational Research , 190(1): 277-291.
  • STEVENS G. 1995. An Approach To Scheduling Aircraft Landing Times Using Genetic Algorithms. Tech. rep.. Department of Computer Science. RMIT University.
  • VADLAMANI S & HOSSEINI S. 2014. A novel heuristic approach for solving aircraft landing problem with single runway. Journal of Air Transport Management, 40: 144-148.
  • VERESNIKOV G, EGOROV N, KULIDA E & LEBEDEV V. 2019a. Methods for Solving of the Aircraft Landing Problem. I. Exact Solution Methods. Automation and Remote Control, 80: 1317-1334.
  • VERESNIKOV G, EGOROV N, KULIDA E & LEBEDEV V. 2019b. Methods for Solving of the Aircraft Landing Problem. II. Approximate Solution Methods. Automation and Remote Control , 80: 1502-1518.

APPENDIX A. RESULTS OBTAINED BY VARYING THE VALUES OF W AND O

In this section we present the detailed computational results obtained for the calibration of the parameters W and O for the ALP.

The RF does not always result in the global optimum for the complete mathematical formulation of the ALP, due to the value of the solution being given through the subproblems obtained with the linear relaxations. These solutions are considered upper bounds for the problem, raising the importance of defining good parameters for the subproblems. The subproblem sizes (W ) and overlap (O) are the essential parameter decisions to obtain good bounds.

First, the value for W was considered less than or equal to 50 because it is the number of aircraft of the instances that were solved with small computational times and returning the optimal value of the solution. Thus, consider W containing 20, 25, 30, 35, 40, 45 and 50 aircraft and varying the O by 5, 10 and 15. In the Tables A1, A2, A3, A4, A5, A6 and A7 we present the results obtained in a time limit of 600.00 seconds.

Table A1
Parameter calibration: W = 20.

Table A2
Parameter calibration: W = 25.

Table A3
Parameter calibration: W = 30.

Table A4
Parameter calibration: W = 35.

Table A5
Parameter calibration: W = 40.

Table A6
Parameter calibration: W = 45.

Table A7
Parameter calibration: W = 50.

For the proper analysis, in Table A8 we compare the values of the solution obtained for each variation of W and O through the difference being the Solution Value of the variation by the Best Solution among all the variations divided by the Solution Value of the variance multiplied by 100. For the proper analysis, in Table A8, we compare the values of the solution obtained for each variation of W and O, considering f as the solution value, which is calculated as follows:

Table A8
Percentage difference of the values of the solutions obtained between the variations of W and O.

Difference(%) = f ( W i O j ) - f ( Best Solution ) f ( i ) × 100 , i = 20, 25, 30, 35, 40, 45 and 50, j = 5, 10 and 15 (20)

We observed that the best variance was with W = 25 and O = 15 among the values tested with an average execution time of 243.76 seconds. In this variation of parameter values, all solutions found are optimal (no instance with time ≥ 600.00 seconds) and only 3 solutions are greater than the best solution among all variations. In Figure A1, the values of these average differences are better observed.

Figure A1
Average percentage difference of the solutions for the variations in the values of W and O.

In Figure A2, we observe the behavior of the variations of O in relation to the average percentage difference with the values of W . In general, higher values for O lead to smaller differences, i.e., better upper bounds. Furthermore, in Figure A3, the average time spent, in seconds, of the values of O according to the variations of W , was analyzed. Noting that higher values for O generate a higher computational cost, on average.

Figure A2
Behavior of the O parameter in relation to the average percentage difference.

Figure A3
Behavior of the O parameter in relation to the average time percentage.

Finally, we add that larger values for W generate fewer but larger subproblems, which become more difficult to optimally solve in a short time. Moreover, larger values of O increase the number of subproblems, requiring more computational effort.

Publication Dates

  • Publication in this collection
    21 Aug 2023
  • Date of issue
    2023

History

  • Received
    25 July 2022
  • Accepted
    02 Feb 2023
Sociedade Brasileira de Pesquisa Operacional Rua Mayrink Veiga, 32 - sala 601 - Centro, 20090-050 Rio de Janeiro RJ - Brasil, Tel.: +55 21 2263-0499, Fax: +55 21 2263-0501 - Rio de Janeiro - RJ - Brazil
E-mail: sobrapo@sobrapo.org.br