2.1 Flight-deck Interval Management (FIM)

Continuous research effort, especially with regard to ASTAR, has led to the development of the Minimal Operating Performance Standards for FIM (MOPS-FIM)⁽⁴⁾, which includes all the important definitions and equations for FIM operation.

The key element for the time-based separation is the spacing error term e(t).

To calculate e(t), the remaining flight time (‘time-to-go’, TTG) of the own aircraft (‘Ownship’, OWN) is compared with that of a pre-assigned target aircraft (‘Traffic-To-Follow’, TTF) to a common waypoint called the Achieve-By-Point (ABP). By adding a nominal spacing target (Δ), e(t) can be calculated as follows:

(1) \begin{equation}e(t) = TT{G_{{\rm{OWN}}}}(t) - \left( {TT{G_{{\rm{TTF}}}}(t) + \Delta} \right)\end{equation}

By adding the current time, TTG can be replaced with the estimated time of arrival (ETA) for each aircraft such that the term within the parentheses (ETA _TTF + △) can be summarised as the required time of arrival (RTA):

(2) \begin{equation}e(t) = ET{A_{OWN}}(t) - RTA(t)\end{equation}

Consequently, a positive error term indicates a delayed arrival whereas a negative one indicates an early arrival. The FIM logic then tries to eliminate any spacing error via speed commands.

Technically, the TTF’s TTG is calculated from its current speed and position obtained via ADS-B, while the ownship’s TTG is obtained from data supplied by the flight management system.

2.2 Interval Management – Speed Planning (IM-SP)

IM-SP is a newly developed logic for FIM, designed to evaluate and re-plan the entire remaining speed schedule to achieve the spacing objective. Modifications are made primarily to existing speed changes (as opposed to the current speed in ASTAR) by changing their target Calibrated Airspeed (CAS) to a higher or lower value, or the location where the speed change is initiated to an earlier or later position (Fig. 1). Additional speed changes are inserted only when deemed favourable or necessary.

Figure 1. Speed change modification principle of IM-SP.^{(Reference Riedel, Takahashi and Itoh12)} Here, AP⁶ is modified, thus affecting AP⁵ for DTG change, or both AP⁵ and AP⁴ for a CAS_TGT change.

To choose the most appropriate speed plan modification, IM-SP uses a two-stage selection process, involving a time-based primary stage and a cost-function-driven secondary stage, as shown in Fig. 2.

Figure 2. Functional selection stages and flow of IM-SP.^{(Reference Riedel, Takahashi and Itoh12)}

The cost for a speed plan modification is based on the characteristics of the modification itself (as further explained in Section 3.3).

Internally, these characteristics are associated with Action Points (APs). APs are important points in the speed schedule. The most prominent APs are the beginning or end of speed changes. Other important APs include the Mach transition point as well as the system initiation and termination points. They are numbered sequentially and store information such as the location (‘Distance-To-Go’, DTG), nominal CAS, target CAS (CAS_TGT) and other scoring-related attributes. Modifications to an AP can affect only subsequent APs, as shown in Fig. 1. The modified AP, subsequently affected APs and unaffected APs are indicated in orange, black and grey, respectively.

2.3 Particle Swarm Optimisation (PSO)

Particle swarm optimisation (PSO) is a Nature-inspired metaheuristic that was first proposed by Kennedy and Eberhard in 1995^{(Reference Kennedy and Eberhart15)}. It is classified as an Evolutionary Algorithm (EA) that is based on swarm or social behaviour, such as the flocking behaviour of birds or schooling behaviour of fish.

PSO works on the principle that the particles in the population move individually through the search space to find the global minimum, where they are driven by both an individual force (‘pursue your own goal’) and a social force (‘follow the leader’).

First, all the particles are distributed over the search space, giving each particle a starting position ( ${\vec p_0}$ ) and an initial velocity ( ${\vec v_0})$ . The particle with the best position, i.e. the current global minimum ( ${\vec g_{Best}}$ ), becomes the leader. Further, each particle memorises its individual best position ( ${\vec p_{Best}}$ ), which is not shared with the other particles.

