Jean-Marc Alliot

Nicolas Durand

Géraud Granger Cena

Cena

FACES: a Free flight Autonomous and Coordinated Embarked Solver

Keywords: Free Flight, autonomous aircraft, conflict resolution, air traffic control

FACES is an autonomous and coordinated embarked (on board) conflict solver for Free Flight airspace. It solves conflict by computing simple manoeuvres that guarantees conflict free trajectories for the next 5 minutes (min). Coordination is ensured by giving sequential manoeuvres to aircraft with a token allocation strategy. FACES can be implemented with the current positioning, broadcasting and flight management technology. Moreover, it is robust to communication or system failure for time up to one or two minutes. FACES was tested with a traffic simulator on busy traffic days over France. Airspace over level 320 was considered as Free Flight. 638 out of 641 conflicts were solved without using vertical manoeuvres. The 3 remaining conflicts could easily be solved with vertical manoeuvres. The mean delay was less than 30 seconds by aircraft manoeuvred, with a max delay of 150 seconds. An airborne implementation of this algorithm can be seriously considered.

.

Introduction

We have all experienced at least once a long wait in an overcrowded air terminal. Reading magazines distributed by airlines during these long hours, we often found that they consider air traffic control as one of the major cause for delays. And it is true that the air traffic control system is becoming saturated. But, if delays due to overloaded airports are easy to understand, it is much harder to comprehend delays due to the en route control system. In fact, if we ask a mathematician to analyze the system in cold blood, it can be proved that the collision probability over flight level 320 is very low for aircraft flying direct routes, especially if some elementary precautions are taken regarding face to face or overtaking conflicts. So, en route control could be considered as expensive (en route charges), inefficient (delays induced) and statistically of very little use. However, if the Free Route and Free Flight concepts are attractive, especially to airlines, we still must consider safety as the first priority, and design new algorithms and systems for these new airspaces(airspace above flight level 320 without control).

The most well known reactive collision avoidance concept is certainly the ACAS/TCAS system. It is a very short term collision avoidance system (less than 60 seconds) that is the last safety filter of an ATC system. Experiments with TCAS to control aircraft on simulated traffic have shown that poor coordination could lead to disastrous 1 Centre d'Etudes de la Navigation Aérienne, 7 av Ed Belin, 31055 Toulouse, France, e-mail : alliot@recherche.enac.fr 2 e-mail : durand@recherche.enac.fr 3 Ecole Nationale de L'Aviation Civile, 7 av Ed Belin, 31055 Toulouse, France, e-mail : granger@recherche.enac.fr situations, especially when more than 2 aircraft are simultaneously involved [START_REF] Bosc | Techniques d'évitement réactif et simulation du trafic aérien[END_REF]. Other simple techniques using repulsive forces [2,[START_REF] Zeghal | A comparison of different approaches based on force fields for coordination among multiple mobiles[END_REF] have also been investigated but drawbacks remain [START_REF] Bosc | Techniques d'évitement réactif et simulation du trafic aérien[END_REF]. This paper presents an algorithm for autonomous embarked conflict resolution with a coordination mechanism. Moreover, it is robust to communication or system failure for times up to one or two minutes. This algorithm could be implemented with current technology (GPS, FMS, ADS-B) at low cost.

First section deals with hypothesis we made and the modelling chosen. The second section presents the ordering strategy. In the third section , the A * algorithm used to optimize a conflict free trajectory for one aircraft is detailed. In fourth section , the basic algorithm is tested with the air traffic simulator CATS [START_REF] Alliot | CATS: a complete air traffic simulator[END_REF] on a heavily loaded traffic day in the French airspace over flight level4 320. Experimental results led us to introduce slight improvements to the basic algorithm also presented in the fourth section ; then we finally discuss results of simulations on the enhanced algorithm.

Modelling 1.Hypothesis

The idea is to build an embarked-that is, an on-board-solver able to compute a manoeuvre each time a conflict is detected with another aircraft in a defined detection area around the aircraft. This solver should continuously (every minute 4 ) guarantee a 5 4 minute conflict free trajectory to each aircraft. This 5 minute conflict free period guarantees that a transient failure of communications would not have a disastrous effect: the system could still restart later on; resolutions would be less optimal, more vertical manoeuvres could be necessary to solve all conflicts, as anticipation would be shorter, but the risk of collision would remain close to zero.

Manoeuvres suggested have to be simple to understand and to execute. No manoeuvre can be given during the first minute (called the quiescent period) in order to give enough time to the solver to compute a solution and inform the pilot (or directly program the FMS). Moreover, only one manoeuvre can be given to one aircraft during a 5 minute time window, and no manoeuvre can start as long as the previous one is not finished.

The algorithm enforces a global resolution order between conflicting aircraft. The general principle is as follows: the aircraft which is first chooses its trajectory without considering other aircraft. Then, the next aircraft in the priority queue takes this trajectory into account, and computes its own, and so on (see the second section).

Airspace over flight level 4 320 is considered as a Free Flight airspace. This area is not a so low density area, especially in France. So it is an excellent test zone for a Free Flight solver. All aircraft entering this airspace are supposed to be separated for 5 minutes when entering the Free Flight zone, and are sent back separated for the next 5 minutes when leaving it. All aircraft entering this airspace have to be Free Flight compliant, i.e.:

• they all have synchronous clocks;

• they are able to receive all broadcast information from other aircraft which are within a 90 nmi zone around them (see part three in this section); • they are all equipped with the FACES solver;

• they are able each minute at the same time to compute, and store their current position, their Free Flight airspace exit point and their predicted trajectory for the next 5 minutes;

• they are able to reliably broadcast the latter information as soon as it has been computed. This information consists of 20 3D-points, one every 15 seconds (in fact, only 16 are needed, those beginning at t=1 min). Extra information is added to the predicted position that indicates its accuracy (the uncertainty model is detailed in the third part in this section). Of course, the more accurate the information, the more efficient detection and resolution. This prediction has to be contractual, i.e. as soon as an aircraft has broadcast the information, it has to keep to this trajectory for the next 5 minutes as long as the solver does not give a manoeuvre. It must be noticed that on exceptional occasions, one aircraft can modify this trajectory, or aircraft not equipped for Free Flight can be accepted in the Free Flight zone. This can also take into account exceptional events such as the failure of one aircraft conflict solver. These aircraft will be given the highest priority number (see second part) and all other aircraft will build their trajectory in order to avoid them. This should be a last resort, as the algorithm might fail if two such aircraft are present at the same time in the same zone.

Manoeuvre modelling

As stated above, time is discretized into 15 seconds5 time steps. As manoeuvres must remain simple to understand and execute, the turning point modelling is chosen in the horizontal plane (see Figure 1). In this article, no manoeuvre is given in the vertical plane 6 .

Figure 1: Turning point modelling.

As shown on Figure 1, a manoeuvre is a heading change of 10, 20 or 30 degrees right or left, it starts at time t 0 , and ends at time t 1 on the original path. As stated above, t 0 (and t 1 ) are always larger than 1 minute.

Uncertainty modelling and 1-to-1 conflict detection

A very simple filter is first applied: only aircraft within a 90 nmi zone are considered as being potential threats. This radius is such that 2 aircraft facing each other at 500 kn cannot be in conflict 7 during the next 5 minutes if they are not in the detection zone of the approaching aircraft.

We then assume that there is an error about the aircraft's future location because of ground speed prediction uncertainties 8 . Uncertainties on climbing and descending rates are even more important 9 . Uncertainties on the future positions of aircraft are all the more important because the prediction is faraway. A 5% of uncertainty on ground speed for an aircraft flying at 480 kn represents 2 nautical miles 5 minutes ahead, and about 6 nautical miles 15 minutes ahead. In the vertical plane, 20% of uncertainty for an aircraft climbing at 2000 feet (ft) per minute represents 2000 ft 5 minutes ahead and 6000 ft 15 minutes ahead.

In the vertical plane, we use a cylindrical modelling (Figure 2). Each aircraft has a maximal altitude and a minimal altitude. To check if two aircraft are in conflict, the minimal altitude of the higher aircraft is compared to the maximal altitude of the lower aircraft. points) with the guaranty that two aircraft can not cross each other without noticing a serious conflict. 6 Vertical manoeuvres were put aside on purpose. They are more difficult to execute, and less comfortable for both pilots and passengers. Results of the fourth part show that they should only be used as a last resort, on the very rare occasions where the solver fails. 7 In this article the separation standards are 6 nmi in the horizontal plane and 1000 ft in the vertical plane. 8 Uncertainties on ground track will not be considered, as they do not increase with time and will be included in the separation standard. 9 The error percentages on vertical and horizontal speed are specific to each aircraft. For example, aircraft with very accurate FMS will have very low percentages. In the horizontal plane, an aircraft is represented by a point at the initial time. The point becomes a line segment in the uncertainty direction (the speed direction here, see Figure 2). The first point of the line ``flies'' at the maximum possible speed, and the last point at the minimum possible speed. When changing direction (t=4), the segment becomes a parallelogram that increases in the speed direction. When changing a second time direction (t=7), the parallelogram becomes a hexagon that increases in the new speed direction. To check the separation standard at time t, we compute the distance between the two polygons modelling the aircraft positions and compare it to the separation standard at each time step of the simulation. It must be noticed that, as only one manoeuvre can be given in a 5 minutes time window, and as no manoeuvre can start as long as the previous one is not finished, the convex can only be a line, a parallelogram or an hexagon.

A classical problem in 1 to 1 conflict detection is symmetry. If aircraft A considers it is in conflict with aircraft B, then B must consider A as a conflicting aircraft. In FACES, broadcasting of positions guarantees that two aircraft that can detect each other share exactly the same information regarding their positions. As detection algorithms are identical, 1 to 1 detection will always be symmetrical.

Ordering strategy 2.1 The coordination problem

Centralized automatic solvers as described by N. Durand [START_REF] Durand | Optimisation de Trajectoires pour la Résolution de Conflits en Route[END_REF] find a global solution to clusters involving many aircraft. Manoeuvres are then given to aircraft simultaneously. An on board solver cannot be based on the same principle: aircraft do not share the same information, as they do not have the same detection zone (limited to 90Nmi). A coordination problem appears and must be solved.

The Free-R [START_REF] Duong | Autonomous Aircraft 1st U.S.A/EUROPE ATM R & D Seminar[END_REF] project uses extended flight rules to solve this problem. The TCAS system uses the transponder code to decide which aircraft has to manoeuvre; giving resolution priorities to aircraft is a way often adopted for solving the coordination problem. A resolution priority order 10 has to be total if we want each aircraft to solve all conflicts when there is more than 2 aircraft. For example, the Visual Flight Rule that gives priority to the aircraft coming from the right does not define a global order if there are more than 2 aircraft simultaneously in conflict. Figure 3 gives an example of conflict involving 3 aircraft for which this priority resolution rule does not define an order because transitivity is not ensured.

