1. INTRODUCTION
Airspace congestion and delay are growing problems both in the terminal airspace of airports witnessing high traffic density and in many en-route sectors. The growth in demand for air transportation increases the workload of air traffic controllers. Airspace capacity is traditionally defined as the maximum traffic that can be controlled under acceptable workload levels (Majumdar et al., Reference Majumdar, Ochieng and Polak2002). The workload of air traffic controllers is a subjective attribution which is individually associated with each controller's perception. The workload is quantified as a function of task demand, which is driven by airspace complexity (Tobaruela et al., Reference Tobaruela, Schuster, Majumdar, Ochieng, Martinez and Hendrickx2014), and varies depending on objective and perceived complexity factors. While directly measurable factors, such as number of aircraft, their changes in altitude, heading or speed, and potential aircraft conflicts, may be classed as objective complexity factors, the factors related to individual differences, such as individual controllers' abilities, age, fatigue and level of experience, may be classed as perceived complexity factors (Djokic et al., Reference Djokic, Lorenz and Fricke2010). In the current air traffic management (ATM) system, airspace is partitioned into small parts called air traffic control sectors and one air traffic controller is assigned to each sector. The horizontal and vertical sector boundaries should be designed to provide a balanced distribution of workload among controllers, in order to divide the responsibilities properly so that the traffic control workload of each sector is maintained within the limits of each controller's capacity (Cao et al., Reference Cao, Zhu, Tian, Chen, Wu and Du2018). Therefore, as the number of aircraft served in an airspace increases, sectorisation of the airspace becomes one of the most important problems in ATM. Eurocontrol (The European Organisation for the Safety of Air Navigation) defines sectorisation as the means of subdividing the totality of control tasks into manageable portions. The main constraints on ATM capacity are airspace limitations and controller workload. The provision of more sectors to overcome these constraints is a finite strategy. The increase in capacity is not proportional to the number of operational sectors available. Reducing the complexity of airspace is the more efficient solution to increase sector productivity and consequently capacity (Eurocontrol, 2003).
Today, operational research techniques are widely used for handling ATM problems. Although airspace sectorisation has an important place among other ATM problems, the studies conducted within this field are quite limited. Recent studies on the sectorisation problem, motivated by the increase in air traffic demand, focus especially on dynamic solutions. Delahaye et al. (Reference Delahaye, Schoenauer and Alliot1998) considered an air transportation network with flows on it inducing a workload spread over the airspace. The problem was defined as a classical graph partitioning problem and was solved by using genetic algorithms. Trandac et al. (Reference Trandac, Baptiste and Duong2003) considered the optimised airspace sectorisation problem with constraints in which a given airspace was partitioned into a number of sectors. They proposed a constraint programming approach to optimise the sectorisation that satisfied the specific constraints. Yousefi (Reference Yousefi2005) developed a methodology for airspace sectorisation based on air traffic controller workload (ATCW). He partitioned the US National Airspace into three layers with different altitude ranges. Each layer was further tiled to hexagonal cells and ATCW was modelled for each cell using various airspace metrics. Yousefi then developed clustering algorithms using optimisation theory to cluster cells and construct sectors. Klein (Reference Klein2005) presented a potential new partitioning mechanism for the National Airspace System (NAS) that utilised a high-resolution hexagonal grid. He described an algorithm that rapidly processed large amounts of traffic data and created potential airspace centre boundaries starting from a selected number of seed locations. Martinez et al. (Reference Martinez, Chatterji, Sun and Bayen2007) described a method for partitioning airspace into smaller regions based on a peak traffic-counts metric. They developed a traffic dependent flow-graph based algorithm for sectorising the airspace. Basu et al. (Reference Basu, Mitchell and Sabhnani2009) modelled the problem of optimal sectorisation as a geometric partition problem with constraints, and developed a precise computational geometric formulation. Yangzhou and Defu (Reference Yangzhou and Defu2014) proposed a new method for dynamic airspace configuration based on a weighted graph model. They developed a graph partitioning algorithm that divides the weighted graph model into sub-graphs. The method attempts to design the sectors with the objective of balancing workload and minimising coordination workload as well as satisfying geometric constraints. Sergeeva et al. (Reference Sergeeva, Delahaye, Mancel and Vidosavljevic2017) introduced a new genetic algorithm in order to optimise airspace configuration. The developed algorithm generates a sequence of sector configurations for one day of operation with minimised controller workload. Gerdes et al. (Reference Gerdes, Temme and Schultz2018) offered a new solution for dynamic airspace sectorisation. Their approach clusters traffic patterns and uses algorithms for optimisation of airspace, focusing on high capacity utilisation through flexible use of airspace, appropriate distribution of task load for air traffic controllers and rapid adaptation to changed operational constraints.
Recent studies connected with the optimisation of airspace sectorisation are generally performed in dynamic environments by using simulation techniques. In the majority of the accessible literature, the heuristic or metaheuristic algorithms have been used in simulation processes. These kinds of algorithms are able to produce appropriate solutions within a short period of time in the following situations:
- reorganisation of airspace sector structure against unexpected situations such as volcanic eruptions and adverse weather conditions,
- rearrangement of existing sector boundaries according to changing traffic flows (Gerdes et al., Reference Gerdes, Temme and Schultz2018).
A heuristic method is a procedure that determines good or near-optimal solutions to an optimisation problem. Heuristic algorithms carry no guarantee that an optimal solution will be found (Eiselt and Sandblom, Reference Eiselt and Sandblom2000). The global optimum solution of the problem can only be obtained by using a suitable mathematical model. The motivation of this study is that the mathematical modelling approach is seldom chosen in the accessed literature for optimising airspace sectorisation.
In this study, we focus on an approach addressing the airspace sectorisation problem through analytical modelling. In the mathematical model proposed, both ATCW and airspace complexity, which are the most important issues of air traffic systems, are taken into account. In this context, initially the ATCW for en-route airspace levels was formulated, considering all the tasks performed by en-route controllers. This formulation is used to measure the total workload of the related airspace. A square grid based method was used to distribute ATCW evenly across the airspace sectors. Real air traffic data was used for the calculations. In the ATCW formulation, the weighting coefficients and task times of different control activities are obtained from a survey conducted with en-route air traffic controllers. The square grid partitioning of airspace is exercised as an alternative technique. Finally, in order to balance the ATCW across the airspace sectors, a multi-objective mixed integer programming (MOMIP) model has been developed. The difficulty related to scalarisation of the proposed multi-objective model is overcome by utilising the conic scalarisation method (CSM) (Kasimbeyli, Reference Kasimbeyli2010). GAMS software is implemented to solve the scalarised problem. GIS software which enables spatial and temporal data analyses was chosen as a four-dimensional airspace analyser.
2. PROBLEM FORMULATION AND THE MATHEMATICAL MODEL
The shape and size of an airspace sector are both functions of the sector's complexity and control tasks. The complexity factors can reduce sector capacity by increasing the ATCW. The complexity factors affecting the controller's functions can be summarised as: airspace factors (number of sectors, shape and size of the sectors, intersection points of traffic flows etc.), traffic factors (number of aircraft, aircraft mix, aircraft movements, interactions between aircraft etc.) and operational factors (presence of severe weather, amount of coordination required, frequency congestion etc.). These characteristics are specific to each airspace sector and change over time (Wang et al., Reference Wang, Vormer, Hu and Duong2013; Tobaruela et al., Reference Tobaruela, Schuster, Majumdar, Ochieng, Martinez and Hendrickx2014; Zhu et al., Reference Zhu, Cao and Cai2017; Rahman et al., Reference Rahman, Sidik, Nazir and Jamil2018).
An airspace sector should be small enough to accommodate sector functions, while providing a balanced workload. At the same time, it should be big enough to accommodate sector functions while not imposing an excessive workload (Eurocontrol, 2003). In this framework, we develop a mathematical model and a solution approach for the airspace sectorisation problem to distribute ATCW among en-route sectors of airspace without exceeding a defined sector capacity and to determine the optimum sector boundaries and locations at a strategic level.
The airspace sectorisation problem can be modelled as a facility location problem, as previously considered in the study by Yousefi (Reference Yousefi2005). Unlike Yousefi's work, in our study the problem was solved with two objective functions instead of three. In addition, a new objective ATCW function was defined in the model and the problem was solved by a different solution method. The workload equation used to calculate ATCW has also been modified. The problem deals with assigning n square units in which different air traffic services are provided, to m potential locations for sectors, in such a way that total workload is distributed among sectors within a defined sector capacity. In this framework, we partition the airspace in square grid units and calculate the workload in each unit, and then we assign each square unit to one of the potential locations. While these square units are clustered to construct optimum sector boundaries, the distance d ij between each square unit i and each sector location j is minimised. All square units within each cluster must be connected to each other.
2.1. The parameters used in the model
2.1.1. Workload of a square unit
Aircraft monitoring (MON), conflict detection and resolution between aircraft (CDR), coordination between sectors (COR), and vertical/horizontal aircraft manoeuvres (ACM) are the main tasks for en-route controllers. Monitoring is directly related to the number of aircraft in a sector and maximum flight time. A conflict is defined as the loss of the minimum required horizontal and vertical separation distance between two aircraft. Coordination includes transfer of control of a flight between air traffic control units (civil or military) and sectors. Changes in flight level, speed, heading and direct route are defined as the types of aircraft manoeuvre (ICAO, 2007). The total workload of the en-route controller (WL) can be determined as the sum of these four tasks, as given in Equation (1).
$$\begin{align}{\rm WL} &= w_{{\rm MON}}(N_{{\rm MON}}\cdot T_{{\rm MON}}) + w_{{\rm CDR}}\left( {\sum\limits_{k = 1}^5 {(N_{{\rm CD}{\rm R}_k}} \cdot C_{{\rm CD}{\rm R}_k}\cdot T_{{\rm CD}{\rm R}_k})} \right) \cr & + w_{{\rm COR}}\left( {\sum\limits_{l = 1}^4 {(N_{{\rm CO}{\rm R}_l}\cdot C_{{\rm CO}{\rm R}_l}\cdot T_{{\rm CO}{\rm R}_l})} } \right) \cr & + w_{{\rm ACM}}\left( {\sum\limits_{h = 1}^4 {(N_{{\rm AC}{\rm M}_h}\cdot C_{{\rm AC}{\rm M}_h}\cdot T_{{\rm AC}{\rm M}_h})} } \right)\end{align} $$where w, N, T and C signify the weighting coefficient of each task, number of aircraft, task time and weighting factor of subtask, respectively; k, l and h also represent the set of indices of subtask. The details of the ATCW measurement equation are defined in the study of Oktal and Yaman (Reference Oktal and Yaman2011).
The calculation of the workload in each square unit necessitates an extensive data analysis. The stay time of an aircraft in a sector changes according to its performance characteristics and the length of the airway. Additionally, different aircraft manoeuvres give rise to differences in the air traffic control service provided in each sector. In this framework, the workload of each square unit in a time period has been calculated by taking into consideration each flight profile occurring in related airspace. Finally, the total airspace workload is calculated by summing the workloads of each square unit as defined in Equation (1).
2.1.2. Estimated reference sector capacity
Even if the number of aircraft controlled by air traffic controllers were the same in each airspace sector, the workload of each controller may change according to the services given in a particular sector (Majumdar et al., Reference Majumdar, Ochieng and Polak2002). Therefore, estimated reference sector capacity (ERSC) is defined as an acceptable level of ATCW in a sector when conflicts (CDR in Equation 1) are put aside. Consequently, ERSC is determined as the sum of monitoring, coordination and aircraft manoeuvre tasks.
2.1.3. The number of sectors
Determining the optimum number of operational sectors is especially significant for airspace design. Within this framework, we calculated the optimum number of operational sectors, NS, by dividing total airspace workload by ERSC, as defined in Equation (2).
$$NS = \displaystyle{{\sum\nolimits_{i = 1}^n {{\rm W}{\rm L}_i} } \over {{\rm ERSC}}}$$The candidate sectors are built up of the square grid units assigned to a defined centre point. For this reason, the fix points, of which the geographical coordinates are known, have to be determined. The aircraft navigate on airways that are designated by radio navigation aids and navigation fixes. The points at which the radio navigation aids and fixes are located were analysed according to their traffic density, and some of them were chosen as the potential centre points of the candidate sectors.
2.1.4. Distance
Distance is an important parameter in the model in order to assign the centre point of each square unit to the potential centre points of the candidate sectors. All the centre points are defined with their geographical coordinates. WGS-84 (World Geodetic Survey-84) as the geodetic datum and cylindrical equirectangular projection are used to calculate the distance between the centre points of the square unit i and the sector j. These Euclidean distances are calculated in the following form by using GIS:
$$\begin{align} d_{ij} =\sqrt {(x_{i1} -x_{j1} )^2+(\,y_{i1} -y_{j1} )^2} \end{align}$$2.2. Mathematical model for airspace sector design
The parameters, the decision variables and the objectives of the proposed model are defined as follows:
Sets and parameters: Let n be the number of square units, m be the number of potential sectors, 
$I = {\rm \{ }1{\rm ,}2{\rm ,}\; \ldots {\rm ,}n{\rm \} }$ be the set of square units, 
$J = {\rm \{ }1{\rm ,}2{\rm ,}\; \ldots {\rm ,}\;m{\rm \} }$ be the set of sectors, N i be the neighbourhood set of square unit i, i ∈ I, NS be the number of open sectors, WLi be the workload for each square unit i, i ∈ I, d ij be the distance between each square unit i and candidate sector centre point j, i ∈ I, j ∈ J.
Decision variables: Let x ij be 1 if square unit i is assigned to sector j and 0 otherwise; i ∈ I, j ∈ J, y j be 1, if sector j is opened and 0 otherwise; j ∈ J, z be the maximum sector workload and u j be the maximum distance, j ∈ J.
Objective functions: We have the following objectives:
The first objective is to minimise the maximum workload per sector.
$$\begin{align} F_1 (x)= \mathop{\min }_x \mathop{\max }_{j\in J} \sum_{i=1}^n {\mbox{WL}_i } \cdot x_{ij} \end{align}$$The minmax objective can be transformed by including an additional decision variable z, which represents the maximum workload:
$$\begin{align} z= \mathop{\max }_{j\in J} \sum_{i=1}^n {\mbox{WL}_i } \cdot x_{ij} \end{align}$$In order to establish this relationship, the following extra constraint must be imposed:
$$\sum_{i=1}^n {\mbox{WL}_i \cdot x_{ij} } \le z\quad \forall j\in J$$ when z is minimised, these constraints ensure that 
$\sum_{i=1}^{n} {\mbox{WL}_{i} \cdot x_{ij} } $ will be less than or equal to z, 
$\forall j\in J$. At the same time, the optimal value of z will be no greater than the maximum of all 
$\sum_{i=1}^{n} {\mbox{WL}_{i} \cdot x_{ij} } $ because z has been minimised. Therefore, the optimal value of z will be both as small as possible and exactly equal to the maximum workload over j.
The second objective is to minimise the distance between each square unit i and each sector location j. This minimisation can be achieved by minimising the sum of distances for all values of i and j as formulated in Equation (7).
$$\begin{align} F_2 (x)=\sum_i^n {\sum_j^m {d_{ij} \cdot x_{ij} } } \end{align}$$Both objective functions are transformed as follows:
$$\begin{align} f_1 (x) &=z \label{eq8} \end{align}$$
$$\begin{align} f_2 (x) &=F_2 (x) \label{eq9} \end{align}$$Under these notifications, we can formulate a MOMIP model for the airspace sectorisation problem described above, in the following form:
$${\rm Min[}f_1(x),f_2(x){\rm ]}$$Subject to
$$\sum\limits_{j = 1}^m {x_{ij}} = 1\quad \forall i\in {\rm \{ }1, \ldots ,n{\rm \} }$$
$$\sum\limits_{i = 1}^n {x_{ij}} \le n \times y_j\quad \forall j\in {\rm \{ }1, \ldots ,m{\rm \} }$$
$$x_{ij} \le \sum\limits_{r\in K_i} {x_{rj}} \quad \forall i,j\;{\rm and}\;K_i = {\rm \{ }s \vert s \,{\rm is}\, {\rm adjacent} \,{\rm square}\, {\rm unit}\, {\rm to}\, i{\rm \} }$$
$$\begin{align} &\sum_{j=1}^m {y_j } =NS \label{eq14} \end{align}$$
$$\sum\limits_i^n {{\rm W}{\rm L}_i} \cdot x_{ij} \le z\quad \forall j\in {\rm \{ 1,} \ldots ,m{\rm \} }$$
$$x_{ij} = {\rm \{ }0,1{\rm \} }\quad \forall i\in {\rm \{ }1, \ldots ,n{\rm \} }\quad \forall j\in {\rm \{ }1, \ldots ,m{\rm \} }$$
$$y_j = {\rm \{ }0,1{\rm \} }\quad \forall j\in {\rm \{ }1, \ldots ,m{\rm \} }$$The constraint sets (11) and (12) guarantee that each square unit must be assigned to only one sector. The constraint set (13) ensures that all the square units must be connected to each other. K i is a neighbourhood set of square unit i. Hence square unit i cannot assign to sector j unless there exists at least one adjacent square unit of i that is already assigned to sector j. The constraint set (14) ensures that the NS sectors must be opened among m potential sectors defined before. The constraint (15), which is the constraint of the first objective, guarantees that the maximum workload must be minimised.
3. SOLUTION APPROACHES
Multi-objective programming problems are generally solved through different scalarisation methods. The problem is transformed into a single-objective optimisation problem involving some parameters or additional constraints. In this study, the CSM developed by Kasımbeyli (Reference Kasimbeyli2010), see also Gasimov (Reference Gasimov2001), is used. Kasımbeyli introduced an explicit class of increasing convex functions which serve for combining different objectives to a single one, putting aside any convexity and boundedness restrictions on objectives and constraints of the problem under consideration. The CSM guarantees to calculate all properly efficient solutions corresponding to the given weights and reference points. The use of the CSM enables the decision maker to calculate some efficient points which cannot be detected by the other scalarisation methods such as weighted sum scalarisation, epsilon constraint, Benson's scalarisation etc. The prominent features of the CSM are mentioned by Kasimbeyli et al. (Reference Kasimbeyli, Ozturk, Kasimbeyli, Yalcin and Erdem2019). The performances of six scalarisation methods used in multi-objective optimisation are also compared in the same study.
Let w = (w 1, w 2) be a vector of weights, a = (a 1, a 2) be a reference point and α be an augmentation parameter 
${\rm \{ }0 \le \alpha \le \min w_1{\rm ,}\;w_2{\rm \} }$. The parameter α is determined by the decision makers in accordance with their preferences. The objective function used for the CSM in the airspace sectorisation model then will be as follows:
$${\rm Min} f(x) = w_1\cdot {\rm (}f_1(x)-a_1{\rm )} + w_2\cdot {\rm (}f_2(x)-a_2{\rm )} + \alpha \cdot {\rm (} \vert f_1(x)-a_1 \vert + \vert f_2(x)-a_2 \vert {\rm )}$$4. CALCULATIONS
The en-route airspace of Turkey has been selected to test the solution approach presented in this paper. Air traffic data and the sector structure of the Turkish airspace from 2007 are used for the analyses because of the difficulty in acquiring actual data suited to this purpose. The working area used in this study is the airspace above FL 245 (flight level of 24,500 feet), which is defined by Eurocontrol as the upper sector of European airspace. Turkish airspace contains 312 waypoints defined by 64 radio navigation aids, 100 reporting points and 148 fixes (DHMI, 2008). Traffic data containing call signs, types, take-off and departure points of the aircraft, the location of the fix points, and arrival time and flight level of the aircraft on these fixes was provided by the General Directorate of Turkish State Airports (DHMI). Air traffic data collected from Turkish airspace in two peak hours in August 2007 were used in the analyses. Related data were converted to appropriate formats and transferred to the GIS environment. Large-scale data analysis was performed to compute the parameters defined in the model. The sector boundaries of Turkish airspace in 2007 are illustrated in Figure 1.
The sector boundaries of Turkish airspace in 2007.
4.1. Calculation of workload of a square unit
In the first step, the traffic density of Turkish airspace above FL 245 was analysed using GIS. A total of 255 flights were detected during the two-hour period in the defined airspace. In the next step, the selected flight levels of Turkish airspace (>FL 245) were partitioned into 471 square grids with 0·5 degree intervals based on geographical coordinates. The border distance of each square unit, given in Figure 2, is about 55·5 km or 30 nautical miles (NM). In the third step, the control tasks and number of aircraft for each square unit were determined by taking all workload variables defined in the workload measurement formulation into consideration. The speed of the aircraft in accordance with the flight level was obtained from the BADA Aircraft Performance Technical Document (Nuic, Reference Nuic2004). The interpolation method was used to determine the speeds for intermediate values of flight levels. Since data on the directions and the velocities of the wind vector was not available, true air speed (TAS) was used instead of ground speed (GS) in the calculations. The speed and the tracking distance of each aircraft were taken into account for the calculation of aircraft staying time in a sector. The longest staying time of any aircraft is taken as the monitoring time in the study. For conflicts, horizontally 10 NM and vertically 1000 feet between aircraft are defined as the protected zone (DHMI, 2008). Lindberg and Värbrand (Reference Lindberg and Värbrand2001) state that 8 NM separation during cruise translates into approximately one minute. Thus, less than one minute separation between aircraft is considered as a conflict. An altitude change greater than 750 feet, an air speed change greater than 10 knots and a heading change greater than 15 degrees are also considered as the subtasks of aircraft manoeuvre (Laudeman et al., Reference Laudeman, Shelden, Branstrom and Brasil1998).
Turkish airspace partitioning and 28 centre points of candidate sectors.
Although the air traffic control services provided by the controllers are generally the same in all types of airspace, the air traffic patterns and the sector characteristics of each part of the airspace may be different. Individual differences between controllers may also affect the measurement of ATCW. In this framework, a survey was administered to 94 en-route air traffic controllers working in Turkish airspace. The mean values of weighting coefficients and task times which are defined in the workload formulation were obtained from the survey. Since the survey results are not the main subject of this study, further details of the survey will not be given here.
4.2. The calculation of ERSC
We presumed that the three main control services defined in Equation (1), excluding conflicts, are given just once to all aircraft in normal traffic conditions. According to the survey responses, the average number of aircraft which could be safely controlled in a sector was found to be 28 aircraft per hour. Each task consists of subtasks as defined in the workload formulation. The reference weighting factors C and the reference times T were calculated from the arithmetic means of the weighting factors and the times of subtasks. The ERSC during two hours was found to be 95·85 units in the case of the control service provided to 56 aircraft (28 aircraft per hour) in an en-route sector.
4.3. Determination of the centre points and the number of candidate sectors
In the study, the sector boundaries are determined by taking into account both air traffic service routes and Area Navigation (RNAV) routes. As mentioned before, the navigation fixes and the radio navigation aids are chosen as the candidate centre points of the sectors. The intersection points with traffic density of more than 2% are taken into account in the model. The locations of these 28 candidate centre points are given in Figure 2.
The total airspace workload during the two peak hours was found to be 508·775 units by summing the workloads experienced in each square unit. In this direction, the optimum number of operational sectors for en-route levels of Turkish airspace was found to be six, using Equation (2).
$${\rm NS} = \displaystyle{{\sum\nolimits_{i = 1}^n {{\rm W}{\rm L}_i} } \over {{\rm ERSC}}} = \displaystyle{{508\cdot 775} \over {95\cdot 85}} = 5\cdot 31\cong 6 \,s{\rm ectors}$$The results of the survey show that the controllers working in Turkish airspace provide traffic services for 45 aircraft on average per hour during heavy traffic periods. In this framework, the workload capacity of each sector is confined to 45 aircraft per hour. The sector capacity for 45 aircraft during the two peak hours was found to be 155 units without any congestion, and hence the maximum workload in a sector has been limited with this value in the model.
4.4. Determination of weightings and reference points in both scalarisation methods
The model, having two objective functions, is transformed into a structure which has a single objective function by using both the weighted sum scalarisation method and the CSM. The sum of the preference weightings used in both scalarisation methods is determined as 
$\sum_{\dot{\!I}=1}^{2} {W_{\dot{\!I}}} =6$.
One of the important points in the CSM is the determination of the reference points for each objective function. First, the model was solved by the weighted sum scalarisation method, and then the optimum values of two objective functions calculated for different values of the weighting coefficients were compared. Finally, the approximate points to these optimum values were chosen for conic scalarisation as the reference points.
5. EXPERIMENTAL RESULTS
The workload parameters of square units and the coordinates for centre points of each candidate site and each square unit were entered into GAMS software in the form of an Excel file in order to increase the rapidity and reliability of the calculations. The results were then visualised with GIS.
The airspace sectorisation problem was scalarised by the weighted sum scalarisation method making use of arbitrary weighting coefficients for two objective functions. The optimum values of two objective functions 
${\rm (}f_1{\rm ,}\;{\rm }f_2{\rm )}$ and six sectors opened (y j) among 28 candidate sectors calculated for different weighting coefficients are shown in Table 1. Since the problem has non-convexity due to integer variables, there may be solutions which are not detected by weighted sum scalarisation. In this framework, the problem is solved by the CSM with different weightings (w), reference points (a) and augmentation parameters (α ). The alternative solutions are listed in Table 2. The CPLEX solver for weighted sum scalarisation and the DICOPT solver for conic scalarisation were used for the calculations. Although the centre points of the sectors are the same in all analyses performed by both solution methods, the values of workload and distance 
${\rm (}f_1,\;{\rm }f_2{\rm )}$ change according to the parameters chosen. As seen from Table 2, two new optimal solutions are obtained from the CSM analyses.
Computational results obtained by using weighted sum scalarisation.
Computational results obtained by using conic scalarisation.
The locations of the six sectors and their centre points (SP4, SP6, SP13, SP18, SP25 and SP28) among 28 candidate sectors obtained by conic scalarisation for 
$\alpha =2{\cdot}9$, the weightings (3, 3) and the reference points (151, 95,000) are shown in Figure 3. This optimal solution is chosen since it creates a more applicable sector structure. The square units used in the model for the determination of sector boundaries are not suitable for practical usage in ATM. Therefore, the sector boundaries are linearised and improved by connecting the centre points of the corresponding unit squares in the GIS environment. By this means, the workload in each grid is equally distributed to adjacent sectors. The new version of the sector boundaries with airways is demonstrated in Figure 4.
The sector boundaries of Turkish airspace obtained from the model.
Linearised version of sector boundaries and airways of Turkish airspace.
Even though the calculated number of sectors is the same as in 2007, the sector boundaries determined from the proposed model show some differences from the 2007 airspace structure. If workloads in the existing sector structure are compared with that proposed, it is found that the western sectors (Istanbul ACC, Ankara South, Ankara West and Istanbul South) of Turkish airspace have higher workload and more traffic complexity than the eastern sectors (East 1 and East 2). As seen in Table 3, the workload of Istanbul ACC especially exceeded the maximum sector capacity of 45 aircraft per hour defined earlier.
The workload distribution according to the existing and the calculated sector boundaries.
In airspace having nonhomogeneous air traffic density like that of Turkey, it is difficult to distribute workloads equally. Otherwise, this condition necessitates opening more and smaller sectors, which increases airspace complexity and operational costs. Therefore, in the proposed mathematical model, the aim is to distribute the total workload across the sectors in a way that reduces traffic complexity and balances sector dimensions as far as possible. The sector workloads of the proposed structure show that the Turkish airspace is divided into two types of zones each having different traffic complexity and workload levels.
6. CONCLUSION
The balancing of the sector workloads strictly is not always feasible, especially in airspace having nonhomogeneous traffic density like that of Turkey. In this type of airspace, the balancing of workloads among sectors may give rise to the composition of sectors having big differences in size. Even though opening more sectors may be a solution for the balanced distribution of workload, this situation may increase sector complexity and personnel and equipment costs. Therefore, with our MOMIP model, we aimed to design sectors having the same horizontal dimensions as much as possible to increase the flexibility of airspace and to decrease sector complexity.
In this study, an airspace sectorisation model aiming to increase flight safety during any air traffic control activity, to reduce the airspace complexity, and finally to balance the ATCW among the sectors, is presented. The mathematical model developed is applied to en-route levels of Turkish airspace (FL 245 and above).
The most important limitation of the study is the difficulty of obtaining a current and detailed data set related to the aircraft manoeuvres occurring in each sector of the airspace during a one- or two-hour period. The proposed model is tested according to the sector structure of Turkish airspace in 2007 because of the difficulties we experienced in obtaining such an extensive data set. If more recent air traffic statistics could be acquired, the model may be applied to evaluate the current sector structure of Turkish airspace.
As mentioned above, dynamic solutions are generally used for the airspace sectorisation problem. Although this kind of solution generates good results in a short computational time for solving the existing problem, it is not certain that the results obtained from these dynamic methods offer optimal solutions. Therefore, obtaining the optimal solution requires the development of a mathematical model and to solve it by using exact solution methods. The most important disadvantage of the mathematical model approach is that the solution time may be longer than the dynamic solutions due to the increase in the number of intersection points where heavy traffic load is experienced.
The analysis results of the proposed mathematical model, which provides the optimum number of sectors and sector boundaries, may be used not only for the medium- and long-term planning and restructuring of airspace, but also for validating the efficiency of the dynamic solutions. Developing new solution methods and new mathematical models that will reduce the solution time will enable mathematical model approaches to be preferred more frequently in airspace sectorisation problems.