All the particles now start to move, and their new velocity ( ${\vec v_{n+1}}$ ) is determined by their inertia ( $\omega$ ) and their individual and social pull, given by Ref. (Reference Shi and Eberhart16) as

(3) \begin{equation}{\vec v_{n+1}} = {\rm{}}\omega \cdot {\vec v_n} + {\rm{}}{c_p} \cdot {r_1} \cdot \left( {{{\vec p}_{Best}} - {{\vec p}_n}} \right) + {c_g} \cdot {r_2} \cdot \left( {{{\vec g}_{Best}} - {{\vec p}_n}} \right),\end{equation}

where c represents the user-defined magnification factors and r random values, the latter of which are newly determined in each iteration.

The particles’ next position ( ${\vec p_{n+1}}$ ) is then given by Ref. (Reference Shi and Eberhart16) as

(4) \begin{equation}{\vec p_{n+1}} = {\rm{}}{\vec p_n} + {\rm{}}{\vec v_{n+1}}\end{equation}

If a particle finds a position better than ${\vec p_{Best}}$ , the new value is memorised in lieu of the former. Further, if a new global minimum is found, ${\vec g_{Best}}$ is updated and the corresponding particle becomes the leader.

PSO ends after a fixed number of iterations or when an abort criterion is met.

2.4 Speed-constrained Multi-objective Particle Swarm Optimisation (SMPSO)

SMPSO is an enhanced version of PSO, precisely OMOPSO^{(Reference Sierra and Coello Coello17)}, a multi-objective PSO algorithm that incorporates additional features (further described below) that have been shown to facilitate a more uniform examination of the search area as well as faster exploration of the Pareto front^{(Reference Nebro, Durillo, GarcÍa-Nieto, Coello Coello, Luna and Alba13)}.

Previous applications of SMPSO in the aerospace context include the optimisation of air routes for hub connections^{(Reference Rahil, Abou El Majd and Bouchoum14)}.

As suggested by its name, SMPSO’s main feature is the speed constraining term that limits an individual particle’s velocity to avoid ‘swarm explosion’, expressing the phenomenon that particles gain high speed, skip unexplored areas and ultimately reach the search space limits.

Based on an idea originally proposed by Clerc and Kennedy^{(Reference Clerc and Kennedy18)}, SMPSO uses the sum ( $\varphi$ ) of the user-defined magnification factors c,

(5) \begin{equation}\varphi = {c_p} + {c_g}\end{equation}

to define a velocity constriction coefficient (χ), which is expressed in a shortened form as

(6) \begin{align}{{\chi }} = \left\{ \begin{array}{ll} \dfrac{2}{2 - \varphi - \sqrt {{\varphi ^2} - 4\varphi } }, &\quad \varphi < 4 \\ 1,&\quad \varphi \le 4 \end{array} \right. .\end{align}

Note that, for SMPSO, c _p and c _g are randomly chosen, and they are limited to the range of [1.5, 2.5].

The new velocity of a particle is then determined by multiplying the velocity of the PSO algorithm (Equation 3) by the constriction coefficient (Equation 6).

Finally, the velocity of the particles in each dimension is limited by half the dimension span:

(7) \begin{equation}{\delta _j} = {{\left( {upper\_limi{t_j} - lower\_limi{t_j}} \right)} \over 2}\end{equation}

Thus, the dimensional speed (v _j) lies within the range [ $-\delta_{j},\ \delta_{j}$ ].

Another important feature of SMPSO is the mutation concept, which is inherited from OMOPSO. Mutation describes the (random) partial change of the particle’s position component values, occurring at a given probability. The added movement improves propagation and lowers the risk of particles becoming trapped in local minima.

3.1 Objective function parameters

A total of four parameters, i.e. one performance, one ecology and two usability factors, are chosen for the objective function. The input argument is IM-SP’s cost function weight vector ( $\vec q$ ), giving the general equation

(8) \begin{equation}\min f({\vec q}) = {\rm{min}}\left( {{f_1}({\vec q}),{f_2}({\vec q}),{f_3}({\vec q}),{f_4}({\vec q})} \right)\end{equation}

For minimisation, some factors have been adapted as described below. The block diagram of IM-SP is shown in Fig. 3.

Figure 3. Block diagram of IM-SP relevant blocks. IM-SP’s logic block is highlighted in blue in the centre. To the left is the ownship’s trajectory as the main input provider; To the right is the flight dynamics block as a post-processor of the speed plan to obtain further performance data. Information flow from left to right: 1) The speed envelope, the ownship’s TTG (TTG_OWN) and the speed plan are output from the ownship’s trajectory. 2) TTG_OWN is compared with the spacing target time (given by the nominal spacing and traffic’s TTG) to obtain the spacing error e(t). 3) From these inputs, IM-SP generates a new speed plan and logs a change history from which the number of total commands and critical commands are obtained. 4) The speed plan is looped back to the ownship trajectory and passed forward to the flight dynamics. 5) From the flight dynamics the new fuel burn and TTG_OWN are calculated. 6) TTG_OWN is looped back to the ownship trajectory and returned as the spacing error correcting variable