Moreover, as aircraft do not share the same information, defining which aircraft will be the first to choose its trajectory is not obvious, even if we have a total priority order. Let's take an example. On Figure 4, four aircraft are in conflict. Aircraft 1 is in conflict with aircraft 2 which is in conflict with aircraft 3 which is in conflict with aircraft 4. However, aircraft 1 can only detect aircraft 2 which detects aircraft 1 and 3. Aircraft 3 detects aircraft 2 and 4 which detects aircraft 3. Let's use a simple priority order between these 4 aircraft as for example (1>2>3>4). In conflict [START_REF] Bosc | Techniques d'évitement réactif et simulation du trafic aérien[END_REF]2), aircraft 2 will have to manoeuvre. In conflict (2,[START_REF] Zeghal | A comparison of different approaches based on force fields for coordination among multiple mobiles[END_REF], aircraft 3 will have to manoeuvre and in conflict [START_REF] Zeghal | A comparison of different approaches based on force fields for coordination among multiple mobiles[END_REF][START_REF] Durand | Optimisation de Trajectoires pour la Résolution de Conflits en Route[END_REF], aircraft 4 will have to manoeuvre. However, aircraft 4 does not detect aircraft 2. Consequently, it can not know that aircraft 3 has to manoeuvre before it can manoeuvre itself. Furthermore, aircraft 3 does not detect aircraft 1 and can not know that aircraft 2 has to manoeuvre before it can manoeuvre itself. New trajectories computed independently can then be inconsistent at the end of the resolution step. So, it is not possible to directly derive the resolution order from a simple priority order based only on transponder code for example.

Building a global resolution order 2.2.1 A token allocation strategy

We have to define a global resolution order such so that each aircraft can know when it can start to build a conflict free trajectory and which aircraft it has to avoid.

We first suppose that a total priority order relation exists on the aircraft population. We can use the simple order based on transponder numbers discussed above (but more elaborate orders will be discussed later). At each 10 An order relation must be anti-symmetrical and transitive. An order relation is total if every pair of individuals can be compared. resolution step, we build a global resolution order from this priority order with the following strategy:

1. First, every aircraft sends its predicted trajectory to its neighbours. Each aircraft is then able to know whether it is conflicting with another aircraft or not for the next five minutes.

2. Each aircraft receives a token from every conflicting aircraft that has a higher priority in its detection zone. Aircraft that are not in conflict never receive any token. In the example of Figure 4, aircraft 2 receives a token from aircraft 1, aircraft 3 receives a token from aircraft 2 and aircraft 4 receives a token from aircraft 3.

Figure 5: Cluster of aircraft in conflict

Tokens are allocated as follows:

• A1 gives a token to: A2, A4, A7; 1 gives the token allocation at the different steps of the resolution.

During step 1, A1 and A5 (0 tokens) choose their trajectories without considering other aircraft in their detection area (those that have at least 1 token). Then, they broadcast their unmodified trajectories and all aircraft that have received a token from them cancel it. A2, A4, and A7 cancel the tokens sent by A1. A7 and A8 cancel the tokens sent by A5. During step 2, A2 (0 token) has no token and modifies its trajectory to solve conflict with A1 (0 token). A3, A4, A6, and A7 cancel the tokens sent by A2. During step 3, A3 (0 token) modifies its trajectory to solve conflict with A2 (0 token), A4 (0 token) modifies its trajectory to solve conflict with A1 (0 token); the new trajectory must not interfere with A2 (0 token). A6 cancels one token sent by A3 and A7 cancels two tokens sent by A3 and A4. During step 4, A6 (0 token) modifies its trajectory to solve conflict with A3 (0 token), the new trajectory must not interfere with A2 (0 token). A7 and A8 cancel one token sent by A6. During step 5, A7 (0 token) modifies its trajectory to solve conflict with A2 (0 token) and A5 (0 token); the new trajectory must not interfere with A1 (0 token), A3 (0 token), A4 (0 token) and A6 (0 token). A8 cancels the token sent by A7. During step 6, A8 (0 token) modifies its trajectory to solve conflict with 5 (0 token); the new trajectory must no interfere with A6 (0 token) and A7 (0 token).

Provability

The allocation-resolution method described above cannot lead to situations where all aircraft would have at least one token or situations where two aircraft detecting each other without any token would have to solve simultaneously. This is guaranteed by the use of a total priority order on aircraft. At each step, an aircraft with no token cannot have any other conflicting aircraft (that has not already solved) with no token in its detection area. In such a case, one of these two aircraft would have given a token to the other. At each step, among the conflicting aircraft that have not already solved, there is one that has the highest priority. This aircraft cannot have any token. It can solve and get back its tokens. The algorithm can be mathematically proved.

Communications

The token allocation method gives a mental picture of how aircraft coordinate themselves. Practically, as stated earlier, no bilateral communication is necessary. If aircraft clocks are synchronous, with the information broadcasted, every aircraft can know how many tokens it has and know its neighbour's trajectories. When an aircraft ends a resolution, it broadcasts its new trajectory. Aircraft that have (virtually) received a token can then cancel it.

It is interesting to give a rough estimate of the capacity needed for the communication channel. As stated previously, 16 ``points'' are broadcast. Each point can be as complex as a regular hexagon (see uncertainty modelling), with a lower and a higher altitude. A regular hexagon can be defined by only 4 (x,y) points; if we consider a GPS resolution of 100 meters, given that the earth circumference is 40000 km, x or y value can be coded with 20 bits. For the sake of simplicity, we will also assume that z values are coded on 20 bits. So, a predicted trajectory will consist of at most 16 _ (4 _ 2 +2) _ 20 = 3200 bits. Some other information has to be broadcast, such as the transponder code, the Free Flight zone exit point, etc. This gives a maximum of 3500 bits.

To implement such a system, the ADS-B messages would have to be extended to include trajectory broadcast. But both STDMA and mode-S would provide enough bandwidth. The 1 Mbits/s capacity of mode-S, even divided by 10, would enable enough re-emissions of messages to guarantee a high level of reliability. Regarding STDMA, it solves the garbling problem and the 4000 slots of 150 bits in 1 minute would give 4 slots of 3000 bits in 1 minute to 50 aircraft in the same zone, twice the capacity needed. Of course, a fine modelling of these systems is required to correctly estimate reliability and availability.

The A * algorithm

As soon as the resolution order is chosen, the problem is to solve a 1 to n conflict problem: we have to find the minimum length trajectory for an aircraft avoiding n already fixed aircraft trajectories, that can be considered as obstacles. This is a classical robotics problem, therefore a classical A * algorithm (see [START_REF]Judea Pearl Heuristics[END_REF]) is used. The A * algorithm finds the shortest path in a tree, given an initial state and a set of final states. It uses a ``best first'' strategy to develop each node. The best node is the one that provides the lowest expected cost.

To use an A * algorithm we need:

• u 0 : the initial state;

• T: the set of terminal states;

• p 1 ,p 2 ,..,p n , the set of rules used to build a leaf from a node; (v ← p i u) means that v is built from u using the rule p i ; • k(u,v): if u and v are two nodes, k is the cost of the arc (u,v). The total cost of a state is the sum of the costs of the arcs used to connect the initial state to the current state. In this application, the cost of an arc is the corresponding trajectory length; • h(u): if u is a node, h is a heuristic that estimates the minimal necessary cost to connect the current state to the final state. The efficiency of the A * algorithm depends on h. For a state u, h(u) tries to approach:

h * (u) = Min (v 1 ,v 2 ,..,v n ) [ k(u,v 1 )+k(v 1 ,v 2 )+... +k(v n ,u t )] with u t a terminal state. A heuristic is perfect if: ∀ u, v, h(u)=h(v) ⇔ h * (u) = h * (v) A heuristic is underestimating if: ∀ u, h(u)≤ h * (u)
The performance of the algorithm strongly depends on the quality of the heuristic chosen. With a perfect heuristic, the first path developed is the optimal path; a underestimating heuristic guarantees that the algorithm always converges to the best solution.

In the present application, the initial state is the state of the solving aircraft at t=1 minute. The terminal states are the possible states of the solving aircraft after 5 minutes of flight or when they have reached their destination.

Each branch of the tree represents a possible trajectory of the solving aircraft. Fortunately, the heuristic function is used to only develop a small part of the tree.

The cost of a path is the trajectory length described by this path. Before starting a manoeuvre, an aircraft is in S 0 state. At each time step, each S 0 state generates 6 S 1 states corresponding to the 6 possible deviations of the trajectory (10, 20, 30 degrees right or left), and 1 S 0 state (the aircraft is not manoeuvred). At each time step, each S 1 state generates one S 1 state (the manoeuvre is extended) and one S 2 state (the aircraft is sent back to its Free Flight zone exit point). Every state generates a terminal S 3 state after 5 minutes or if the aircraft has reached its destination.

The cost function k(u,v) measures the distance between the position of the aircraft at node u (time step t tu , state S su ) and the position of the aircraft at node v (time step t tv , state S sv ). If a conflict occurs between node u and node v, the value k(u,v) is widely increased so that the corresponding branch is no longer developed.

The heuristic function h(u) is here the direct distance between node u and the Free Flight exit point (destination) of the solving aircraft. This heuristic is clearly a underestimating one, which guarantees that the optimal solution will always be found.

Figure 6: Example of manoeuvre optimization. Figure 6 details an example of a search in a tree. In this example, the solving aircraft has no manoeuvre at time t=1 minute (t=1 minute corresponds to the beginning of the optimization), it is in S 0 state. The first node chosen is the node that has the lowest estimated cost (h+k). In the example, the lowest estimated cost (100) is found if no manoeuvre is chosen at time 1 minute (step 1). At time 1:15 the lowest estimated cost is still 100 if no manoeuvre is chosen (step 2). At time step 1:30 the lowest estimated cost is 1000 because of a conflict occurring between 1:30 and 1:45 whatever the manoeuvre chosen is. Consequently the branch that has the lowest estimated cost is then developed. Here, at time 1:15 a manoeuvre of 30 degrees left gives an estimated cost of 105 (step 3). The solving aircraft is moved to S 1 state. At the next step, the lowest estimated cost is found if the manoeuvre is extended (107). The alternative choice (S 2 state: turn back to the destination) generates a conflict (step 4). At time 1:45, the lowest estimated cost is found if the aircraft turns back to the destination (S 2 state, step 5). At t=2:00, only one state (S 2 ) can be generated and its estimated cost is still the lowest (step 6).