3.1.1 Final spacing error

The final spacing error indicates the deviation between the Actual Time of Arrival (ATA) and the RTA at the ABP, or the measuring gate; thus, it is the primary indicator of IM performance. Ideally, this value is 0s, i.e., neither too early nor too late; thus, the partial objective function is

(9) \begin{equation}{f_1}({\vec q}) = \left| {e(t)} \right|\end{equation}

In a realistic environment, FIM operation is deemed successful if 95% of all arrivals are within ±10s of the RTA^{(4,Reference De Gelder, Bussink, Knapen and in ‘t Veld19)} . For easier comparison with Ref. (Reference Riedel, Takahashi and Itoh12), ±5s is used as the target range here.

3.1.2 Fuel burn

The fuel burn is used as an ecological indicator of FIM operation. In general, the lower the fuel burn, the lower the cost. Therefore,

(10) \begin{equation}{f_2}({\vec q}) = Fuel\ Burn.\end{equation}

For the reference profile, a fuel burn of 1600.4kg was estimated. However, as FIM operation requires deviation from the reference profile, a slight increase in fuel burn must be expected.

3.1.3 Total number of commands

The number of commands represents how often the speed is changed; thus, in a federated system environment, it represents the required number of interactions by the crew. Therefore, it is the first indication of the workload imposed by the system on the crew. In general, fewer commands are favourable; however, too few commands could result in undesirably high speed change magnitudes. Therefore, ideally, the number of speed commands (n) will be equal to the number of changes in the reference profile (n _Ref):

(11) \begin{equation}{f_3}({\vec q}) = \left| {n - {n_{Ref}}} \right|\end{equation}

Here, the reference profile has six commands.

3.1.4 Number of critical commands

As MOPS-FIM expects a command to take effect within 11s of its annunciation (including latency and engine delay)⁽⁴⁾, the crew is required to immediately recognise and comply with the command. Otherwise, there is a risk that the FIM logic’s intentions might not be met. To reduce this risk, IM-SP was designed to avoid commands on short notice, unless absolutely necessary (Subsection 2.2.). These commands are further referred to as ‘critical commands’ (n _Crit) and include all commands with lead times shorter than 15s. Ideally, no critical command is required during IM-SP-based FIM operation. Thus,

(12) \begin{equation}{f_4}({\vec q}) = {n_{{\rm{Crit}}}}\end{equation}

3.2 Other non-objective function parameters

In addition to the objective function parameters, four other parameters, two indicating usability and two that are related to the speed profile, are evaluated in this paper to show the side effects of the cost function setting as well as for clearer comparison with previous studies.

3.3 Cost function elements and weights

As shown in Fig. 2, IM-SP includes a cost-function-based selection stage involving the scoring vector ( $\vec s$ ). This vector consists of five elements: two primary selection attributes (s _AEM and s _TTG) and three secondary penalty attributes (s _TTR, s _APD and s _Type).