Generally, many different paths are developed and the depth of the tree is 16 (4 minutes). In this applications, the solution is given in less than 5 seconds on a Pentium II 300, even for the biggest 1-to-35 conflict.

Experimental results and improvements

1-20°S 1 -10°S 0 S1 +10°S 1 +20°S 1 +30°S1 +30°S 1 +20°S 1 +10°S 0 S1 -10°S 1-20°S 1 -30°S1 -30°S 1-20°S 1 -10°S 0 S1 +10°S 1 +20°S 1 +30°S 2 S1 30°S1 +30°S 2 S2
6 t=1:00 t=1:15 t=1:30 t=1:45 t=2:00 t=2:15

The CATS simulator

The algorithm was tested on the CATS [START_REF] Alliot | CATS: a complete air traffic simulator[END_REF] simulator. The core of the CATS system is an en route traffic simulation engine. It is based on a discrete, fixed time slice execution model: the position and speed of aircraft are computed at fixed time steps, usually every 5, 10 or 15 seconds.

Aircraft performances are in tabulated form describing ground speed, vertical speed, and fuel burn as a function of altitude, aircraft type and flight segment (cruise, climb or descent.) Two main datasets for aircraft flight performance are used:

• CAUTRA / ENAC aircraft performance tables, extracted from the French flight data processing system;

• Base of aircraft data (BADA) performance summary tables derived from the total energy model of EUROCONTROL. 69 different aircraft types are described. Synonym aircraft are used to model the rest of the fleet. The Airbus A320 (EA32) is used as default aircraft. In the further applications, aircraft use direct routes to their destination. The separation standard used is 6 nautical miles in the horizontal plane and 1000 ft vertically11 . Conflicts were not solved under flight level 320, and a delay t e was added when necessary for aircraft entering the Free Flight zone in order to separate them on entry points. Uncertainties on speed (either vertical or horizontal) were set to minimal values.

Results using an arbitrary order

Results presented in this part are obtained with the 6381 flight plans of the 21 st of June 1996 with no regulation. The Free Flight zone defined is the airspace above flight level 320. The allocation-resolution strategy described earlier is repeated every minute and the trajectory prediction is done on the next 5 minutes.

2763 aircraft enter the Free Flight zone. 641 conflicts are detected in this zone during the day.

The algorithm requires the definition of a total order among aircraft. The simplest order that can be chosen is an arbitrary order based on the CAUTRA number of each aircraft. This kind of order is already used in the TCAS system.

The A * algorithm is called 2781 times, which represents a mean of 4 times per conflict. At the end of the simulation, 104 conflicts remain. The simulations lasts only a few minutes.

The failing cases were investigated and most of the time the following situation appeared: an aircraft that has a low priority order starts a manoeuvre. A few minutes later, this aircraft is involved in a new conflict while its manoeuvre is not finished. As it still has a low priority order (its CAUTRA number has not changed), it has to solve the conflict. It cannot find a conflict free manoeuvre because the started manoeuvre can not be called into question.

This situation led to the following conclusion: an aircraft that has already started a manoeuvre should have a higher priority order than an aircraft that has not started any manoeuvre. This rule was added to the previous order definition and the simulation was executed again.

Manoeuvre free priority order

The new total priority order defined in this section is the following: an aircraft that is manoeuvre free has a lower priority order than an aircraft that has already started a manoeuvre. The CAUTRA number is used to compare two manoeuvre free aircraft or two manoeuvred aircraft.

With this new priority order, the A * algorithm is called 2654 times. At the end of the simulation, 29 conflicts remain.

Most of the remaining conflicts are due to the horizon effect12 already observed by N. Durand in [START_REF] Durand | Optimal Resolution of En Route Conflicts[END_REF]. The criteria optimized by the A * is the trajectory length between the initial state and the terminal state at the end of the time window (5 minutes). For converging aircraft, the criteria is very often minimized by postponing the conflict out of the time window instead of solving it. Therefore, the solving aircraft is given a heading that tries to move it away from the other aircraft. Figure 7 gives an example of such a conflict: aircraft 3467 is first turned left to move away from aircraft 3660. 1 minute later, aircraft 3660 becomes the solving aircraft and turns right to move away from aircraft 3467. 1 minute later, as both aircraft are not manoeuvre free, the conflict can not be solved. To solve this problem the cost criteria used by the A * algorithm must be changed to take into account the efficiency of a manoeuvre. Therefore, when two aircraft must cross, the manoeuvre that enforces crossing must have a lower cost than the manoeuvre that postpones the crossing. This new criteria is included in the A * algorithm.

Manoeuvre efficiency criteria

In this section, the order previously defined is used and the cost function is modified to take into account the efficiency of the manoeuvre as defined above. The simulation is executed once again. Three conflicts remain unsolved, but they can easily be solved by a very simple vertical manoeuvre.

Delays

897 manoeuvres are given during the day to 367 aircraft (2.44 manoeuvres per aircraft). 13.28% of the aircraft flying in the Free Flight zone are manoeuvred. 73 aircraft are delayed when entering the Free Flight zone, because they are not conflict-free at that moment (that is 2.64% of the traffic in the Free Flight zone). The maximum manoeuvre lasts 2:30 minutes.

Table 3 gives the number of steps aircraft have to wait because they have been given tokens. In the most complex conflict, one aircraft has to wait for 7 resolutions before it can choose its trajectory. So, regarding delays, the performance of the algorithm is very good.

Number of Steps

Unsolved conflicts and priority order

There are 3 aircraft in each remaining unsolved conflict. These conflicts appear because the order between aircraft is not well chosen. F. Médioni [START_REF] Médioni | Méthode d'optimisation pour l'évitement aérien : systèmes centralisés[END_REF] in his PhD thesis showed that a very good solution or no solution at all could be found for a simple conflict involving only 3 aircraft, depending on the order chosen. This situation is described in Figures 8,9 and 10. Moreover, it looks extremely difficult to devise an algorithm that would find the best possible order without seriously increasing the complexity of the global algorithm and the necessary capacity of the communication medium. Embarked conflict solvers which have only a partial information on the global situation will almost certainly remain suboptimal, while centralized conflict solvers are able to find the global optimal solutions. However, this may not be a too serious concern in the upper airspace: the simulation above shows that this algorithm is almost always able to solve conflicts, even with situations as complex as the one presented on on Figure 11 where 35 aircraft are involved, while delays remain small.

Conclusion

We have demonstrated in this paper that an efficient on board algorithm for Free Flight conflict resolution can be designed and implemented. This algorithm has the following advantages:

• Compared to a centralized automated system, the development of such a system could be relatively low cost. Most assumptions are quite weak: synchronous clocks are already available with GPS, FMS are now elaborate enough to provide the information needed for trajectory prediction in the next 5 minutes, ADS-B, or even better STDMA, provides sufficient capacity for communications; the 1 to n resolution algorithm is simple to implement and has already been widely used for similar problems in robotics; computing power needed fits in a standard PC computer. • Compared to rule based systems, the algorithm is mathematically provable, and the simulation above shows that it would be efficient in upper airspace, even when the density is quite high, and even with quite a large separation standard: without vertical manoeuvres, only 3 conflicts out of 641 remain. 13% of the aircraft are manoeuvred and the mean delay is 27 seconds by aircraft effectively manoeuvred, with a maximum delay of 150 seconds.

• Compared to purely reactive systems [2], which usually require constant changes in headings, the manoeuvre model is classical and easy to implement. Further, and this is the main point to stress, as trajectories are guaranteed conflict free for at least 5 minutes, a transient failure of communications would 4956 160 329= 450 1347 131 330= 470 not have a disastrous effect: the system could still restart later on; resolutions would be less optimal, more vertical manoeuvres could be necessary to solve all conflicts, as anticipation would be shorter, but the risk of collision would remain insignificant.

• The system could be progressively put into service by first defining Free Flight airspace over oceanic areas, and gradually extending them. This would help in solving the classical transition problem from the current system to a partially automated one. The main issue that has to be refined is the value of the quiescent time window (set to one minute in our simulation). During this period, each aircraft has to build its predicted trajectory and broadcast it. Then the whole loop of the algorithm has to be executed, with resolutions, new trajectories broadcast, etc. It is still unclear if one minute will be enough. Simulations show that one minute is a correct value when communications are instantaneous. A fine modelling13 of the communication medium is needed to confirm or deny this value. However, simulations could be conducted with larger values to test the behaviour of the algorithm.

Simulations have also to be conducted using larger uncertainty parameters. As trajectory prediction is done for the next 5 minutes only, results should not be very different. We will also test the solver with higher traffic to find out its limits regarding density. Vertical manoeuvres will be added as last resort to correctly complete the resolution step.

Figure 2 :

 2 Figure 2: Modelling of speed uncertainties.

Figure 3 :

 3 Figure 3: 3 aircraft conflict.

Figure 4

 4 Figure 4: 4 aircraft cluster.

Figure 7 :

 7 Figure 7: Second step: aircraft 3660 turns right

Figure 9 :

 9 Figure 8: order (1,2,3), optimal solution

Figure 10 :

 10 Figure 10: order (2,1,3), no solution.

Figure 11 :

 11 Figure 11: 35 aircraft in the detection zone.

Table 1 :

 1 Token allocation at the different steps of resolution

	• A2 gives a
	token to: A3,
	A4, A6, A7;
	• A3 gives a
	token to: A6,
	A7;
	• A4 gives a
	token to: A7;
	• A5 gives a
	token to: A7,
	A8;
	• A6 gives a
	token to: A7,
	A8
	• A7 gives a

Table 2 :

 2 Delays

	3467 191
	330= 470
	3660 160
	330= 480

Table 3 :

 3 Number of waiting aircraft.

These parameters can be modified. This first study does not discuss the opportunity of increasing or decreasing these values.