For each scoring element, a corresponding weight (q _AEM, q _TTG etc.) exists in the weight vector ( $\vec q$ ). The final score (S) for each solution (i) is then obtained as the scalar product of the two vectors:

(13) \begin{equation}{S_i} = \overrightarrow {{s_i}} \cdot \vec q\end{equation}

Each individual weight can assume a value between 0 and 1; however, the primary attributes are dependent on each other, giving the constraint q _AEM + q _TTG = 1. Therefore, if q _AEM is 1, q _TTG will be 0 and vice versa. All other weights are independent of each other; thus, the total number of inputs can be reduced to four. Therefore, the optimisation problem involves four input variables and four output objectives.

Table 1 Weight factor settings by profile

The following subsections outline each element’s function and motivation. A more detailed explanation including the corresponding equations can be found in Ref. (Reference Riedel, Takahashi and Itoh12). The function weights used in the initial study to benchmark the optimised weights are given in Table 1.

3.3.4 Action point distance (s _APD, q _APD)

The action point distance is introduced to allow for a sufficient distance between APs to avoid a high frequency of speed commands. In particular, this function penalises speed changes that commence less than 5NM apart from each other.

Figure 4. Horizontal and vertical routing of the ownship^{(Reference Riedel, Takahashi and Itoh12)}.

4.1 Aircraft and calculation model

All the performance calculations in this study were based on EUROCONTROL’s BADA model, version 3.12⁽²⁰⁾ and MOPS-FIM⁽³⁾ with standard atmospheric conditions. The aircraft model used for the ownship and traffic was a Boeing 787-8 at reference mass.

4.2 Flight path

The horizontal and vertical flight paths are shown in Fig. 4. The flight originates at waypoint SMOLT and continues to the KAIHO arrival toward the ILS X approach for Runway 34 L at Tokyo International Airport (RJTT). Vertically, the profile is initiated at 38,000ft (FL380) and laid out as a continuous descent approach (CDA) with a fixed geometric flight path angle of −2.2°^{(Reference Itoh, Fukushima, Hirabayashi, Wickramasinghe and Toratani21,Reference Itoh, Wickramasinghe, Hirabayashi and Fukushima22)} . The angle is kept from the top of descent until glide slope capture, from where the descent is continued with −3.0°. Speed constraints were applied below 10,000ft (250Kn), at KAIHO (180Kn) and at the Final Approach Fix (FAF, 150Kn).

4.3 Scenarios and data sets

The original study on IM-SP consists of 147 simulations based on 7 fixed offset scenarios ranging from −30s to +30s in 10-s intervals as well as 5 error patterns with 28 different settings, i.e. 2 amplifications (high or low) and 2 directions (positive or negative) for each of the 7 offsets.

In this study, one new error pattern was taken from the SPICA simulator^{(Reference Riedel, Takahashi, Tatsukawa and Itoh23)}, increasing the total number of simulations to 175.

The new pattern is used for various examples throughout this paper (Sections 5.1 and 6.1), and it was included because of its demanding characteristics, introducing a nearly linear progressing error of 10s within the last 30NM.

Besides the full data set (175 scenarios), a reduced optimisation set was added for computation time considerations. This set contains only 18 scenarios, 6 fixed offsets and 6 error patterns in both directions, but only with high amplification and without any initial offset.

4.4 SMPSO

The SMPSO source code used in this study is based on the Multi-Objective Evolutionary Algorithm (MOEA) Framework 2.12 by Dave Hadka^{(Reference Hadka24)}. In this implementation of SMPSO, the values of c _p and c _g are randomised for every particle in every iteration and the probability of a mutation is 1 to 6. In total, 81 particles were simulated for 50 steps each, giving a total 4050 calculations.

5.1. Arrival expedition margin versus time-to-go

The increase in q _AEM, and conversely the decrease in q _TTG, has a significant effect on all the attributes. In particular, an increasing trend is observed in the final spacing error, number of commands, number of accelerations and RPD. By contrast, a decreasing trend is observed in the fuel burn, speed brake use time and deviation from nominal CAS.