This value is not chosen at random. With 15 s time steps, detection can be made only on these points (and not on the segments between these

Then, each conflicting aircraft with no token solves conflicts with every aircraft in its detection zone that has no token. It does not take into account aircraft that have one or more tokens.

When this trajectory has been computed, the aircraft broadcasts its new trajectory; all aircraft which have received a token from this aircraft take this new trajectory into account, and cancel the token received from this aircraft.

Steps 2 and 3 are repeated until no token remains.On Figure4, aircraft 1 does not change its trajectory as it has no zero token aircraft in its detection zone. It broadcasts its new (unmodified) trajectory. So aircraft 2 takes this trajectory into account and cancels the token received from aircraft 1. Then aircraft 2 solves its conflict with aircraft 1 as aircraft 1 is in its detection zone and has no token. Aircraft 3 takes this trajectory into account, cancels token from aircraft 2 and solves the conflict with aircraft 2. Aircraft 4 then applies the same procedure.

The 6 Nmi value is certainly quite high for an embarked solver. We do think that it could seriously be reduced, regarding GPS and FMS precision. This would even improve airspace capacity and resolution efficiency.

This effect has been observed and discussed for years now in other tree search algorithms: game algorithms.

This modelling would also be very useful to find out the exact availability, capacity and error rate of the channel.

We are conscious that the whole system depends on the reliability and availability of transmissions. Requirements on the bandwidth are low enough to enable multiple emissions of messages. But error correlations would have to be considered. Specific data and results on these issues do not appear to be presently available. However, we believe that an airborne implementation of this algorithm can be seriously considered.

Acronyms


Jean-Marc Alliot
email: alliot@recherche.enac.fr

Nicolas Durand
email: durand@recherche.enac.fr

Géraud Granger
email: granger@recherche.enac.fr

A statistical analysis of the influence of vertical and ground speed errors on conflict probe

Keywords: air traffic control, conflict probe, statistical analysis

Conflict probes will be important components of air traffic control tools in the next years. They appear in almost every project (CINCAT, ERATO, HIPS, PHARE, and so on). To be useful, they have to fulfill two goals : reliability (they have to detect all conflicts) and efficiency (they must minimize the number of false alarms). Conflict probes rely on trajectory prediction, and their reliability and efficiency highly depend on the accuracy of trajectory prediction. In this paper, we present a quick mathematical insight of the influence of ground speed errors on trajectory prediction, results of arithmetic simulations on real traffic both on ground and vertical speed errors, and a statistical analysis of these results to model the influence of vertical and ground speed errors on conflict probe.

A mathematical overview

In this section, we present quickly some mathematical results regarding two aircraft conflict.

Two aircraft conflict at constant speed

Figure 1 show a classical two aircraft conflict. Aircraft on the lower segment fly at speed v 1 , and aircraft on the upper segment fly at speed v 2 . The angle of incidence is α.

We will use the auxiliary variables r = v 2 /v 1 , and D the separation standard. Let's suppose that we have an aircraft p 1 on the lower segment at a distance l 1 of the crossing point. We want to know which interval on the

❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❛ ❡ ❡ ✪ ✪ ✟ ✟ ✛ ✲ α v 1 v 2 l 1
Figure 1: Two aircraft conflict upper segment will contain conflicting aircraft with this one. Let's assume that an aircraft p 2 is at distance l 2 of the crossing point. Then we have (with a classical orthonormal referential):

x 1 = v 1 t -l 1 y 1 = 0 x 2 = cos(α)(v 2 t -l 2 ) y 2 = sin(α)(l 2 -v 2 t)
If the two aircraft are to be in conflict, there must exist t such that the following inequality is satisfied :

(x 1 -x 2 ) 2 + (y 1 -y 2 ) 2 ≤ D 2
This is a second degree inequality in t. We notice that for t → +∞ or t → -∞, the inequality is not satisfied. So, it will only be satisfied if the equation:

(x 1 -x 2 ) 2 + (y 1 -y 2 ) 2 -D 2 = 0
has at least one root. If the discriminant of the above equation

∆ = D 2 (v 2 1 -2v 1 v 2 cos(α) + v 2 2 ) -sin(α) 2 (l 2 v 1 -l 1 v 2 ) 2
is negative, then aircraft p 2 will never be in conflict with p 1 ; if the discriminant is positive, they will be in conflict during [t 1 , t 2 ] where t 1 and t 2 are the roots of the equation.

The discriminant is itself a polynomial of degree 2 in l 2 . Thus, the distances l 2 which satisfy ∆ >= 0 belong to a single interval, and the extremal points of the interval are the roots of the equation ∆ = 0

r 1 = l 1 v 2 v 1 + D 1 + ( v2 v1 ) 2 -2( v2 v1
) cos(α) sin(α)

(1)

r 2 = l 1 v 2 v 1 -D 1 + ( v2 v1 ) 2 -2( v2 v1 ) cos(α) sin(α) (2)
However, we are only interested in the length of the segment, which is the difference between the two roots. We do not take into account the degenerate cases α = 0 (takeover) and α = π (facing aircraft). It comes:

L = 2D 1 + (v 2 /v 1 ) 2 -2(v 2 /v 1 )cos(α) sin(α)
We take r = v 2 /v 1 as the ratio of speed and we express the length of the segment in number of separation standard:

L D = f (r, α) = 2 1 + r 2 -2r cos(α) sin(α)
The function doesn't depend on l 1 (it could have been expected) and the minimum is reached for r = cos(α), and we have f (cos(α), α) = 2. This is independent of α and perfectly normal: it represents aircraft p 2 at a distance -D < d < D of the crossing point when aircraft p 1 is exactly at the crossing point. We represent on figure 2 a contour plot of the above function. The x axis is r (the ratio of v 2 /v 1 ); we use values of r ranging from 0.5 to 1.5. The y axis is the angle of incidence α in degrees, from 20 degrees to 150 degrees. We do not represent values above 150 degrees or below 20 degrees, as they exhibit pathological behavior. These curves are not exactly new, as quite similar ones already appeared in documents describing the Gentle Strict (GS) algorithm of the AREA project[NFC + 83, Nie89b, Nie89a]!

The darkest part of the contour plot represents segment lengths ranging from 2 to 2.5 separation standard. Each line represent one half separation standard more.

One should naively expect the function to be symmetrical with f (r, α) = f (1/r, α), as the problem looks symmetrical: aircraft on the upper segment should "see" the same number of conflicts with aircraft on the lower segment than aircraft on the lower segment with aircraft on the upper segment. Actually, we have:

f ( 1 r , α) = 1 r f (r, α).
This can be intuitively understood on a simple example.

Let's suppose that aircraft on the upper and the lower segment are equally spaced by exactly one separation standard (maximal rate). Each aircraft on the lower segment will detect roughly f (r, α) conflicts; in a given time period T , n 1 = T × ( v1 D ) aircraft will pass on the lower segment and the total number of conflicts detected by aircraft on the lower segment will be n 1 f (r, α). During the same time period, aircraft on the upper segment will detect n 2 f ( 1 r , α) conflicts. These two numbers have to be equal, so:

n 1 f (r, α) = n 2 f ( 1 r , α), or f (r, α) = n2 n1 f ( 1 r , α)
. By replacing n 1 and n 2 by their values:

f (r, α) = T × ( v2 D ) T × ( v1 D ) f ( 1 r , α) = v 2 v 1 f ( 1 r , α) = rf ( 1 r , α)

Two aircraft conflict with uncertainties

We suppose now that speeds are not exactly known, but that they belong to a given interval, with a percentage of error. We have then:

x 1 = v 1 (1 + e 1 )t -l 1 y 1 = 0 x 2 = cos(α)(v 2 (1 + e 2 )t -l 2 ) y 2 = sin(α)(l 2 -v 2 (1 + e 2 )t)
By taking the following variable: r = v2 v1 , l = l1 D and k = 1+e2 1+e1 , we just have to substitute in equations 1 and 2 to get the new extremal points of the segment:

r 1 D = lkr + 1 + r 2 k 2 -2rk cos(α) sin(α) r 2 D = lkr - 1 + r 2 k 2 -2rk cos(α) sin(α)
If we suppose that both e 1 and e 2 belongs to the same interval [-e, e], we would like to plot the length of the segment as a function of e:

L D (α, l, r, e) = max k∈[ 1-e 1+e , 1+e 1-e ] r 1 D - min k∈[ 1-e 1+e , 1+e 1-e ]
r 2 D

Let's consider the maximization of r 1 (k). It is a constrained maximization problem, so the extremum will be either on the hull of the convex (here at one of the points 1-e 1+e or 1+e

1-e ) or at an extremal point in the convex. This last point will be found by solving the equation:

∂r 1 ∂k = l r + k r 2 -r cos(α) sin α 1 + k 2 r 2 -2 k r cos(α) = 0
This equation has only one root in k:

k r = -l 1 -l 2 sin 2 α sin 2 α + (1 -l 2 sin 2 α) cos α r(1 -l 2 sin 2 α)
and this root only exists when l < 1 sin α . The second order derivative is:

∂ 2 r 1 ∂k 2 = r 2 sin(α) (1 + k 2 r 2 -2 k r cos(α)) 3 2 > 0
We have also:

∂r 1 ∂k (k = 0) = r(l - 1 sin α ) lim k→∞ ∂r 1 ∂k (k) = r(l + 1 sin α )
Then we can conclude that the extremum of the function will never be reached inside the interval, but only at the upper or lower bound. Then:

l > 1 sin α : max k∈[ 1-e 1+e , 1+e 1-e ] r1 D = r1( 1+e 1-e ) D l < 1 sin α : max k∈[ 1-e 1+e , 1+e 1-e ] r1 D = max r1( 1+e 1-e ) D , r1( 1-e 1+e ) D
An extremely similar discussion arises from the minimization of r 2 . The conclusion is: We can now plot the length of the segment L D (α, l, r, e). However, it is much more interesting here to represent the ratio of the number of conflicts detected (which is proportional to L) to the number of conflicts which would be detected if there was no uncertainty. So we do plot the function L(α,l,r,e) L(α,l,r,0) on figure 3. The x-axis is the uncertainty e with values from 0 to 30%; the y-axis is l, and is also the time to the crossing point in minutes if we make the approximation that an aircraft flies one separation standard in one minute; values range from 3 to 10 minutes to the crossing point. We set r = 1 and α = π/4. The darkest part (on the left) is the zone where we have from 0 to 10% conflict more. Each line crossed towards the lightest zone adds 10% of conflict.