Remarkably, q _AEM values lower than 0.25 show fewer commands than the reference profile and have only a minor impact on the fuel burn and speed brake use. However, beyond 0.25, a significant change is observed in these values, especially RPD.

Considering the development of the final spacing error and commands, it can be easily assumed that (if fuel burn was of no concern), a q _AEM of 0, i.e. a q _TTG of 1, would render the best results for the optimisation problem. However, further investigation of the generated profiles indicates an operationally undesirable phenomenon called ‘backloading’, as shown in Fig. 9.

Figure 5. Overall results dependent on q _AEM (x-axis), q _TTG given by 1-q _AEM.

Figure 6. Overall results dependent on q _TTR.

Figure 7. Overall results dependent on q _APD.

Figure 8. Overall results dependent on q _Type.

Figure 9. Example of the backloading phenomenon. Upper graph: The ownship’s speed profile within the speed envelope Lower graph: The external spacing error originating from the TTF’s trajectory The x-axis shows the DTG of the ownship, aligned for both graphs, the graph’s left border (0NM) indicates the runway threshold or final position; the right border (130NM) the initial position, at which FIM was commenced. Chronologically, the graph is read right-to-left: 1) At 130NM (initial position), a spacing error of +15s is indicated, meaning that the ownship needs to expedite its arrival by 15s, compared with the originally estimated time. 2) From 130NM until 80NM the error gradually increases to +30s. 3) From 80NM to approximately 30NM the error mostly remains at +30s. 4) From 30NM to 0NM (final position) the error reduces to +20s. This shows that, in hindsight, the ownship would have needed to only shorten its remaining flight time by 20s (compared with the initial estimation) to arrive with no spacing error. Further, it is seen that the TTF slowed down during the last 30NM; Thus, if the ownship had already compensated for +30s at this time, it would now be required to delay its arrival by 10s.

Owing to the characteristics of s _TTG, a setting of q _AEM = 0 (q _TTG = 1) biases the modifications to be summarised and applied to the speed changes during the later flight phases. In the case of an initially high and further developing spacing error (the arrival must be further expedited), the system will try to assume the maximum speed profile, or if it is not sufficient to compensate the error, it will add a single late acceleration that is consequently shorter but has a higher CAS value.

In the example given below, this would cause an initial acceleration from 310Kn to 340Kn at a DTG of approximately 85NM, followed by another acceleration to 347Kn, which is kept for a short distance of 5NM.

The above-mentioned observation was made for q _AEM values below 0.3, and it coincides with the sudden rise in RPD as well as the reduction in the fuel burn and speed brake use time, as shown in Fig. 5.

Another problem that arises with such behaviour is that the profile may be saturated, i.e. it may become the maximum speed profile, once the initial deceleration is commenced, and it might not be able to deal with further increasing error. A detailed discussion on this issue is provided in Section 6.1.

Consequently, solutions that incorporate a q _AEM of less than 0.3 are disregarded.

5.2 Time-to-react

Figure 6 shows that the TTR penalty portion can suppress commands with short lead times (here, from 2.8 to 1) and reduce the overall commands and accelerations, as intended by its design. Furthermore, a minor reduction in the final spacing error, which is minimised for a q _TTR of 0.33, can be observed. Beyond 0.5, all of the above-mentioned attributes stabilise, while none of the others are affected significantly.

5.3 Action point distance

The primary purpose of the APD penalty is to avoid many commands over short distances. While this is not reflected by any of the attributes directly, in Fig. 7 for q _APD values starting from 0.25, a sudden increase can be observed in RPD and the maximum nominal CAS deviation, accompanied by a reduction in fuel burn and speed brake use, which supports the bias of the penalty toward fewer but higher-magnitude speed changes, compared with many but small speed changes.

Other noteworthy points include the maximum final spacing error (0.34) and the (local) maxima for total and critical commands for a q _APD of 0. Further, a marginal increase is observed in the total commands for higher values of q _APD.

5.4 Action point type

The type penalty was introduced to prevent additional speed commands if the same error reduction can be realised by modifications to existing speed changes. As shown in Fig. 8, this is successfully achieved, reducing the total number of commands from 8 to 6 over the range of q _Type, with only a minor impact on the final spacing error. Furthermore, it can be seen that the penalty effectively reduces the RPD at q _Type values of 0.22 to 0.56. However, a distinct increase in fuel burn is observed for q _Type values greater than 0.25.

5.5 Single parameter analysis summary

Two observations can be made from all the graphs.

First, the penalty portions (s _TTR, s _APD and s _Type) perform as intended, i.e. they have only a minor effect on the final spacing error compared with s _AEM and s _TTG, and they augment the final speed profile for improved characteristics (fewer overall and critical commands).

Second, a reduction in fuel burn is mainly at the expense of more (critical) commands or a higher final spacing error, and vice versa. In other words, the objectives are in conflict with each other and thus require a trade-off solution.

Nevertheless, as the above-mentioned analysis involves only one parameter at a time, potential correlations are not reflected. For example, naturally, a larger number of total commands increases the probability of one of them being a critical command; however, the inverse assumption (many critical commands imply a large number of total commands) cannot be made. Therefore, this analysis should serve only as an indicator of the quality of each weight’s influence and the result to be expected from the multi-objective optimisation described in the next section.

6.1 Speed profile characteristics

Figures 13 and 14 show two examples for comparing the resulting speed profiles generated by T_Opt and F_Opt. Figure 13 shows an example of an initially positive and increasing spacing error, identical to the one shown in Fig. 9.

Figure 13. Profile comparison for both settings with an initial positive error.

Figure 14. Profile comparison for both settings with an initial negative error.

The profile generated by T_Opt uses fewer commands and commences its initial acceleration at a later time (DTG, 94NM), albeit with higher speed magnitudes.

By contrast, the profile generated by F_Opt reacts to the initial spacing error with an earlier acceleration (119NM), resulting in more but smaller speed steps and giving a lower top speed that contributes to the lower fuel burn. In the latter half (<60NM), the speed profiles are nearly identical, with a small difference for the final or next-to-final deceleration.

In the second example (Fig. 14), the same error profile is used, but in the opposite direction with a corresponding negative initial spacing error. Both settings react with the same initial deceleration. However, F_Opt adjusts to the further development of the spacing error by expediting the deceleration from 250Kn to 220Kn, while the profile generated by T_Opt remains on the nominal profile and reduces the next deceleration from 220Kn to 180Kn to a lower target CAS.

In summary, it can be seen that F_Opt, using a higher q_AEM, reacts earlier to changes in the spacing error, while T_Opt reacts later.

Depending on the progression of the error, this makes T_Opt more susceptible to errors progressing in the same direction, e.g. a positive error that increases, than F_Opt. This is also reflected by the final spacing error distribution of the full data set shown in Section 6.2.

By contrast, for an error that changes in direction or cancels itself, e.g. a positive error that decreases (or vice versa, as shown in Figs. 13 and 14), T_Opt reacts later and has more time to combine or revert previous changes, thus resulting in fewer overall commands.

However, in rare cases, i.e. if all other options have been exhausted, the late reaction also increases the potential for reversals, i.e., a deceleration immediately followed by an acceleration (as seen in Fig. 14), to compensate for the change in the spacing error.

6.2 Full data set results

The results for the full data set, including all 175 scenarios, are shown as boxplots in Figs. 15 and 16 and summarised in Table 3. The boxplots are drawn according to the Tukey convention and grouped by the weight factor setting: Original (Orig.), T_Opt and F_Opt. The outliers are indicated with red ‘+’ marks, the mean values are indicated by a black ‘x’ mark and reference or target values, where applicable, are indicated by a grey dashed line.

Figure 15. Box plots for each parameter by setting.

Figure 16. Box plots for each parameter by setting (continued).

Table 3 Median, mean and standard deviation values for each parameter by setting