l > 1 sin α : min k∈[ 1-e 1+e , 1+e 1-e ] r2 D = r2( 1-e 1+e ) D l < 1 sin α : min k∈[ 1-e 1+e , 1+e 1-e ] r2 D = min r2( 1+e 1-e ) D , r2 (
For example, if we detect conflicts 5 minutes before the crossing point, with an angle of convergence of 45 degrees, a speed ratio of 1 and an uncertainty of 15%, we will detect 50% more conflicts than if uncertainty was 0.

Notice that the time used is not the time before the start of the conflict but before the crossing point is reached. The difference is quite small for large values of α, but large for small values of α.

If we suppose that we detect conflicts far enough from the crossing point (a reasonable assumption), the maximum of r 1 is r 1 ( 1+e 1-e ) and the minimum of r 2 is r 2 ( 1-e 1+e ). Then we have:

L D (α, l, r, e) = r 1 ( 1+e 1-e ) D - r 2 ( 1-e 1+e ) D
We can develop this equation:

L D = 4 e r l 1 -e 2 + 1 sin(α) 1 + (1 -e) 2 r 2 (1 + e) 2 - 2 (1 -e) r cos(α) 1 + e + 1 + (1 + e) 2 r 2 (1 -e) 2 - 2 (1 + e) r cos(α) 1 -e
It is interesting to notice that the second term doesn't contain l and that for large values of l, this function is equivalent to 4 e r l 1-e 2 . For values of e not too large, we can write the order 1 Taylor series in e:

L D (α, r, l, e) = 4 l r e + 2 1 + r 2 -2 r cos(α) sin(α)
and find back in the second term the length of the segment when e = 0 that was computed in the previous section. For that case, the increase in the number of conflicts is:

L(α, r, l, e) L(α, r, l, 0) = 1 + 2 r sin(α) 1 + r 2 -2 r cos(α) l e = 1 + C(r, α) l e
This last expression is simple and interesting. The 1 represents the "standard" number of conflicts, and the second term is the increase in percentage. Figure 4 is a contour plot of C(r, α) for r ∈ [0.5, 1.5] and α ∈ [20, 150] degrees. The white zone corresponds to 1.9 < C(r, α) < 2, and each line crossed towards the darkest zone decreases C(r, α) by 0.1. The general shape of the curve is not surprising: for very small (≃ 0) or very high (≃ 180) values of α, conflicts are almost certain, so the percentage of "false alarms" remains small. Regarding the speed ratio, aircraft having quite similar speeds are intuitively more prone to give false alarms than aircraft with clearly different speeds. On a simple numerical example (r = 1, α = 30), we have:

L(α, r, l, e) L(α, r, l, 0) ≃ 1 + 2 l e
where l can be considered as the time in minutes to the crossing point (with, again, the approximation that an aircraft fly one separation standard in one minute), and e is the uncertainty. for very small values of e (e < 0.05, i.e. 5%), the second term remains less than 1, even for l up to 10 minutes. However, for larger values of e, the function looks like 2 l e, a direct proportionality to l and e.

The last two curves (figure 3 and4) are extremely interesting. Acceptance of conflict probe systems depends highly on their efficiency, and this efficiency can be evaluated as the percentage of "false alarms". The theoretical model shows that, to be able to detect conflicts soon and efficiently, we must be sure that the uncertainties on speed are as small as possible, and certainly less than 5%.

A statistical analysis 2.1 Introduction

The mathematical model makes lot of assumptions (constant speeds for example), many of them being quite unrealistic. Moreover, we haven't taken into account the vertical plane. We are now going to try to find a statistical model of the number of conflicts detected depending on the value of vertical and ground speed errors and of the value of the anticipation (the time window used for conflict detection).

To get data for our statistical analysis, simulations were used. The following subsections describe the experimental protocol.

The Air Traffic simulator

We used the OPAS [DAB97, DABM97] simulator (sometimes known as CATS) for simulations. It uses a tabulated model for aircraft performances: ground speed, vertical speed, and fuel burn are functions of altitude, aircraft type and flight segment (cruise, climb or descent.) The main dataset for aircraft flight performance is the base of aircraft data (BADA) performance summary tables derived from the total energy model of EUROCONTROL. 69 different aircraft types are described. Synonym aircraft are used to model the rest of the fleet. The Airbus A320 (EA32) is used as default aircraft.

Aircraft follow classical routes (from way-point to waypoint). The flight model is simple: an aircraft first climbs up to its RFL, then remains leveled till its top of descent, then descends to its destination.

Aircraft fly with a timestep that can be chosen at the start of the simulation. The timestep is always chosen in order to guarantee that two aircraft face to face flying at 500 kts could not cross without being closer than one standard separation at at least one timestep. For most of our simulation, we use a 15s timestep.

Flight plans are data of the COURAGE system, an archiving system of the operational French CAUTRA Air Traffic Control system. We have been using initial flight plans (without regulation), and we have used one of the module of the OPAS simulator to give slots to aircraft.

Conflict detection

Trajectory prediction is done each three minutes by a simulation of a given duration inside the global simulation. This duration is what we call the anticipation.

We assume during each of these detection simulations that there is an error about the aircraft future location because of ground and vertical speed prediction uncertainties.

Then, an aircraft is represented by a point at the initial time of the conflict detection window. But the point becomes a line segment in the uncertainty direction (the speed direction here, see figure 5). The first point of the line "flies" at the maximum possible speed, and the last point at the minimum possible speed. These maximal and minimal speeds depend of course on the uncertainty chosen: for 5% uncertainty on ground speed, the first point will fly at a speed of 1.05 v and the last point at 0.95 v, if v is the nominal speed of the aircraft.

When changing direction on a beacon, the heading of the line segment "fastest point" changes as described on figure 5.

To check the standard separation at time t, we compute the distance between the two line segments modeling the aircraft positions and compare it to the standard separation at each timestep of the simulation.

In the vertical plane, we use a cylindrical modeling (figure 5). Each aircraft has a mean altitude, a maximal altitude and a minimal altitude. To check if two aircraft are in conflict, the minimal altitude of the higher aircraft is compared to the maximal altitude of the lower aircraft. The separation standard used is 6 nautical miles in the horizon- Conflicts detected can be merged: if a conflict is detected a time t 1 , and detected again three minutes later, the two conflicts are only considered as one. We maintain a hash table during the whole simulation to determine which conflicts are to be merged. When using direct routes, the model is slightly modified (see figure 6) and much more simple. As aircraft never change direction, the aircraft is, at the beginning of the time window, a simple point, which becomes a line segment whose size grows with the uncertainty on speed. The modeling in the vertical plane is identical to the one used for standard routes.

t=0 t=1 t=2 t=3 t=5 t=6 t=4 ✁ ✁ ✁ ✁ ✁ ✁ ✁ ✁ ✁ ✁ ✁ ✁ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✂ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ✄ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ☎ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✆ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✝ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✞ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✟ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✠ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ✡ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛ ☛

Beacons

Results

We ran 540 simulations with anticipation ranging from 4 minutes to 12 minutes, vertical speed error ranging from 1 to 40% and ground speed error ranging from 2 to 16%. In table 1, we give a very short sample of data found by simulations.

t w is the duration of the time window in minutes, e g is the horizontal uncertainty and e v is the vertical uncertainty. N d is the number of conflicts measured for these values for direct routes and N s for standard routes.

It must be noticed (these results are not shown in the table) that when t w = 0, N = 1726 for every value of e g and e v (when no predictions are made, there can't be more conflicts), and that N = 1726 also for all values of t w when both e g and e v are equal to 0 (without uncertainties, conflicts can be detected from the beginning).

Statistical modeling (direct routes)

The problem when doing statistical analysis is to fit the data into the models that can be built and evaluated. In the first two subsections, we describe simple multiplicative models whose parameters can be easily found with linear regression (by using the logarithm), while in the third section we present a more complex model, requiring more elaborate methods to find the parameters.

A simple multiplicative model

In order to be able to use simple statistic methods, we tried the following modeling:

N = N 0 + a 0 t a1 w (1 + e g ) a2 (1 + e v ) a3
t w , e g , e v and N are the variables described in the subsection above, while N 0 is the number of conflicts when t w = e g = e v = 0.

Using logarithms, we get:

log(N -N 0 ) = log(a 0 ) + a 1 log(t w ) + a 2 log(1 + e g ) + a 3 log(1 + e v )
Estimate SE Tstat PValue log(a 0 ) 5.068 0.0278 182 0 a 1 0.917 0.0122 75 0 a 2 10.78 0.1225 88 0 a 3 2.075 0.0444 47 0 R 2 = 0.946689

Table 2: N = N 0 + a 0 t a1 w (1 + e g ) a2 (1 + e v ) a3
We can then use a linear regression on all our data samples. These results are summarized in table 2. Results look good if we consider the statistical estimators presented [START_REF] Saporta | Probabilités, analyse des données et statistique[END_REF]: good R 2 , high T s for all variables with a low probability value. However, this model works on the logarithm of all variables. An other standard, but more interesting (and intuitive), estimator of the validity of the model is:

S e = i=no i=1 (N i -Ni ) 2 n o
where n o is the total number of observations, N i is the number of conflicts observed, and Ni is the number of conflicts predicted by the model. S e can be "interpreted" as a quadratic mean of the differences between observed values and predicted values.

For that model, we have S e = 487, which is not a so small value if we consider that the actual number of conflicts range from 1726 to 14216.

Another indicator is:

s e = i=no i=1 (1 -Ni Ni ) 2 n o
The S e indicator has the disadvantage of giving more weights to large values, while this one estimates a global "percentage" of error. For this model, we have s e = 0.12, a quite large value again. Moreover, this model has a major problem. With both e g and e v equals to zero, we should get N = N 0 , but this is not the case.

Another multiplicative model

Another model is:

N = N 0 + a 0 t a1 w e a2 g e a3 v
Using logarithms, we get:

log(N -N 0 ) = log(a 0 )

A multiplicative and additive model

The two models above both exhibit congenital problems.

We would like to have a model which gives N = N 0 only when t w = 0 with any value of e g or e v , or when e g = 0 and e v = 0 with any value of t w . The simplest model having these two features is:

N = N 0 + a 0 t a1 w (e a2 g + a 3 e a3 v )
But this model is both multiplicative and additive: we can not use standard regression methods to find the optimal values of the parameters. The usual method is then to try to find the set of parameters minimizing the S e or the s e functions described above. We finally decided to optimize s e .

We used both global optimization (interval arithmetic and branch and bound) and local optimization (simplex) to find the optimal set of parameters [START_REF] Hansen | Global optimizations using interval analysis[END_REF]. We finally get (when dividing by N 0 ): For that set of parameters and that model, we have s e = 0.048, an excellent value, much lower than the values found for the other two models.

N N 0 = 1 + t 0.
When e v = 0, the increase of the number of conflict given by this model is almost identical to the formula found in the theoretical part of this paper. Exponents of t w and e g are extremely close to 1.

We notice that the exponent of e v is also almost equal to 1, but this term has a lesser contribution, only 1/7 of the e g contribution. We can try to explain this result: aircraft are stable during a large part of their flight, so the vertical uncertainty only has an influence during a smaller part of the flight than the horizontal uncertainty.

Statistical modeling (standard routes)

In this section, we consider the problems on standard routes. We only compute results for the best model, i.e. the additive and multiplicative one.

A model for standard routes

We now apply the following model:

N = N 0 + a 0 t a1 w (e a2 g + a 3 e a3 v )
to the results found on standard routes. We used exactly the same methods than above to find the parameters of the model. Then, we have: We have here s e = 0.058, an excellent fit again. We are still quite close to the theoretical model, but we must note some differences: the exponents are slightly larger than 1 (closer to 1.2), but the most striking fact is the value of the constant of the ground speed error. While it remains the same for vertical speed error (around 0.5), it is more than twice the direct route model value for ground speed error (8.05 instead of 3.66). This is easy to understand: as aircraft fly on routes, there are much more problems of takeover and face to face.

N N 0 = 1 + t 1.

Conclusion

In this paper, we have presented a theoretical model and a statistical analysis of the influence of the vertical and ground speed errors on conflict probe.

The formulas which summarizes this study can be written approximatively:

N d N 0
= 1 + t w (3.5 e g + 0.5 e v )

N s N 0 = 1 + t 1.2 w 8 e 1.2 g + 0.5 e 1.2 v where N N0 represents the ratio of conflicts detected over conflicts really happening (N d is the formula for direct routes, and N s the formula for standard routes), t w being the prediction anticipation, e g being the ground speed error (in percentage) and e v being the vertical speed error (in percentage also).

This model shows that, if we want to detect conflicts 10 minutes before they appear, and accept to detect twice the actual number of conflicts, we would need, for example to have e g < 0.014 = 1.4% and e v < 0.1 = 10% on direct routes; these values are currently out of reach of ground based trajectory prediction systems without using FMS informations.

The study also shows that improving significantly trajectory prediction will drastically reduce the number of conflicts detected, thus giving the controller a significant increase in his comfort, and probably a significant increase also in sector capacity.

These might be tracks that could be followed in the near future.

Biography

Jean-Marc Alliot graduated from the Ecole Polytechnique de Paris in 1986 and from the Ecole Nationale de l'Aviation Civile (ENAC) in 1988. He also holds a Ph.D. in Computer Science (1992). He is currently in charge of the global optimization laboratory of CENA and ENAC in Toulouse.

Nicolas Durand graduated from the Ecole Polytechnique de Paris in 1990 and from the Ecole Nationale de l'Aviation Civile (ENAC) in 1992. He has been a design engineer at the Centre d'Etudes de la Navigation Aérienne (CENA) since 1992 and holds a Ph.D. in computer Science (1996).

Géraud Granger graduated from the Ecole Nationale de l'Aviation Civile (ENAC) in 1998. He is completing a Ph.D. in Computer Science at the Ecole Polytechnique de Paris.

Figure 2 :

 2 Figure 2: Two aircraft conflict

Figure 3 :

 3 Figure 3: percentage of additional conflicts (r=1, α = π/4)

Figure 4 :

 4 Figure 4: C(r, α)

Figure 6 :

 6 Figure 5: Modeling of speed uncertainties (standard routes).

Table 1 :

 1 Results of simulations

	t w	e g	e v	N d	N s
	0 0.00 0.00	1726	2758
	4 0.02 0.01	2211	4001
	4 0.02 0.40	3323	6264
	4 0.16 0.01	4501 13653
	4 0.16 0.40	6006 17847
	12 0.02 0.01	2937	6474
	12 0.16 0.40 14216 68200

Table 3 :

 3 +a 1 log(t w )+a 2 log(e g )+a 3 log(e v ) N = N 0 + a 0 t a1 w e a2 g e a3 v S e and s e gives S e = 453 and s e = 0.12, values extremely close to the one found for the other multiplicative model. This model corrects the problem above, but a new one arises: with this formula, having e g = 0 or e v = 0 will give N = N 0 , which is, of course, incorrect.

	Linear regression results are summarized in table 3. There
	again statistical results look good, but the computation of


Jean-Baptiste Gotteland
email: gotteland@recherche.enac.frdurand@tls.cena.fr

Nicolas Durand

Jean-Marc Alliot
email: alliot@dgac.frpage@cena.fr

Erwan Page

Aircraft Ground Traffic Optimization

Air traffic growth and especially hubs development cause new significant congestion and ground delays on major airports.

Accurate models of airport traffic prediction can provide new tools to assist ground controllers in choosing the best taxiways and the most adapted holding points for aircraft. Such tools could also be used by airport designers to evaluate possible improvements on airport configurations and airport structure.

In this paper, a ground traffic simulation tool is proposed and applied to Roissy Charles De Gaulle and Orly airports. A global optimization method using genetic algorithms is compared to a 1-to-n strategy to minimize time spent between gate and runway, while respecting aircraft separation and runway capacity.

In order to compare the efficiency of the different optimization methods, simulations are carried out on a one day traffic sample, and ground delay due to holding points or taxiway lengthening is correlated to the traffic density on the airport.

Introduction

Traffic delay due to airport congestion and ground operations becomes more and more penalizing in the total gate-to-gate flight cycle. This phenomenon can be in a large part attributed to recent hubs development, as all departures and arrivals are tending to be scheduled at the same time. Moreover, many ATC problems and environmental inefficiencies can appear as a result of taxi queueing and take-off time uncertainty. As airport designers are in charge to build new taxiways to reduce congestion and improve ground operations, ground simulation tools become essential to validate their choices before realization.

Even if most research projects are concentrated on decision making tools for airspace controllers and do not consider ground operations utilities, highly detailed models of airport operations already exist, such as SIMMOD1 or TAAM2 . They can be useful to evaluate qualitatively the relative effects of various airport improvements. The DP3 project ([IDA 98]) focuses on improving the performance of departure operations. Finally, a component of the TARMAC4 project focuses on the ATC-related traffic planning systems for airport movements.

In this paper, a ground traffic simulation tool with a conflict resolution module is introduced and tested on a one day traffic sample on Roissy Charles De Gaulle and Orly airports. Different optimization strategies are used to find the best trajectory and the most adapted holding points for taxiing aircraft. The goal is to minimize the time spent from gate to take-off or from landing to gate, respecting the separation with other aircraft and the runway capacity. During the optimization process, actual one way taxiways, operational airport configurations and speed uncertainty are considered.

Problem modeling

The problem is to find, for each aircraft, an optimal path from its gate to a given runway take-off position or from its runway exit to its gate position, respecting a given separation between aircraft.

An optimal path can have different definitions: for example, the length of the path or the total taxiing time. At the same topic, holding on a taxiway can be more or less penalizing than increasing the length of the path or holding at the gate position. Therefore a global optimum criteria will have to be defined in the following. However, the purpose of this article is not to discuss the choice of such criteria, which can be refined without modifying the algorithm itself, considering many different factors related to the airport geometry, the traffic, or airlines preferences... By the way, it is quite difficult to predict with a good accuracy the future positions of aircraft on taxiways. First of all, the exact departure time is generally known only a few minutes in advance (many factors can cause delays), and the exact landing time depends on the runway sequencing. Hence, the proposed model should take into account speed uncertainty and must be regularly updated with real aircraft positions.

Airport structure

An airport is described by its gates, taxiways and runways. Different kind of taxiways can be differenced: ✁ Gate specific access (entries, forward exits or push-backs), characterized by a very low speed; ✁ Runways access (entries and exits), containing the actual holding points before take-off and exit points after landing with specific speed limitations; ✁ Taxiways intersecting runways, with access re- strictions ;

✁ Simple taxiways, where speed limitations is mod- eled as a function of the turning rate (figure 3).

Connections between taxiways are limited (it is not always possible to proceed from a taxiway to another, even if they are intersecting). The airport description specify usable taxiways connections.

Thus, the airport is defined by a graph: links represent connections between taxiways whereas nodes are taxiways themselves, gate positions, and landing or take-off points. The cost from a taxiway node to its connected nodes is the time spent to proceed via this taxiway, taking into account speed limitations due to this taxiway. The cost from the other nodes (gates and runway positions) to their connected nodes is null.

Figure 1 represents the graphs of Roissy and Orly airports. These graphs are obviously connected. Classic graph algorithms can be used to compute alternative paths for aircraft: ✁ A Recursive Enumeration algorithm [START_REF] Victor | Computing the K Shortest Paths : A New Algorithm and an Experimental Comparison[END_REF] using the Dijkstra's result can compute the ✆ best paths from a given node to another.

✁ A Branch and Bound algorithm [START_REF] Horst | Global Optimization, Deterministic Approaches[END_REF] can com- pute all alternate paths lengthening the best path less than a given distance or time.

Aircraft model

Aircraft are described by their flight-plan (ident, departure or arrival time, gate position, requested runway, eventually their CFMU slot...), their wake turbulence category (low, medium or high) and their take off or landing distance.

In order to perform conflict detection, a model for aircraft separation is defined. This model takes into account runways area, 90 meters away from each side of the runway (or 150 meters away on bad weather conditions). On these area, aircraft are considered on the runway even if they are not taking off or landing.

Aircraft separation model is defined as follows :

✁ aircraft on gate position are separated with all other aircraft.

✁ The distance between two taxiing aircraft must never be lower than 60 meters.

✁ No more than one aircraft at a time can take off or land on a given runway.

✁ A time separation of 1, 2 or 3 minutes (depending on the aircraft category) is necessary after a take off to clear next take off or landing from wake turbulence.

✁ When an aircraft is proceeding for take off or land- ing on a given runway, other aircraft can be taxiing on the same runway area only if they are behind the proceeding one.

Speed uncertainty

Speed uncertainty is modeled as a fixed percentage of the initial defined speed (which is function of procedures and turning rate). Therefore, an aircraft is considered to occupy multiple possible positions at a given time.

Separation is ensured if all of the possible aircraft positions are separated from others, as defined before.

When an aircraft is following an other one, its speed uncertainty will be reduced as the pilot won't go faster than the first one.

Speed uncertainty reduces the validity period of predictions. Thus, simulations with speed uncertainty will be carried out with a lower time window (see 2.5).

Aircraft maneuvers

In order to minimize the total delay and to ensure separations, the path of an aircraft can be modified and aircraft can hold position at the gate, on taxiway or queue at the holding point before take off.

Thus, a ground control order is described by : ✁ The path that the aircraft must follow, choosen among the computed possible paths for the aircraft;

✁ Eventually, the holding position ✝ on this path and the time ✞ until which the aircraft must hold on.

In order to perform acceptable maneuvers, only one holding order should be given to the pilot at a time, and proposed alternative paths should not lead an aircraft to use the same taxiway twice.

With such a holding model (hold at position ✝ un- til time ✞ ) uncertainties defined before can be reduced, while referencing a precise holding position and a precise end holding time (see figure 4).

Simulation model

As the aircraft future positions and movements are not known with a good accuracy, it is necessary to regularly update the situation, every ✟ minutes for example. By the same time, looking a long period ahead is not possible as predictions are not good enough.

Consequently a time window ✠ ☛✡ ✌☞ ✍✟ is defined. Only aircraft taxiing in the time window will be considered. The time window will be shifted every ✟ minutes, the problem reconsidered and a new optimization performed (see figure 5).

At each simulation step (every ✟ minutes), traffic prediction is performed for the next ✠ ☛✡ minutes and pairs of conflicting aircraft positions are extracted. Conflict resolution for this simulation step consist in choosing for each aircraft a path among the given set of possible paths and an optional holding point and time to ensure separations.

Global optimum criteria

In the current version, the global criteria to minimize is defined by the total rolling time (including queueing for runway time), added to the time spent in lengthened trajectory. With this definition, lengthening trajectory appears to be twice more penalizing than holding position.

A*: 1-to-n strategy

In this strategy, aircraft are sorted and considered one after the other.

The optimization problem is reduced to one aircraft: the algorithm must find the best path and/or the best holding point for the aircraft, taking into account the trajectories of the other aircraft already considered. In this point of view, first considered aircraft have priority on last considered aircraft.

Graph modeling

The 1-to-n strategy for an aircraft can be modeled as a graph exploration problem :

✁ A node

✁ Each non conflicting node has two sons :

-The first son is the next position in the same path ✝ ✎ at time ✞ ✒✑ ✔✓ (the aircraft go forward). The cost to reach this son is ✓ .

-The second son is the same position at time ✞ ✕✑ ✖✓ (the aircraft holds position at time ✞ ). The cost to reach this son is ✗ , as a delay is given to the aircraft.

✁ The root nodes are the first position on each path ✝ of the aircraft at current time ✞ ✙✘ . ✁ The terminal nodes are the ones describing a non conflicting position of the aircraft at time ✞ ✙✘ ✚✑ ✛✠ ☛✡ .

An A* algorithm can easily find the best solution for the aircraft.

Sorting method

As last considered aircraft are extremely penalized (they must avoid all first considered aircraft) the way to sort aircraft is a determining factor.

A simple way to assign priority levels is to consider the flight-plan transmission time to the ground controllers.

This option seams the most realistic as ground controllers can hardly take into account an aircraft without its flight-plan. In the simulation context, this is equivalent with sorting aircraft by their departure or arrival time.

However, this option must be refined :

✁ As landing aircraft can't hold position before exit- ing runway, their priority level must be higher than all taking off aircraft.

✁ Queueing for runway aircraft should be sorted in their queue order.

In order to satisfy these principles, a time ✠ ✢✜ is af- fected to each aircraft as a function of its begining time ✠ ✏✘ and its remaining time ✞ ✤✣ : ✠ ✢✜ ☎✥ ✦✠ ✧✘ ★✑ ✩✞ ✤✣ for departures, ✠ ✜ ✥ ✦✠ ✘ ✫✪ ✓ ✭✬ ✯✮ ✭✰ ✏✱ for arrivals. Aircraft are sorted by increasing values of ✠ ✜ .

Genetic Algorithms

In these strategies, classical Genetic Algorithms and Evolutionary Computation principles such as described in the literature [START_REF] Goldberg | Genetic Algorithms in Search, Optimization and Machine Learning[END_REF][START_REF] Michalewicz | Genetic algorithms + Data Structures = Evolution Programs[END_REF] are used. The algorithm is used every ✟ minutes on the problem defined in sec- tion 2.5. Two strategies are developed : in the first one, the algorithm finds a path and an optional holding position for each aircraft. In the secund one, the genetic algorithm finds a path and a priority level for each aircraft, and an A* algorithm is used to compute the resulting trajectories.

Data structure

During each optimization process, each aircraft trajectory is described by its own parameters:

✁ The first strategy needs ✲ numbers (✳ , ✝ , ✞ ) for each aircraft: ✳ is the number of the path, ✝ and ✞ the evantual holding position for the aircraft (if ✝ is reached after ✞ , the aircraft does not stop) as detailed in section 2.4.

✁ The secund strategy needs ✗ numbers (✳ , ✝ ✯✱ ✭✴ ✙✮ ) for each aircraft: ✳ is the number of the path and ✝ ✯✱ ✭✴ ✙✮ its priority level.

Fitness function

For the two strategies, the fitness function must ensure that a solution without any conflict is always better than a solution with a conflict. Consequently it was decided that the fitness of a solution with a conflict should be less than ✵ ✶ and the fitness of a solution without any con- flict more than ✵ ✶ . Thus, for a solution with ✳ ✸✷ remaining conflicts,

✹ ✥ ✓ ✓ ★✑ ✩✳ ✷

For a solution without any conflict,

✹ ✥ ✓ ✗ ✑ ✓ ✗ ✺✑ ✼✻ ✾✽ ✎ ❀✿ ✵ ❁ ✎ ✑ ❃❂ ✎
where ❁ ✎ is the delay of aircraft ✴ and ❂ ✎ the time spent by aircraft ✴ in lengthened trajectory.

Crossover operator

The conflict resolution problem is partially separable as defined in [START_REF] Durand | Genetic crossover operator for partially separable functions[END_REF][START_REF] Durand | Automatic aircraft conflict resolution using genetic algorithms[END_REF]. In order to increase the probability of producing children with a better fitness than their parents, principles applied in [START_REF] Durand | Genetic crossover operator for partially separable functions[END_REF] were applied. For each aircraft ✴ of a population element, a local fitness ✹ ✎ is defined as :

✁ for an aircraft with ✳ ✷ ☞ ❅❄ conflicts, ✹ ✎ ❆✥ ✖✓ ❇❄ ❈❄ ✒❄ ❊❉ ★✳ ✸✷ ;

✁ for a non conflicting aircraft ✹ ✎ ❆✥ ❁ ✎ ❋✑ ❃❂ •✎ .
The crossover operator is presented on the figure 6. First two population elements are randomly chosen. For each parent

Mutation operator

For each candidate to mutation, parameters of an aircraft having one of the worst local fitness are modified. The crossover and mutation operators are quite deterministic at the beginning as many conflicts have to be to solved. They focus on making feasible solutions. When solutions without conflict come in the population, they become less deterministic.

Clusters

In order to lower the complexity of the problem as often as possible, a transitive closure is applied on conflicting aircraft pairs and gives the different clusters of conflicting aircraft [START_REF] Durand | Automatic aircraft conflict resolution using genetic algorithms[END_REF]. The different clusters will be solved independently at first. If the resolution of two clusters creates new conflicting positions between them, the two clusters are unified and the resultant cluster is solved.

Sharing

The problem is very combinatorial and may have many local optima. In order to prevent the algorithm from a premature convergence, the sharing process introduced by Yin and Germay [START_REF] Yin | A fast genetic algorithm with sharing scheme using cluster analysis methods in multimodal function optimization[END_REF] is used. The complexity of this sharing process has the great advantage to be in ✳ ◆▼ ❖ ❈◗ ✏| ❙✳ ✸❚ (instead of ✳ ✶ for classical sharing) if ✳ is the size of the population.

A distance between two chromosomes must be defined to implement a sharing process. Defining a distance between two sets of ❯ trajectories is not very simple. In the experiments, the following distance is introduced:

❱ | ✂ ❳❲ ❍ (❚ )✥ ✻ ✽ ✎ ✿ ✵ ❬ ❂ •❭ ☛( ✙❂ •) ☛( ❬ ❯ ❂ ❭ ( (resp ❂ ) ( ) is the ✴ ❵❴ ❜❛ aircraft path length of chromo- some ✂
(resp ❍ ). As the paths are sorted according to their length, the distance increases with the difference of lengths.

Ending criteria

As time to solve a problem is limited, the number of generations is limited: as long as no available solution is found, the number of generation is limited to ❝ ❈❄ . The algorithm is stopped ✗ ❞❄ generations after the first acceptable solution (with no remaining conflict) is found.

Experimental results

Simulations

Simulations are carried out with real flight plans of Roissy Charles De Gaulle and Orly airports on a complete day (May ✓ ❇❡ ❈❢ ❞❴ ❜❛ 1999).

Three strategies are compared :

✁ in the "1-to-n method", aircraft are sorted as de- scribed in 3.2. They keep the same priority level during all the simulation and an A* algorithm finds the best solution.

✁ in the "Global method", a genetic algorithm finds a path and an optional holding position for each aircraft in order to minimize the global criteria described in 2.6.

✁ in the "Mixed global method", a genetic algo- rithm finds a path and a priority level for each aircraft and the fitness function is computed by an A* algorithm applied on sorted aircraft.

Simulations parameters:

Time window : ✠ ✡ ✥ ❣❝ ✐❤ •✳

Comparing the strategies

As Roissy and Orly simulation results has given the same relative conclusions about the 3 strategies efficiency, figures related in this article only concerns Roissy Charles De Gaulle airport. Figure 7 gives the mean delay as a function of the number of aircraft moving on the taxiways for the different methods. When number of aircraft increases, the mixed method appears to be the best one.

Figure 8 gives for the different strategies the number of aircraft simultaneously moving as a function of time. It appears that the mixed method keeps a lower number of moving aircraft during heavy time periods: a good resolution of ground traffic conflicts allows to decrease delay and then leads better situations with less moving aircraft.

Genetic algorithm efficiency

In order to observe the GA efficiency, figure 9 gives the number of generations required by the GA as a function of time for the two GA strategies.

For the global method, the different peaks which appears at 7, 10 am, 1 and 7 pm are the traffic peaks.

For the mixed method, the global optimum is always found with a few number of generations : sorting aircraft by evolutive priority levels seams to be very efficient as far as ground conflicts resolution is concerned.

Conclusion and further work

A preliminary work has shown that it was possible to build a taxiway adviser in order to optimize the ground traffic on busy airports such as Roissy Charles de Gaulle and Orly. It can be noticed that the modeling was easily improved with new runways on Roissy Charles De Gaulle, different speeds, uncertainties on speeds etc... without changing the algorithm itself. Genetic Algorithms are very efficient on the problem as they search the global optimum of the problem whereas a deterministic algorithm such as a 1-to-n strategy causes more delay.

Further work will concentrate in improving the global criteria for Genetic Algorithms, taking into account for example take off sequencing needs of approach sectors or priority levels for slotted departures.

Figure 1 : 1000 Figure 2 :

 110002 Figure 1: Roissy airport graph -Example of shortest and alternate paths

Figure 3 :

 3 Figure 3: Speed limitation as a function of turning rate ✁ An ✂ ☎✄ algorithm [Pea84] can compute the best path and the corresponding minimal time spent between two given nodes (gate and runway entry for example).

✁

  A Dijkstra algorithm[START_REF] Ravindra | Network Flows, Theory, Algorithms and Applications[END_REF] can compute best paths and corresponding minimal time spent from a given node to every other node.

  Figure 4: Uncertainty reduction on holding points

  of the graph is a position in a path ✝ ✏✎ of the aircraft at time ✞ . ✁ An heuristic function for this node is the minimal remaining time to reach the end of the path. ✁ If a node represents a conflicting position with al- ready considered aircraft, it has no son.

✂

  and ❍ , fitness ✂ ✎ and ❍ ✎ of aircraft ✴ are compared. If ✂ ✎ ✚■ ❍ ✎ , the children will take aircraft ✴ of parent ✂ . If ❍ ✎ ❏■ ✂ ✎ , the children will take aircraft ✴ of parent ❍ . If ✂ ✎ ▲❑ ❍ ✎ children randomly choose aircraft ✂ ✎ or ❍ ✎ or even a combination of ✂ ✎ and ❍ ✎ .

  Figure 6: Crossover operator

Figure 7 :Figure 8 :

 78 Figure 7: Mean delay as a function of the number of moving aircraft.

Figure 9 :

 9 Figure 9: Number of generations as a function of time

SIMulation MODel (FAA)

Total Airspace and Airport Modeler(Preston Group)

Departure Planer

Taxi and Ramp Management And Control (DLR) 1

Biography

Jean-Marc Alliot graduated from the Ecole Polytechnique de Paris and from the Ecole Nationale de l'Aviation Civile (ENAC). He also holds a Ph.D. in Computer Science (1992). He is currently in charge of the global optimization laboratory of CENA and ENAC in Toulouse.

Nicolas Durand graduated from the Ecole Polytechnique de Paris and from the Ecole Nationale de l'Aviation Civile (ENAC). He has been a design engineer at the Centre d'Etudes de la Navigation Aérienne (CENA) since 1992 and holds a Ph.D. in computer Science (1996).

Erwan Page graduated from the Ecole Nationale de l'Aviation Civile (ENAC). He has been in the design of Flow Management and Decision Support tools for French ATC for 15 years and is head of the "Aéroport, Tour et zone Terminale" (ATT) division at CENA/Paris.

Jean-Baptiste Gotteland graduated from the Ecole Nationale de l'Aviation Civile (ENAC). He is completing a Ph.D. in Computer Science at the Institut National Polytechnique de Toulouse.


T Miquel

F Mora Camino Version

A NONLINEAR CONTROL APPROACH FOR RELATIVE GUIDANCE OF COMMERCI AL AIRCRAFT 7 JULY 2003 CENA/SAS/NT03-768

SUMMARY

The delegation to the flight crew of tasks currently performed by the air traffic controllers provides a new perspective to significantly increase air traffic control capacity. Then, merging and station keeping operations have appeared to be of interest even if they rise new technical challenges. Indeed, these unusual maneuvers in the field of Civil Aviation seem difficult to be performed manually and since in this case they will result in an increase of the flight crew workload, new on-board automated functions should be developed. This paper investigates the design of an auto pilot mode dedicated to merging and station keeping maneuvers behind a leading aircraft. The proposed relative flight guidance law is based on a non linear sliding mode contr ol technique. Simulation results are displayed and discussed.

KEY WORDS

Relative Guidance, Flight Control Laws, Non Linear Control.

INTRODUCTION

The increase of air traffic urges the design of new devices to improve significantly air traffic control capacity [START_REF]CARA/ASAS Activity 5 description of a first package of GS/AS applications[END_REF] . So, to meet this challenge, new concepts such as "relative guidance" are presently of interest to contribute to the delegation of tasks to the crew [START_REF] Agelii | Flight deck simulations of station keeping[END_REF]. The relative guidance function should provide an aircraft the ability to establish and to ma intain a given time or distance behind a designated aircraft (see figure 1). Nowadays no automatic control mode is available on-board civil aircraft to perform this task and few studies have been published in this field [START_REF] Ivanescu | Impact of ADS-B link characteristics on the performances of in -trail following aircraft[END_REF][START_REF] Giulietti | Au tonomous Formation Flight[END_REF]. This communication focuses on the design of a nonlinear controller which enables an aircraft to perform merging maneuvers and to maintain station keeping behind a designated aircraft.

The inves tigated approach is based on the inversion of a simplified non linear model of the relative dynamics coupled with a sliding mode control technique. Preliminary evaluations are performed using simulation models of two wide body aircraft.

STATE SPACE REPRESENTATION OF RELATIVE DYNAMICS

The range and the relative bearing of the two aircraft can be derived from the along track distance TK and the cross track distance XTK by:

2 2 XTK TK + = ρ and ) / ( TK XTK arctg + = ψ µ (1) while )) ) cos( /( ) ) sin( (( y x W V W V arctg + + = ψ ψ χ (2)
These relations translate the variables manipulated by air traffic controllers and pilots (i.e.

"merge 5 NM behind, 1 NM left") into variables processed by the relative guidance device.

The state vector retained here to describe the relative dynamics of the two aircraft is chosen as follows:

[ ] ψ µ ρ , , , V T X = (3)
This leads to the following non linear affine state space representation of the relative guidance kinematics, where u denotes the control vector : CENA/SAS/NT03-768/T. MIQUEL, F. MORA CAMINO

VERSION 1.0 4 u X g X f X ) ( ) ( + = & (4) with [ ] ϕ , c T V u = (5)
where

            - - - - - - - = 0 / / ) sin( ) sin( ) cos( ) cos( ) ( V L L L L V V V V V X f τ ρ µ ψ µ ψ µ ψ µ ψ and             = V g X g V / 0 0 0 0 / 1 0 0 ) ( τ (6)
The output vector, representative of the relative position of the aircraft, is chosen as :

[ ] µ ρ , T y = (7)
The control objective retained in this study is to lead the track variables to the respective desired values TK d and XTK d while the track angle of the trailing aircraft χ becomes equal to the one of the leader. Through relations [START_REF]CARA/ASAS Activity 5 description of a first package of GS/AS applications[END_REF], this objective can be fully expressed with variables ρ, µ and χ. The resulting transient and permanent regimes should be compatible with the handling qualities of the aircraft and the required degree of comfort for a passenger aircraft.

THE RELATIVE GUIDANCE CONTROL LAW

Following non linear control theory [START_REF] Isidori | Nonlinear control systems[END_REF], the relative degree of theses outputs is obtained through their derivation with respect to time until terms of the input vector appears. Here we get:

( ) ( ) u X X y ⋅ ∆ + ∆ =       = 0 µ ρ (8) with ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )                   - - - - ( ) ( ) ( ) ( ) X X X l 0 1 ∆ ∆ - = - (13)
At this point, to fully determine the control law, a transient dynamics must be specified for the output vector y. In another work [START_REF] Miquel | Relative guidance of aircraft: analysis and comparison of two control approaches[END_REF], the authors have evaluated the case in which the objective of the non linear inverse control law is to get linear decoupled output dynamics.

This former approach has shown some strong limitations: high solicitation of the actuators leading to their saturation, generation of complex relative maneuvers and poor robustness with respect to modeling errors. Then the current study investigates the use of a modified sliding mode controller which leads the outputs to adopt progressively first order linear dynamics. Indeed, sliding mode control has been recognized to be a technique which can provide robustness in the control of non linear systems [START_REF] Mora -Camino | Zero gravity atmospheric flight by robust nonlinear inverse dynamics[END_REF].

For the derivation of a sliding mode controller, let us select two switching surfaces as follows: where k1 and k2 are positive numbers, and ρd and µd are the desired range and relative bearing derived from the desired along and cross track distances. Observe that their derivatives take the form:

( ) ( )    - + = - + = d d k s k s µ µ µ ρ ρ ρ
      + =       µ ρ 2 1 2 1 k k y s s (16)
Then, a classical Lyapunov function [START_REF] Slotine | Sliding controller design for nonlinear systems[END_REF] can be introduced:

2 2 2 1 2 1 2 1 s s W + = (17)
It is clear that if dW/dt is forced to remain negative, ρ and µ will be inclined to follow the above (relations ( 14)) first order sliding mode dynamics towards ρd and µd . For that, let us choose a control law such that the derivatives of the switching surfaces satisfy:

( ) ( )    - = - = 2 22 2 1 11 1 s sat a s s sat a s & & (18)
where a ii are positive parameters and the saturation functions are given by :

( )      > ∀ ≤ ∀ = max max max ) ( i i i i i i i i s s s s sgn s s s s sat (19)
Then the relative guidance control law adopts the following expression:

( ) ( ) ( ) ( ) X l k s sat a k s sat a X u +       - ⋅ - - ⋅ - - ∆ = µ ρ & & 2 2 22 1 1 11 1 (20)

CASE STUDY

Different scenarios have been des igned in order to evaluate the performances of the proposed control law. One of them is displayed here. This scenario considers an approach and includes the following features over a time period of 700 seconds for two wide body aircraft:

The leading aircraft starts at x 0 = 0 NM, y 0 = 0 NM, with an initial conventional airspeed of 200kts and a heading of 90 degrees. The bank angle of the leading aircraft remains at zero, except between 220 sec and 320 sec, when the leader changes its heading by about 175 de grees (near a half circle) through setting at 20 the controlled bank angle . Furthermore, the controlled conventional airspeed of the leading aircraft is set to 160 kts for t > 600 sec, which

Figure - 1

 1 Figure -1 In-bound traffic at Roissy Charles -De-Gaulle

Figure- 2

 2 Figure-2 Relative guidance parameters

  deceleration. The trailing aircraft starts at x0 = -8 NM, y0 = +4 NM, with initial conventional airspeed and heading of 200 kts and 90 degrees respectively. The requested separation for the trailing aircraft is 90 seconds behind the leading aircraft. The whole maneuver is supposed to take place within a wind of 20 kts blowing from North. During the maneuver, the inputs have been limited to ± 20° for ϕ with a maximum rate of ±5°/sec and to values between 140 kts and 250 kts for V c with a maximum rate of ± 1 kts/sec. The time constant τV has been chosen equal to 40 seconds. The values of the parameters of the control law have been taken as follows: the leading and the trailing aircraft in the horizontal plane as well as the normal load factor are represented hereafter :

Figure- 3 :Figure- 4 : 8 Figure- 5 :Figure- 6 :

 34856 Figure-3 : Aircraft movement in the horizontal plane (axes in NM)Figure-4 : Load factor(minus 1) versus time (in sec)The evolutions with time of current delay and range between aircraft are shown in the

Figure- 7 :Figure- 8 :

 78 Figure-7 : Controlled and actual airspeed (in kts) versus time (in sec) Figure-8 : Controlled and actual bank angle (in degrees) versus time (in sec)

& & & & & &

---------------------

& & & & &

Note that in the above relations, the time rates of ψL and VL have been neglected. Here the relative degrees of the outputs are both equal to two and the whole dynamics is covered by their dynamics. Then, there is no ground left for inner dynamics at this level of detail for the considered model of the relative flight dynamics. Moreover, as the rank of the matrix ∆(x) is in normal conditions (ρ > ρ min > 0) equal to 2 , this last relation can be inverted to get a non linear control law such as:

ACKNOWLEGMENTS

The authors wish to thank Mr Francis Casaux from CENA/EUROCONTROL-Brétigny and Professor Jean Lévine from Ecole Nationale Supérieure des Mines de Paris for helpful discussions during this research.