Performance profiling for predictive train schedules

Journal of Rail Transport Planning & Management 4 (2014) 98–114

Contents lists available at ScienceDirect

Journal of Rail Transport Planning & Management journal homepage: www.elsevier.com/locate/jrtpm

Performance profiling for predictive train schedules R.L. Burdett ⇑, E. Kozan Decision Science Discipline, Queensland University of Technology, Australia

a r t i c l e

i n f o

Article history: Received 29 August 2014 Revised 10 November 2014 Accepted 17 November 2014 Available online 4 December 2014 Keywords: Schedule analysis Sensitivity analysis Profile of performance Robustness Railways

a b s t r a c t In this article an alternate sensitivity analysis is proposed for train schedules. It characterises the schedules robustness or lack thereof and provides unique profiles of performance for different sources of delay and for different values of delay. An approach like this is necessary because train schedules are only a prediction of what will actually happen. They can perform poorly with respect to a variety of performance metrics, when deviations and other delays occur, if for instance they can even be implemented, and as originally intended. The information provided by this analytical approach is beneficial because it can be used as part of a proactive scheduling approach to alter a schedule in advance or to identify suitable courses of action for specific ‘‘bad behaviour’’. Furthermore this information may be used to quantify the cost of delay. The effect of sectional running time (SRT) deviations and additional dwell time in particular were quantified for three railway schedule performance measures. The key features of this approach were demonstrated in a case study. Ó 2014 Elsevier Ltd. All rights reserved.

1. Introduction 1.1. Background A train schedule is a plan of all train movements that are supposed to occur in a railway system over a given period of time. Unfortunately operational conditions regularly cause actual events to deviate quite considerably from the schedule. The schedule is therefore only a prediction of what will happen. For environments such as this, the creation of a robust schedule is of considerable importance. The definition of robustness however varies from person to person and is viewed in different ways. For example in Salido et al. (2008) robustness is the capability ‘‘to absorb delays more readily’’, whereas in Cacchiani and Toth (2012) it is described in terms of whether a schedule ‘‘avoids delay propagation as much as possible’’. In Burdett and Kozan (2014) a schedule is robust if for specific key performance indicators (kpi) the schedule’s performance can be shown to be insensitive to delays of a specified level. Robustness can be further categorised in regard to the sensitivity of the schedule to both single delays and multiple delays, where a single delay refers to an incident to a single train, and multiple delays refer to separate simultaneously occurring incidents to more than one train. A schedule that is robust for multiple delays can be thought of as being more robust whereas a schedule that is robust to single delays is less and perhaps only partially robust. ⇑ Corresponding author. http://dx.doi.org/10.1016/j.jrtpm.2014.11.002 2210-9706/Ó 2014 Elsevier Ltd. All rights reserved.

In practice a ‘‘passive reactive approach’’ is typically taken and schedules are executed and conflicts and infeasibilities are dealt with as they occur, for example by rescheduling strategies. For schedules that are applied regularly (i.e. passenger services), the same delays and conflicts may occur and must be repeatedly resolved. With this in mind schedules should be analysed before implementation and courses of action that resolve infeasibilities and reduce degradation should be determined in advance. In other words, what can go wrong and how badly it can go wrong should be determined beforehand as well as what to do about it if it does. This is a ‘‘proactive scheduling’’ methodology instead of a ‘‘passivereactive’’ approach. This strategy may even lead to permanent changes in the schedule at the planning stage. A proactive scheduling approach can be implemented in one of two ways. Corrective strategies can be devised for any and all bad behaviour or they can be devised for specific expected bad behaviour. The second alternative obviously requires significantly less effort. 1.2. The research objective This article introduces an alternate way to characterise the strengths and weaknesses of a train schedule and its volatility and robustness to delays. The purpose of our characterisation is also to provide useful information that will allow an existing train schedule to be improved, for example as part of a proactive scheduling approach. In particular we specifically introduce the concept of performance profiling and the application of a sensitivity analysis for train schedules.

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

The main purpose of this article is hence to describe how best to obtain these profiles and to report on the nature of these profiles for different performance metrics. Furthermore this article explores whether profiles of performance can be used and are beneficial. Concrete details will be provided in due course in Section 2 however the main idea behind this approach is quite straightforward and is now described. The purpose of the sensitivity analysis is to identify the impact of delays to different ‘‘individual’’ parts of the schedule. For example each operation in the schedule is separately analysed and a unique profile of the schedules performance is created and stored. The advantage of this approach is that it is conceptually simple and straight forward. It is in essence just a sensitivity analysis of the processing times within the schedule. In theory the operations that have the worst profiles (i.e. that result in the greatest deterioration in schedule performance) should first be concentrated upon in some type of schedule refinement process. That topic is however for another paper. To be perfectly clear, there is no single output profile, but many profiles. The development of a single profile or metric for robustness, that summarises all the obtained profiles, has not been considered for a variety of reasons. The foremost reason is that the likelihood of specific operations being delayed is probabilistic and open to argument. This approach is not an approximation and does not consider a set of randomly generated incidents; the impact of all possible values of delays is quantified. The term ‘‘delay’’ is used to characterise additional dwell time and sectional running time deviations. The approach assumes that the original schedule is not altered and recovery techniques are not applied. Those actions however could also be integrated, perhaps in another stage of analysis. This approach is to be used offline to improve planning. However to our knowledge there is no reason why it cannot be used online to help with real-time control. This approach is a first step towards addressing the more complex situation involving multiple delays. In practice delays to trains hardly ever come in isolation, however the effect of single delays is still useful, as that information is indicative of the impact of delays to specific trains. An analysis of multiple delays is beneficial because it shows the overall effect, but its drawback is that it hides the exact impact of specific trains. In this article a disjunctive graph model of train operations is used as a delay propagation and evaluation mechanism. Other alternatives exist like simulation. To our knowledge simulation is most often the primary tool to analyse the effect of different types/ magnitudes of delay. For a set of likely delays and deviations, the expected performance of the schedule can be reported. The results however are only as good as the ‘‘incidence’’ generator that attempts to approximate the frequency and occurrence of realworld events. Even so, only a small subset and not all possibilities can be considered. The approach in this paper is anticipated as a superior, more accurate and more computationally efficient strategy to replace, enhance or supplement a standard simulation approach. Existing simulation models for railways are commonly based upon the perturbation of train arrival and departure times at different locations. That is, the analysis is directly applied to the train schedule and either train conflicts are counted or randomly resolved to give some idea of the overall performance. Simulation has not been performed in terms of the sequencing representation to our knowledge. Three performance measures were considered in this article. They are schedule time horizon or makespan (cmax), total train delay (ttd) and total time window violation (twv). The total time window objective is used to characterise adherence to the original schedule where the scheduled entry and exit times are use to define strict time windows for each operation. These performance measures were chosen because they are the predominant way of summarising the most important features of a train schedule. In

99

this article profiles are separately created for each performance measure. However in reality there is a relationship between the different performance metrics. For example, changes to the schedule that reduce the makespan can also simultaneously reduce the ttd and the twv. As this article is a proof of concept, we have first chosen to analyse each performance measure independently. In the future a combined approach should be investigated. A multi-objective approach is certainly feasible for instance. Search strategies are required to facilitate the sensitivity analysis. These are also addressed in Section 2. Section 3 follows with a discussion of how to interpret and use the results of the sensitivity analysis. In Section 4 the information from the sensitivity analysis is used to define costs of delays for railways. In addition the expected impact of an operation delay is considered in Section 5. A case study is provided in Section 6 for demonstrative purposes. Finally conclusions and research directions are given in Section 7. 1.3. Related work This article is a companion to Burdett and Kozan (2014) in which a sensitivity analysis was proposed to determine the affected operations after a single delay. The information provided by that analysis provides another measure of schedule robustness and also provides control information that can be used when delays occur in practice. Upon completion, the sets of affected operations can also be used to define the impact of all delays without further disjunctive graph evaluations. The techniques developed in Burdett and Kozan (2014) and in the current paper form a complete accounting of the effect of single delays. For example they describe the impact on performance metrics, on other operations, and on schedule feasibility. This article and our companion paper utilised the hybrid job shop scheduling framework of Burdett and Kozan (2010a,b). That framework was used to characterise train operations and to evaluate train schedules. In that framework sections of rail are viewed as machines, passing loops (sidings) are viewed as capacitated buffer machines, train journeys are jobs, and operations are train movements on particular sections of rail. Sectional running times and dwell times of trains contribute to the processing time of train operations as does the length of the train by way of a time lag. The schedule is represented by machine sequences. The schedule is therefore a temporal realisation of the sequences and is obtained by evaluating a disjunctive graph using a non-delay scheduling policy. The disjunctive graph is simply a directed graph that contains nodes and arcs. Each node refers to an operation. Each arc refers to track occupation precedence or else precedence between different operations of the same train. Both nodes and arcs have numerical values attached to them, typically called weights. The nondelay assumption prepares for unexpected events and ensures that trains are scheduled as early as possible and are not restricted from entering a section if it is feasible to do so (i.e. unforced idle time is not allowed). The train scheduling model and approach of Liu and Kozan (2009, 2011) is a variation of our approach that treats sidings as parallel machines. It could also be used as the basis of this approach. Similarly the alternative graph model of D’Ariano et al. (2007) could be used. Train scheduling can also be used for capacity level identification in railways. A review of capacity determination methods can be found in Burdett and Kozan (2006). A review of sensitivity analysis and robustness related research was provided in Burdett and Kozan (2012, 2014). In summary, the problem of determining criticality and sensitivity in activity networks has been considered (Elmaghraby, 2000). The degradation of the performance of a static scheduling algorithm in the presence of online disturbances has been devised (Penz et al., 2001). Leon et al. (1994) has analytically defined robustness measures as a function of expected delay for the case where only a single disrup-

100

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

tion occurs. For job shop problems the calculation of a ‘‘stability’’ measure to determine the credibility of an optimal schedule has been developed (Sotskov et al., 1997). In that approach the question is to determine such changes in the processing times of operations as do not destroy the optimality of the schedule. The stability radius denotes the largest quantity of independent variations of the processing times. Within the train scheduling field, this topic has received increasing levels of attention. Cacchiani and Toth (2012) have overviewed train scheduling with particular emphasis on those that provide robustness, and advantages and drawbacks were identified. Takeuchi et al. (2007) considered the qualitative and quantitative characteristics of robust train schedules. They advocated an approach based upon the simulation of incidences (i.e. delays) and the evaluation of a robustness index. Their robustness index is defined as the expected value of the total increase in overall passenger disutility, which is a measure of the reduced level of transportation service to passengers. Their robustness index however is only associated with the affect to passengers and is not generic. Hence it is not applicable (as is) to the vast assortment of train scheduling applications/scenarios and to any other scheduling kpi. Salido et al. (2008) proposed two formulas for comparing schedule robustness. They are based upon the topology of the railway infrastructure and buffer times. These equations however only allow a comparison to be made between two schedules, i.e. they report which schedule is more robust. They do not return a measure for the actual robustness of the schedule per se. D’Ariano et al. (2008) have introduced the idea of a ‘‘flexible schedule’’ the purpose of which is plan fewer trains in practice, in order to provide a greater freedom for the recovery of disturbances (i.e. by dispatchers). That idea was compared with the practice of inserting buffer time, which also provides schedule robustness. Greedy heuristics and a branch and bound approach were developed to create flexible schedules. This article concluded that flexible schedules are better than rigid ones. Fischetti and Monaci (2009) considered optimisation problems in which the exact value of input data are uncertain. A new modelling framework called ‘‘Light Robustness’’ was devised. It couples robust optimization with stochastic programming. The new approach results in reduced computational effort and time and is more flexible than the over conservative robust optimization approach that assumes each parameter takes (i.e. at the same time) its worst value. For demonstrative purposes a train scheduling problem was considered. A separate simulator was developed to measure the reaction and performance of a schedule to different delays and deviations. Hence robustness was measured in terms of the level of adjustment required (i.e. the difference between the scheduled and actual times). Comprehensive details of the train scheduling approach and the numerical investigation are described in Fischetti et al. (2009). In that article four stochastic methods were compared for improving the robustness of a nominal schedule. A simulation engine was used as a secondary device to verify solution robustness. Liebchen et al. (2009) proposed a new paradigm for robust train timetabling. In essence the new paradigm labelled ‘‘recoverable robustness’’ considers the joint optimisation of a plan and a strategy for limited recovery. This approach consists of a planning phase and a recovery phase. The recovery phase solves a feasibility problem and attempts to maintain the original schedule with a defined level of computational effort. Liebchen et al. (2010) have also considered the creation of robust train schedules and have attempted to bridge the gap between delay resistant timetabling and delay management by considering two difficult optimization problems. An alternative objective criterion was used and resistance to delay was measured by simulating individual delays and solving the corresponding delay management problem.

Railway schedule stability analysis has also been performed using ‘‘max-plus’’ algebra. For example Goverde (2007) proposed an approach to evaluate timetable realizability and stability using max-plus spectral analysis. Critical path algorithms were also used to quantify robustness. Train schedules were described as stable if a delay can be compensated for by time reserves elsewhere in the schedule. Goverde (2010) also modelled a railway system as a discrete event dynamic system using max-plus algebra. In particular models and algorithms were proposed to compute the propagation of initial delays over a periodic railway schedule. A timed event graph (i.e. petri net) representation was exploited and different types and magnitudes of delay were considered. In that approach delay propagation was minimized using buffer times. This allows the schedule to be recovered sooner. This articles approach has similar objectives to that of Goverde (2007, 2010). The nature of the two approaches however is vastly different. In addition, only a sensitivity analysis has been developed here. A recovery approach was not developed as that is a re-scheduling approach and the point of this article is to evaluate the schedule as it was originally intended.

2. Details of the sensitivity analysis This section describes how train schedules are analysed. The terminology used in this section and later ones will be introduced as needed, however a full listing can be found in Appendix A. In this article a machine scheduling framework has been utilised. This means that each train (i) is viewed as a job, each section (j) is a machine, and each train movement is classified as an operation oi,k with processing time pi,k and speed of travel vi,k. The index k refers to the position or stage within the train’s itinerary where the last operation of a job is denoted by Ki. The processing time is predominantly the sectional running time (SRT) or else the section occupancy time (SOT). As each train has length, a time lag denoted by lagi,k is added to the processing time of the operation. That lag defines the time difference separating the front of the train and the rear at each position on the railway. Without loss of generality the length of a train is hereby denoted as Ltrain and the i length of a section by Lsect . The superscripts are used to differentij ate between lengths of trains and lengths of sections, and sub scripts refer to particular trains or sections. As the ratio Ltrain =Lsect i j is used in a number of derivations and formulae, it is denoted by Li;j or more generally as Lj when the train type is not specified. Each operation is scheduled and its start and end timings are denoted by (bi,k, ei,k). The essence of the proposed approach is to determine for each train operation oi,k a complete account of delay. More specifically for different amounts and intervals of delay, the effect on a given performance metric (i.e. objective function) should be defined. The outcome is a complete profile of performance or in mathematical terms a function yobj i;k . Of particular interest is how the objective criterion is altered by different levels of delay. The value of delay that causes the schedule to become infeasible (i.e. the point of infeasibility), and the value of delay that causes a deterioration in the objective function (i.e. the point of increase) are sought. In a machine scheduling framework, delaying the entry time of a train increases the exit time on the previous section. Likewise delaying the exit time increases the entry time on the successor section. Increasing the processing time of a train operation by Dt in a schedule may have a different effect to delaying the exit time of the operation by Dt time. This is because the time for the rear of the train to depart a section is different in each case. After the occurrence of a delay, the original schedule is either feasible or infeasible and the objective value (performance measure) has increased, remained static or it has decreased. Therefore

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

for a single operation delay there are six possible results. After a delay, the schedule must be recomputed using the disjunctive graph and the non delay scheduling policy. In the context of train scheduling, infeasibility is associated with a cycle (i.e. precedence impossibility), an overloaded passing loop (i.e. a buffer overload) or both. In this analysis, infeasibility is only manifested as buffer overloads, as delays and SRT deviations do not cause precedence impossibilities. Readers should also note that buffer overloads are a consequence of the model and are not an actual reality. For delays in a specified interval [tl, tu] the analysis returns a vector of two tuples of the form (interval, status) where interval and status are themselves 2 tuples. For example the result of an analysis could be the following:

ððt l ; t 1 Þ; ðfeas; no deteriorateÞÞ; ð½t1 ; t 2 Þ; ðfeas; deterioratesÞÞ; ð½t2 ; t u Þ; ðinfeas; deterioratesÞÞ The schedule may be called stable during intervals where the objective function is static. Similarly a schedule may also be defined as stable during intervals where feasibility is maintained and unstable when it becomes infeasible. To determine the effect of a single delay Dt the following steps are taken: For the last operation of a job, a delayed departure from the section is equivalent to a predefined dwell. The dwell parameter of oi,k is incremented by Dt as is the value xi,k which is the weight of the associated node in the disjunctive graph. For all other operations of a train, a release time is given to the successor operation. More specifically, entry to the next section is delayed by assigning the new release time: rlti,k+1 = bi,k+1 + Dt. Finally the time lags of the train are recomputed. The determine the effect of a SRT deviation the following steps are taken: The processing time parameter and the node weight of the associated operation are both incremented by Dt, i.e. pi,k + = Dt and xi,k + = Dt. The time lags are then recomputed. If the operation is in the last position of the associated machine sequence then the arc to the sink node must be updated with lagi,k. Any disjunctive arc emanating from the operation node must also be updated. Burdett and Kozan (2010a,b should be consulted for more details about the exact details of disjunctive arc calculations. In order to apply the sensitivity analysis two conditions must be met. Lemma 1. (Dt > 0) _ ((k = Ki ^ dwelli,k + Dt P 0) _ (k < Ki ^ rlti,k+1 + Dt P 0)) for delays.

Proof. All positive delays are valid. For negative delays one of two conditions must be satisfied. The first case is for the last operation of a job and the requirement is that the delay nullifies (cancels out) an existing dwell. The second case is for operations that have a successor operation and the requirement is that the delay nullifies an existing release time for the successor operation. In other words there is a release time that can be ‘‘undone’’. h Lemma 2. pi,k + Dt > 0 for deviations in sectional running. Proof. It is impossible to travel across a section in negative time. Hence operations must have a positive ‘‘resultant’’ processing time. h In both cases the disjunctive graph must then be re-evaluated in order to obtain a new schedule. The computations can be reduced using the topological ordering of the operations. For example in a traditional job shop the operations that occur before the delayed operation in the topological ordering do not need to be rescheduled. For train scheduling problems this is also true. However unlike classical job shops, some other operations may need to be

101

re-scheduled even though their associated node in the disjunctive graph occurs before the ‘‘altered’’ region in the topological order. This is because operation exit times are determined from a later operation of the same job. In other words if this later operation is altered (i.e. the longest path value changes) then the predecessor operations exit time is affected. Determining all values of Dt and their effects is more difficult and computationally demanding. The selection of upper and lower  limits for delay (i.e. denoted by tþ lim and t lim ) has a significant effect on the sensitivity analysis. The number of function (i.e. disjunctive graph) evaluations that are necessary increases as these values are increased. The amount of positive delay that is relevant for an operation can be bounded by the schedules original horizon and the operations position within the schedule. For example that value is: cmax  (ei,k  delayi,k). Hence an operation that is performed towards the end of the schedule has less opportunity to be delayed than an operation towards the beginning. It is however possible to þ define some other bounds, say Dt  i;k , Dt i;k , which are to be used in the sensitivity analysis. Consequently tþ lim may be defined in the followþ ing way: t þ lim ¼ maxðcmax  ei;k þ delayi;k ; Dt i;k Þ. For negative delays a reduction greater than or equal to the original SRT cannot be  entertained. Therefore, t lim ¼ maxðpi;k ; Dt i;k Þ. It may also be necessary to restrict reductions in SRT because the associated speed is too high for the type of train, i.e. pi;k < ðLsect =v min j i;k Þ < Dt < 0 and sect t =v min lim ¼< ðLj i;k Þ. For obvious reasons the amount of negative delay is typically much smaller than the amount of positive delay.

2.1. The makespan objective 2.1.1. Positive delays Numerical investigations have generally shown that an infeasible train schedule does not become feasible again by further delaying an operation. Consequently the smallest number of intervals that results from a sensitivity analysis is one and this occurs when a schedule becomes immediately infeasible as a result of an operation delay. It has also been found that once deterioration occurs in the makespan, then there is no other ‘‘later’’ interval where deterioration does not occur. Generally the affect on the objective function value is not valuable once the schedule becomes infeasible. The total number of is at most three. The three forms to  intervals  þ þ be identified are: ð 0; t þ lim Þ; ðð0; t 1 Þ; ½t 1 ; t lim Þ, ðð0; t 1 Þ; ½t 1 ; t 2 ; ½t 2 ; t lim Þ. Two values of time (i.e. t1, t2) need to be found. It is assumed that the schedule is initially feasible and has an objective value of Z. In   the interval 0; tþ lim , delays and their effects can be characterised for this objective criterion according to Table 1. The seven cases are defined analytically in Table 2. For a specific delay, the feasibility of the solution is denoted by the Boolean parameter f. The associated objective function value is denoted by z. The subscripts are used to describe specific situations. For example: fu  f(tu), zu = z(tu). By default tl = e and tu ¼ tþ lim and e is a specified tolerance. For the makespan objective, the effect of delays and deviations is a piecewise linear function that is strictly increasing after the point of increase, i.e.

Table 1 Effects of delay. Case

Feasibility

Makespan

1 2 3 4 5

Always Always Always Always

Feasible until t þ inf

No increase Continuously increasing – Increases from time tþ Inc: Static while feasible

6

Feasible until t þ inf

Continuously increasing while feasible

7

Feasible until t þ inf

þ Static to tþ Inc: . Continuously increasing until t inf

feasible feasible infeasible feasible

102

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

Table 2 Summary of sensitivity analysis results.

Table 4 Summary of sensitivity analysis results.

Case Condition

Intervals

Interval feasibility

Increase in cmax

Case Condition

1

fl ^ fu ^ (zu = zl)

 ð0; tþ lim

Yes

No

2

fl ^ fu ^ (zl > Z)

ð0; tþ  lim

Yes

Yes

3

:f l ^ :f u fl ^ fu ^ (zu > zl)

ð0; tþ  lim

No

NA

þ þ ð0; tþ Inc: Þ; ½t Inc: ; t lim 

(Yes, Yes)

(No, Yes)

1 2 3 4

(Yes, No)

(No, –)

5

(Yes, No)

(Yes, –)

6

(No, Yes, –)

7

4 6

fl ^ fu ^ (tinc P tinf) ð0; tþ Þ; ½t þ ; tþ  inf inf lim fl ^ fu ^ (zl > Z) ð0; tþ Þ; ½t þ ; tþ  inf inf lim

7

fl ^ fu ^ (tinc < tinf)

5

þ þ þ þ (Yes, Yes, ð0; tþ Inc: Þ; ½t Inc: ; t inf Þ; ½t inf ; t lim  No)

  zðDtÞ  zðDt þ 1Þ8Dt 2 t Inc: ; tþ lim . Once the objective starts increasing then it continues to do so. The proposed sensitivity analysis can be computationally demanding as the number of operations increases. An analysis of those operations more prone to deviations and delay could possibly be concentrated upon. This would significantly reduce the computational overhead while still providing a reasonable sensitivity analysis. 2.1.2. Negative delays In this section reductions (negative deviations) that result in earliness, are considered. In other words for all deviations   Dt 2 tlim ; 0 where t lim < 0, the effect on the schedule is identified. Secondary points of infeasibility and increase are now required. For the makespan objective, a deviation results in no effect whatsoever, an immediate improvement and then no improvement, or an immediate and continual improvement. Negative dwells can occur, but only if there is an original dwell that can be undone. From numerical investigations it has been found that for Dt < 0, if there is going to be any reduction (i.e. improvement) in makespan, then it will occur straight away. For example the reduction of a critical operation will result in an immediate improvement, whereas the reduction of a non critical operation will have no effect. In the interval ðt lim ; 0Þ SRT deviations and their effects can be characterised according to Table 3. The seven cases are defined analytically in Table 4. For this analysis tu = e and tl ¼ t lim . The application of the sensitivity analysis on different case studies has shown that fewer graph evaluations are required to determine the effect of negative deviations as opposed to positive deviations. The reason for this is that the majority of the effects can be classified by Case 1 to 4 which do not require points of infeasibility, and hence more graph evaluations to be determined. The results of separate positive and negative analysis can be combined to give nine general scenarios and these are shown in Fig. 1. The yaxis describes the makespan and the x-axis refers to the value of delay. These functions could alternatively be plotted in terms of the change in makespan or even the relative change. Instead of

Table 3 Effects of delay. Case

Feasibility

Makespan

1 2 3 4 5

Always feasible Always feasible Always infeasible Always feasible Feasible until t inf

6 7

Feasible until t inf Feasible until t inf

No decrease Continuously decreasing – Stops decreasing from time t  Inc: Continuously decreasing while feasible No decrease while feasible Continuously decreasing to t Inc: . Static until infeasibility

Intervals

Interval Decrease feasibility in cmax

ðt lim ; 0Þ ðt lim ; 0Þ ðt ; 0Þ    lim  t lim ; t Inc: ; ðt Inc: ; 0Þ     fl ^ fu ^ (zl < zu < Z) t ; t inf ; ðt inf ; 0Þ    lim fl ^ fu ^ (zu = Z) ; t t inf ; ðt inf ; 0Þ  lim         :f l ^ f u ^ ðtinf < tInc: Þ tlim ; t inf ; t  inf ; t Inc: ; ðt Inc: ; 0Þ

fl ^ fu ^ (zu = zl) fl ^ fu ^ (zl < Z) :f l ^ :f u fl ^ fu ^ (zl < zu)

(Yes) (Yes) (No) (Yes, Yes)

(No) (Yes) (–) (No, Yes)

(No, Yes)

(–, Yes)

(No, Yes)

(–, No)

(No, Yes, Yes)

(No, No, Yes)

performing two separate sensitivity analysis, a single combined analysis could be also be performed. Each of the nine cases shown in Fig. 1 could be separately identified. Points of infeasibility are separately determined from this combined analysis. There are one to four intervals in each case that are static, increasing or decreasing. The cases may be described in the following way: Case (a), (b), (d) & (i): zl < Z < zu Case (c) & (h): zl < Z = zu

Case (e) & (g): zl = Z < zu Case (f): zl = Z = zu

It is not difficult to distinguish between (c) and (h); the gradient of the line between 0 and t lim can be used to obtain an intermediate point that is then evaluated and compared with the associated value on the line segment. To distinguish between (e) and (g) the same process is applied; the gradient of the line between 0 and tþ lim is used to obtain an intermediate point that is then evaluated and compared with the associated value on the line segment. To distinguish between (a), (b), (d) and (i) two projections must be performed. In summary two to four objective function evaluations are needed per operation to identify the correct case; two for the static case (f), three for case (c), (e), (h) and (g), and four for the remainder. 2.1.3. Extension Numerical investigations have shown that the gradient of the slope, from the point of increase to the upper bound, is one or is close to one. In other words for every second of delay, the makespan increases by approximately one second. Therefore, the point of increase can be exactly determined without need for a 1D search. In order to do so, a number of lemma’s are required, as shown below. Lemma 3. For deviations in sectional running time and for additional delays, tþ ¼ t u  ðzu  zl Þ=ku where tu ¼ t þ , the gradient is inc lim ku ¼ ðzu  zn Þ=ðt u  tn Þ, and tn is a chosen point such that tþ < tn < tu . inc Proof. See Fig. 2. Fig. 2b shows the position of the point of increase in relation to a possible gradient of one. h Lemma 4. For deviations in sectional running time the gradient of the function representing the change in makespan is no more than ð1 þ Lj Þ and no less than Lj where index j refers to the kth section traversed by train i (i.e. j = mi,k). Proof. The velocity and time lag are determined as follows:

v i ¼ Lsect =pi;k j

and lag i ¼ Ltrain =v i ¼ Lj  pi;k

ð1Þ

For an increase in SRT of one minute, the velocity decreases and the time lag increases:

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

103

Fig. 1. Complete set of functions for the makespan objective.

(a)

(b) Fig. 2. Determining the point of increase.

0 v i ¼ Lsect =ðpi;k þ 1Þ and lag i ¼ Ltrain =v i ¼ Lj ðpi;k þ 1Þ j

ð2Þ

The difference between successive time lags is hence: 0

lag i  lag i ¼ Lj  pi;k  Lj ðpi;k þ 1Þ ¼ Lj

ð3Þ

In summary, for a SRT deviation of one minute the exit time of the operation increases by 1 þ Lj . The change in velocity also causes a change of Lj in the time lag of the predecessor operation. In other words the predecessor operation will also finish at a later time. Any increase in makespan is manifested from one of these two sources. If it is from the former, then the makespan will increase by 1 þ Lj , if it is the latter then the makespan will increase by Lj . h Lemma 5. The actual gradient lies between Lj and 1 þ Lj . The point of increase is determined using the maximum possible gradient and a test point denoted by tn in the following way: t þ Inc: ¼ t u þ e  ðzu  zl Þ=ku where ku ¼ ðzu  zn Þ=ðt u  tn Þ and t n ¼ tu  ðzu  zl Þ=ð1 þ Lj Þ. tþ Inc:

Proof. If other gradients were used then tn could occur before and this would give an invalid ku value. The point of increase would then be invalid too. h In summary, the effect of the time lag, in determining when the makespan is affected, or not, arises. An operation that does not directly increase the makespan after an SRT deviation typically results in a gradient of one. These operations delay succeeding

operations which directly increase the makespan one unit of time for every unit of deviation (delay). Operations which directly increase typically have a gradient greater than one. Slightly different behaviour however is often exhibited by the last operation of a job, or buffer operations. The slack time of an operation also indicates when the makespan will increase if it is right shifted. For example the operations start time can be postponed by this amount before the makespan increases. For the makespan objective the slack time may also be used to roughly identify when an SRT deviation increases it. The slack time is not the ‘‘point of increase’’ exactly. The reason for this is that the operations duration increases as a result of an SRT deviation. Hence its completion time will not match the completion time signified by the slack value. The slack time is only valid for an operation of a specific duration otherwise it is invalid. 2.2. The total train delay objective For the total train delay objective and a non-delay scheduling policy, a positive SRT deviation may result in a decrease in total train delay. This is in complete contrast to the makespan objective which could only stay as is or increase as a result of an SRT deviation. In general the total delay may stay static, increase, decrease, or decrease and then increase. However the total delay cannot increase and then either stabilise or decrease. In the interval   þ ð0; tþ lim ; ðt lim ; 0 or ðt lim ; t lim Þ, delays and their effects can be characterised according to Table 5.

104

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

Table 5 Effects of delay. Case

TTD

1 2 3 4 5 6

Remains static Continuously increasing Continuously decreasing Static until some point, then continuously increases Decreases until some point and then remains static Decreases until some point, then remains static, then continually increases Decreases until some point then continually increases

7

Case 1 typically occurs when an operation has no successor. Consequently any delay (or deviation) has no affect on any other part of the schedule. In Case 2 a delay immediately delays other operations. Case 3 is a part of Cases 5–7 and only occurs in isolation for a small interval of time. The situation in Case 4 occurs when an operation is not directly adjacent to its successor operation in the sequence, i.e. the entry time of the successor operation is strictly greater than the exit time of the operation. Once this separation time is used up then the successor operation is delayed and only then does the total delay increase. In Case 6 and 7 a decrease in total delay occurs because a blocked operation is delayed within its enforced ‘‘blocking’’ time. Once the delay exceeds the blocking time then other operations are immediately affected. In the intervals where the objective function is increasing or decreasing, there may be several separate segments each with its own gradient. The function of delay is therefore piecewise linear. The reason for this is that within specific sub intervals the effect is constant and for every unit of delay the objective function increases proportionally. The gradient of the decrease is usually quite ‘‘shallow’’. Experimentation has shown that often the gradient (i.e. on the left) is negative one if the function is not static there. Furthermore the gradient of the decrease typically does not change, i.e. it is constant till the point tsta or the point of increase. The gradient of the increase varies quite considerably, i.e. from shallow to steep. Evidently the number of affected

operations affects the gradient. However an exact expression to link the two has yet to be devised. For the total train delay objective determining which case occurs is more difficult. The approach for the makespan objective, which was to compare Z, zl and zu, is now insufficient. To correctly identify what is going on, the gradients at each boundary must also be determined. The gradients at the two boundaries respectively are: kl ¼ ðzm  zl Þ=e and ku ¼ ðzu  zn Þ=e, and require that the following points be evaluated: (tm = tl + e, zm), (tn = tu - e, zn) and (tu, zu). Table 6 shows the conditions for correctly identifying the cases for positive, negative, or positive and negative delays. To completely describe Cases 4–7 extra work is however required to identify the minimum value and the point or interval where this occurs. No shortcut exists to do this like there was for the makespan objective. A 1D Golden Section search could be used in the worst case but a computationally superior ‘‘projection’’ method has been developed. The new approach which can be found in Appendix B, is specifically designed to find the minimum of a piecewise linear, non differentiable, convex uni modal function. The algorithm is similar to the Newton–Raphson method that finds the root of an equation. It differs in that two projections are performed at each step. This is shown in Fig. 3a. The gradients are also determined numerically. In more detail: The search starts from each boundary. The line on the left has a negative gradient and the line on the right has a positive gradient. A new point is determined every iteration by finding the intersection of the two currently defined lines. The point is evaluated and the gradient is determined. If the objective value equals the projected objective value then the minimum has been found and there is no static interval. If the gradient is zero then the minimum has also been found but there is a static interval around that point. Otherwise one of the existing lines is discarded and replaced with another. If the gradient is positive the line on the right is replaced. If the gradient is negative, the line on the left is replaced. The number of steps required is at most the number of line segments in the function. The number of function evaluations is at most twice the number of line segments, and this is due to the need to calculate the gradient at each point.

Table 6 Summary of sensitivity analysis results for positive and negative delays. Case

Condition

Intervals (For positive delay)

Intervals (For negative delay)

Total train delay

1

ðkl ¼ 0Þ ^ ðku ¼ 0Þ

Þ ð0; tþ lim

ðt lim ; 0Þ

(Static)

2

ðkl > 0Þ ^ ðku > 0Þ

Þ ð0; tþ lim

ðt lim ; 0Þ

(Incr)

3

ðkl < 0Þ ^ ðku < 0Þ

ð0; tþ Þ lim

ðt lim ; 0Þ

(Decr)

4

ðkl ¼ 0Þ ^ ðku > 0Þ

ð0; tInc: Þðt Inc: ; tþ Þ lim

ðt lim ; t Inc: Þðt Inc: ; 0Þ

(Static, Incr)

5

ðkl < 0Þ ^ ðku ¼ 0Þ

ð0; tsta Þðt sta ; tþ Þ lim

ðt lim ; t sta Þðt sta ; 0Þ

(Decr, Static)

6

ðkl < 0Þ ^ ðku > 0Þ

Þ ð0; tsta Þðt sta ; tInc: Þðt Inc: ; tþ lim

ðt lim ; t sta Þðt sta ; t Inc: Þðt Inc: ; 0Þ

(Decr, Static, Incr)

7

ðkl < 0Þ ^ ðku > 0Þ

Þ ð0; tInc: Þðt Inc: ; tþ lim

ðt lim ; t Inc: Þðt Inc: ; 0Þ

(Decr, Incr)

Fig. 3. Demonstration of projection method: (a) complete function (b) magnified region.

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

105

Fig. 4. Complete set of functions for the total train delay objective.

Determining the interval where the function is minimal is not trivial. If both of the ‘‘bounding’’ line segments are identified as being directly adjacent (i.e. to the static line segment) then the interval can easily be determined as the point of intersection between two different pairs of lines. However this does not always happen (see Fig. 3b) as the projection method tends to bypass some line segments. In order to exactly pinpoint the boundary of the static interval, the projection method should be applied twice with either the left or right line segment fixed as occurs when Case 4 and 5 is being considered. The results of separate positive and negative analysis can be easily combined to give 11 general scenarios and these are shown in Fig. 4. These functions could alternatively be plotted in terms of the change in total train delay or even the relative change. There are one to three types of intervals in each. There are still seven general shapes which the functions can take. Therefore an analysis over the complete positive and negative spectrum is more computationally beneficial. Numerical investigations have confirmed this and fewer function evaluations were necessary. A potential requirement of this analysis is perhaps the explicit identification of the interval and gradient of each line segment. The relevancy or more accurately the usability of this information is questionable. It is foreseeable that several variations of the projection method could be used to accomplish this, though the exact details of such an approach are outside the scope of this paper. Numerical investigations have shown that the gradient of adjacent line segments was quite close and a difference of one was typical. 2.3. The total time window violation objective Deviations from the original schedule may also be of interest and should be included in the sensitivity analysis. In order to measure this strict time windows may be defined for each operation according to their current entry and exit times. A similar and perhaps contrived situation occurs when previous time windows have been defined in which case the previous actions are unnecessary.

Unlike previous objectives, the schedule has a total time window violation of zero and this value cannot be improved upon. In other words the function is always increasing as a result of delay or deviation. This is because of knock on effects which cause more and more operations to be affected. Similarly as the amount of deviation increases in the negative sense, then the total time window violations also increase. There are no intervals where the objective value is static. This is because a time window violation is always incurred by the operation that is delayed. Four scenarios typically result and these are roughly shown in Fig. 5. The y-axis describes values of the twv. These functions linearly increase (decrease) or the function is piecewise linearly increasing (decreasing) on either side of the y axis. It should be pointed out that case (c) and (d) are not convex; more specifically the function is not convex for negative delays. It should also be mentioned that for positive delays the gradient is continually increasing but this is not always true for negative delays. For example in case (c) and (d) the gradient becomes shallower. The reason for this is that at some point the operation ‘‘earliness’’ ceases to affect other operations. In order to identify which of the four cases occurs, one must determine if the function is linear or piecewise linear, for positive and negative delays. A projection from each boundary is sufficient to determine linearity. In particular the intersection of the projected line and the zero total time window violation line should give a y-axis intercept of zero if the function is linear; otherwise it must be piecewise linear. Algebraically cu ¼ zu  ku t u ¼ 0 implies linearity for positive delays, while cl ¼ zl  kl tl ¼ 0 implies linearity for negative delays, where kl ¼ ðzl  zn Þ=ðt l  tn Þ, ku ¼ ðzu  zm Þ= ðt u  t m Þ, tm = tu/2 and tn = tl/2. The gradient of the projected line is determined by evaluating two other points tm and tn. Four objective function evaluations are therefore necessary per operation. If these points are too close to the boundaries, then numerical issues such as rounding error could cause the wrong case to be identified. The midpoint was consequently chosen. Numerical testing has shown this to be a good choice. Once the case has been identified, detailed information about the number of piecewise linear segments and the gradient of each

Fig. 5. Complete set of functions for the twv objective.

106

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

segment is a potential requirement for part of case (b) and (c) and all of all case (d). The relevancy or more accurately the usability of this information however is uncertain. The projection method may be altered so that all piecewise linear segments are determined. For example the projection method can be called on any sub intervals to identify line segments that have been ‘‘skipped’’. Points of infeasibility are still of interest and must be determined for this objective.

middle point is not feasible, then tm is the point of infeasibility. Otherwise the point of infeasibility must occur at the upper boundary, i.e. tu. To correctly find the point of infeasibility when Dt < 0, the left and right boundaries for the Bi-Section 1D search should be reversed; otherwise the following alterations are required:

Ifð:f m Þthenðt l ; f l Þ :¼ ðt m ; f m Þ; elseðt u ; f u Þ :¼ ðt m ; f m Þ; ==Increase left or right boundary

2.4. Search strategies The proposed sensitivity analysis requires that the value of delay, that causes a change in the objective value, or makes the schedule infeasible, be determined. For the makespan objective the point of increase can be determined analytically but this is not so for the other objectives and is certainly untrue for determining the points of infeasibility which occur independently of the objective criterion. To ascertain these points, separate one dimensional search must be performed. The problem is similar to that of finding the minimum of a function. Several numerical strategies may be used such as Golden Search, Fibonacci Search, Quadratic Search and Cubic Searches. All of these strategies however must be modified as they are primarily designed to find the minimum of uni-modal functions. It should also be noted that the Quadratic and Cubic Searches are particularly unsuitable for finding the point of infeasibility because there is no gradient information. A Bi-Section (Mid-Pt) search approach was selected instead of the others. It has been found to be superior to the Golden Section and Fibonacci search strategy because the interval of uncertainty is reduced the most at each step (i.e. it halves). Therefore the smallest number of objective function evaluations, are necessary. The Bi-Section method is primarily used to find the roots of a function. It generates a new middle point at each step, which is then evaluated. Based upon the comparison of the three objective function values, the new interval of uncertainty is (tl, tm) or (tm, tu). The signs at each boundary must be the opposite in order for the root to be correctly bounded, i.e. sign(zl) = sign(zu). For tolerance e the relative difference between successive midpoints and the relative difference between successive midpoint objective values is used as the termination criteria. The Bi Section approach may not be suitable for finding the minimum because three points are insufficient to tell on which side of the middle point the minimum will occur. Three points are only sufficient if zl < zm < zu or zu < zm < zl. The full details of the two variant Bi Section algorithms can be found in Appendix B. The new variants differ from the original algorithm by the way the new interval is selected at each step and by how the point of infeasibility or increase is identified once the termination condition has been met. Normally the minimum would occur at the midpoint tm. Makespan objective & points of increase: The function is static or increasing from left to right. If the value at the midpoint is greater than the value at the left hand boundary, then the point of increase must be to the left of the midpoint (i.e. decrease upper boundary). Otherwise the point of increase is to the right of the middle point (i.e. increase lower bound). Upon termination the point of increase is identified as the middle point tm if zm > zl. Otherwise it is the upper boundary tu, i.e. when zu > zm. Points of infeasibility: For Dt > 0 the function is feasible on the left and infeasible on the right. If the middle point is feasible then the point of infeasibility occurs to the right (i.e. increase lower boundary) otherwise it occurs to the left (i.e. decrease upper boundary). Upon termination the point of infeasibility could equally occur at any of the three points. To determine which one: If the lower boundary is not feasible then tl is the point of infeasibility. Otherwise if the left boundary is feasible and the

Determining the point of infeasibility here is different because the objective function is infeasible on the left and feasible on the right; previously it was feasible on the left and infeasible on the right. 3. Interpreting the results This section discusses how best to use the information provided by the sensitivity analysis. 3.1. Points of infeasibility The point(s) of infeasibility are particularly important as they show that some type of corrective action like rescheduling would be necessary in order to keep the schedule valid. Ideally while the schedule is being executed, ‘‘traffic controllers’’ should ensure that specific operations do not get delayed beyond levels that result in infeasibility. Determining points of infeasibility may be performed separately from the rest of the sensitivity analysis. It is also invariant to the objective criteria and therefore this part of the analysis need not be repeated for different objective criteria. The analysis described in aforementioned sections can be augmented by determining some additional information. For example, where the infeasibility occurs in regards to ‘‘elapsed schedule time’’ has not been provided yet. Similarly the extent of the infeasibility in terms of the number of buffer occupancy violation (BOV), or the buffer machine upon which the BOV will occur, could be recorded. This information should be added as an outcome of the analysis as it can be quite useful. BOV are however specific to the hybrid job-shop framework which we utilised, and is not present in other frameworks. The actual infeasibility occurs between the exit time of the ‘‘delayed’’ operation and the maximum completion time of the entire schedule, i.e. tactual 2 ðei;k þ t inf ; cmaxÞ. If the schedule inf becomes infeasible immediately or soon after the delay (i.e. tactual  ei;k is small) then the operation may be viewed as critical, inf though not critical in the usual sense which is that the operation lies on a critical path of some activity network. In addition the schedule may also be viewed as ‘‘locally unstable’’. Otherwise there is some room to manoeuvre and rescheduling (i.e. recovery activities) for example could be applied. To be viewed as unstable a schedule must be intolerant to any potential delays or at least a great number of potential delays. The degree of instability could be indicated by the number of critical operations, where criticality is defined as in the above paragraph. The number of critical operations could be used to compare competing schedules. As a rule of thumb it can be assumed that the level of infeasibility increases or stays static as the delay increases. It could be useful to identify the level of infeasibility at the point of infeasibility and at the upper (lower) bound. This would then give a reasonable indication of how the level of infeasibility is changing and what amount of rescheduling effort is required. For example the effort to resolve one BOV is relatively small in comparison to the effort of resolving say ten. An exact accounting of the level of infeasibility could be determined if necessary (i.e. a time versus number of BOV

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

chart could be plotted) but would require a large computational overhead if done for all operations. A way of altering / correcting the schedule so that the predicted infeasibility does not occur and the resultant schedule is not degraded, is outside the scope of this paper and is a source of future research. 3.2. Points of increase The points of increase (or change) are also important. For the makespan objective the points of increase define where the objective is no longer static and where delays start to have a detrimental effect. In other words the interval where delay is acceptable is indicated. For the total train delay objective the point of increase has a similar usage, however it is also used to show where the function stops decreasing and starts increasing, i.e. where it is minimal. For the total time window violation objective there is no unique point of increase. For the case where some deterioration in the objective is acceptable, the points of increase (change) can be used to pinpoint the exact value of delay that causes this effect. Fig. 6a shows a typical example. For the makespan objective this value of delay can be determined most easily because there is only one shape of function other than a completely static one, and the function is purely linear. For the total train delay objective there are more cases to consider. For the total train delay and total time window violation objective criteria the main difficulty occurs when there are several piecewise linear segments. A typical example is shown in Fig. 6b. Without this information it is hard to pinpoint the delay exactly. The proposed projection method would have to be applied to determine this value exactly. 3.3. Comparing cases In theory the obtained profiles can be ordered (listed) from best to worst or vice versa. This is necessary in a proactive corrective approach, as most effort should be directed at the worst operations. Operations that have been identified as having a large impact on the performance of a schedule should be ‘‘buffered’’. For example, additional slack should be introduced with a magnitude proportional to a specified tolerance that describes the acceptable change in the performance of the original schedule. For example if the delay causes a deviation of d in the objective and the tolerance  d is less than d minutes say, then buffering should be proportional to d   d. In the case of the makespan objective, and a gradient of one for the increase, the buffer size will be exactly d   d. For operations identified with the same profile type (i.e. case), the comparison is simplest. The operation with the greatest period of no change is best. Otherwise the operation with the smallest

107

gradient is best. Otherwise an operation is best if the schedule remains feasible longest when that operation is delayed. For comparing operations with different profiles it should first be noted that cases with guaranteed feasibility are better than those that become infeasible immediately. Otherwise it is like ‘‘comparing apples with oranges’’. However a case with a static period is better than one without. For cases that become infeasible, a larger point of infeasibility is better however the presence of a static interval should take priority in the comparison. When comparing cases with a point of infeasibility to those without: a comparison should first be made with respect to the presence of a static period, and then with respect to the gradient and duration of the linearly increasing interval. 3.4. Schedule refinement The sensitivity analysis provides profiles of performance which are then sorted and ranked. The operations that are most volatile and contribute to the greatest reduction in the schedules performance can hence be identified. In order to improve the schedule two courses of action could be taken. The first is to buffer volatile operations. For example, additional separation time should be added between train movements. Buffering will of course change the schedules timings and cause deteriorations in the performance of the overall schedule. In other words there is no improvement in the performance metric itself. The reliability of the corrected schedule will however be improved. If predicted delays do not eventuate then some of the original timings could re-instated and the schedules performance could be improved; but those changes would have to be performed dynamically as the schedule is implemented. The only other course of action is to actually change the sequences of the train movements which would completely change the timings of many trains within the schedule. This analysis does not currently identify how to change those sequences. 4. Quantifying the cost of delay A train access charging methodology was developed in Kozan and Burdett (2005). A similar approach is developed here that utilises the performance profiles explained in preceding sections. It should be noted that railway systems cost significant sums of money to build, operate, maintain and extend. Hence to recover these costs, train services must be charged. When train services are delayed positively or negatively, there can be significant financial repercussions. The cost of delay can be quantified for each objective criteria obj e {cmax, ttd, twv} and for each profile yobj . i;k Furthermore a more generic cost of delay can be defined that includes all three, though any one of these may be removed as desired. The details of the cost functions are now discussed.

Fig. 6. Determining the delay for a prescribed level of acceptable deterioration.

108

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

4.1. Cost of Changes in Capacity utilisation The profile ycmax is required here and allows a complete profile i;k of cost C1i;k to be defined for operation oi,k and all values of delay Dt. b and there is some If a delay changes the makespan from C to C function U1 for this charge, then the cost of delay could be defined in the following way:

C1i;k ðDt; C;

b ¼ U1 ð C b  CÞ CÞ

ð4Þ

If this amount is positive then a loss of income occurs; unless it is recovered through a surcharge to the offending train service or else split equitably amongst all train services running at the time or in the schedule. If this amount is negative then the additional income could be refunded or else retained. If there is a fixed charge 1 f then U (Dz) = Dz  f. If overhead costs are linked to the loss of capacity utilisation then the cost of delay could be defined in the following way:

b ¼ C1i;k ðDt; C; CÞ

b  CÞXT ðC T

i:e: U1 ðDzÞ ¼ Dz  XT =T

ð5Þ

 ¼ XT =T In this equation the overhead cost per unit of time is X and the overhead for an interval of time T is X T. This equation is based upon the observation that if the makespan changes then  to C bX  . The difference in the overhead cost also changes from C X  is constant, it is the same as these values is the loss or gain. As X a fixed charge. 4.2. Cost of additional train delay The profile yttd i;k is required here and allows a complete profile of cost C2i;k to be defined for operation oi,k and all values of delay Dt. A delayed train may cause delays to other train services. If a delay to b and operation oi,k of Dt changes the total train delay from D to D function U2 for this charge, then the cost could be defined in the following way:

b ¼ U2 ð D b  DÞ C2i;k ðDt; D; DÞ

ð6Þ

For a fixed charge f then U2(Dz) = Dz  f. If the cost of train delays was dependent on the train or service type (say gi) then there could be several U2 functions in use. Consequently the following cost function could be defined:

b ¼ C2i;k ðDt; D; DÞ

X

^ d Þ U2gu ðd u;k u;k

ð7Þ

8ou;k jd^u;k >du;k

^ signify the predicted and actual delay The variables di,k and d i;k to operation oi,k respectively. The actual delay for each operation however is not currently recorded in the sensitivity analysis. The total train delay per train or service type (i.e. Di) however could be more easily recorded. In that case the cost function could be redefined in the following way:

C2i;k ðDtÞ ¼

X b u  Du Þ U2u ð D

ð8Þ

8u

4.3. Cost of schedule non adherence The profile ytwv i;k is required here and allows a complete profile of cost C3i;k to be defined for operation oi,k and all values of delay Dt. If a delay to operation oi,k of Dt changes the total time window viob and there is a function U3 for this charge, then lation from V to V the cost could be defined in the following way:

b Þ ¼ U3 ð V b  VÞ C3i;k ðDt; V; V

ð9Þ

For a fixed charge f then U3(Dz) = Dz  f. It is possible to charge different amounts for different train or service type. The cost function could therefore be defined as:

C3i;k ðDtÞ ¼

X

b u  V uÞ U3u ð V

ð10Þ

8u

The two previous equations penalise any deviations from the original schedule. It is however possible to separate earliness and lateness and to charge them differently. If there is a fixed charge or function U3 for earliness and a fixed charge or function U4 for lateness then the cost could be defined in the following way:

C3i;k ðDtÞ ¼

X

X

^  b ÞÞ þ ðU3 ðb u;k u;k

8ou;k 2CB

ðU4 ð^eu;k  eu;k ÞÞ

ð11Þ

8ou;k 2CE

^ – b g and CE ¼ fo j^e – e g. The variables where CB ¼ fou;k jb u;k u;k u;k u;k u;k ^ bi,k and bi;k signify the predicted and actual start (i.e. entry) time for operation oi,k respectively. Similarly the variables ei,k and ^ei;k signify the predicted and actual finish (i.e. exit) time for operation oi,k respectively. If the cost of earliness and lateness was train type specific then there could be several U3 and U4 functions in use:

C3i;k ðDtÞ ¼

X

^  b ÞÞ þ ðU3gu ðb u;k u;k

8ou;k 2CB

X

ðU4gu ð^eu;k  eu;k ÞÞ

ð12Þ

8ou;k 2CE

The actual entry and exit times for each operation however are not currently recorded in the sensitivity analysis. The total earliness (i.e. earlyi) and lateness (i.e. latei) could more easily be extracted from the total time window violation and recorded. Therefore the cost function could be redefined in the following way:

C3i;k ðDtÞ ¼

 X d u  early Þ þ U4 ð late d u  lateu Þ U3u ð early u u

ð13Þ

8u

4.4. Cost functions In theory there are many possibilities for the functions Un (where n = 1, 2, 3, 4). In the context of railways they should be quite simple. Exhibiting different behaviour over different intervals of delay however is a very strong possibility. The shape of the resulting cost functions can also be quite diverse. The particular shape is determined by the interaction of the performance profile n yobj i;k and the charge function U . In general the following signs should be satisfied:

Un ðDzÞ > 0 8Dz > 0; n 2 f1;2;3; 4g; Un ðDzÞ < 0 8Dz < 0; n 2 f1;2g Un ðDzÞ < ðor >Þ0 8Dz < 0; n 2 f3;4g ð14Þ 2

A decrease in total train delay is always beneficial so U can return negative values. In some situations early entry and exit is also desirable and U3 and U4 may also be allowed to return negative values. As the values of Dz are typically increasing linearly as the value of delay increases, the output of function Un, and hence the value of cost, will also be increasing linearly. All that is required is to scale the input value appropriately and to possibly add some other fixed cost/charge. Hence: Un(Dz) = aDz + b. A function like Un(Dz) = a  expbDz + c, is also a possibility for sufficiently small values of |b|. If |b| is too large then the cost would increase (or decrease) too dramatically, and this makes little practical sense in a railway context. To use a function like this the question of why the cost is permitted to change more greatly as the delay increases, also requires further investigation. In the context of charges however, a function like this could be defined in order to

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

109

provide incentives for reducing delays in practice. A step wise function also makes sense because costs and charges may be deemed fixed over specific intervals. A comparison of these different functions is shown in Fig. 7b, c, d for the performance profile shown in Fig. 7a.

þ max running time deviations t  ¼ pmin i;k  pi;k and t ¼ pi;k  pi;k . In the event that no value is more likely than any other, a uniform distribution could be utilised. The probability that a particular delay in min interval [t-, t+] occurs is 1/(t+ - t-)  1/ ðpmax i;k  pi;k Þ. The expected impact for each operation is therefore as follows:

5. Quantifying the expected impact of delay

 Z tþ h i 1 E Iobj Zð:; tÞdt ¼ Z þ i;k tþ  t t

The expected impact of an operation delay is quantified in this section. In order to do this Z is defined as the original objective value. In addition Z is defined as the function that returns the objective value after operation oi,k is delayed by time D t, i.e. Zðoi;k ; DtÞ. It is also abbreviated as Zð:; DtÞ below. This abbreviation is also used in other appropriate places. Which type of delay is being considered is known in advance. The impact of an operation delay is given functionally as Iðoi;k ; tÞ which returns the value Zð:; DtÞ  Z. Variable Iobj i;k is similarly defined to record the value for a particular objective criterion. The expected impact of an operation delay may be quantified in the following way:

E½Iobj i;k  ¼

Z

tþ

Ið:; tÞPð:; tÞdt

ð15Þ

An analytical result can be derived for the integral if the profiles of performance determined from the sensitivity analysis are available,      i.e. Zð:; t  Þ; Zð:; t þ Þ; Z :; t þ Inc: ; Z :; t Inc: . This expression also requires þ þ the following values: t  ; t  Inc: ; t Inc: ; t ;. It should be noted that þ t ; t are the values of delays on each side of the x axis respecInc: Inc: tively, where a change in gradient occurs. For the makespan objective:

 Z t Z 0 h i Inc: 1 E Icmax Z ð :; t Þdt þ Zð:; tÞdt ¼ Z þ i;k tþ  t t t Inc: ! Z tþ Z tþ Inc: Zð:; t Þdt þ Zð:; t Þdt þ 0

t

where the extent of the possible delays is over interval [t-, t+] and (D., t) is a suitable probability density function. For sectional

ð16Þ

where

Z

b

Zdt ¼ a

ð17Þ

tþ Inc:

Zð:;aÞðb  aÞ if Zð:; aÞ ¼ Zð:; bÞ ba ðZð:;bÞ þ Zð:;aÞÞ if Zð:; aÞ – Zð:bÞ 2

Fig. 7. Comparing cost functions.

Fig. 8. The considered train schedule with analysis results overlayed. Trains are denoted as J0–J9 and cases are denoted by ‘C’.

ð18Þ

110

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

Proof. The function is piecewise linear and the gradient of each segment is positive, i.e. each segment is static or increasing in the positive direction. The area under each line segment is rectangular and or triangular and can be determined with basic equations. h As there are nine cases for the makespan objective (see Section 2.1), specific equations could be derived for each. Removal of the minus Z term also results in an equation for the expected completion time of the schedule. For the total train delay objective there are eleven cases (see Section 2.2) and negative gradients can occur. For the total time window violation objective there could be more than three or four line segments in a profile. Therefore the expected impact is further generalised in the following way:

h i E Itwv ¼ Z þ i;k

! X Z tiþ1 1 Zð:; tÞdt tn  t 1 i2½1;n ti

ð19Þ

where t n ¼ t þ and t 1 ¼ t  . More computational effort is required to identify the line segments. For most other probability distributions the following integral must be calculated numerically.

Z h i E Itwv ¼ i;k

tþ

t

Zð:; t ÞPð:; tÞdt  Z

Z

tþ

Pð:; tÞdt

ð20Þ

t

6. Case Study For demonstrative purposes a schedule consisting of ten trains was analysed. It is shown in Fig. 8 and has been obtained from Burdett and Kozan (2014). Only the most important parts of our approach are shown here; a full approach would consider more than one objective. The case study is part of the important North Coast Line in Queensland, Australia. This line runs from Brisbane to Cairns and is approximately 1600 km long in entirety. The portion considered here is approximately 60 km and involves 10 sections and 9 sidings. The line is bi-directional. Each siding is one kilometre long and the section lengths (in km) are as follows: 4.75, 1.13, 3.56, 4.72, 9.93, 2.9, 5.25, 2.48, 8.93, 5.55. Four trains travel at 80 km/h, three at 100 km/h and three at 120 km/h. Each train has 19 operations; hence our analysis consists of 190 components. Note that operations are numbered from 0 to 18, where O0 refers to movement on section S0 and O18 on S18. An analysis for the makespan objective and for increases to the SRT’s are now discussed. The makespan of the original schedule is 88.82 min. The complete profiles for the operations of three of the ten trains are given in Appendix C. In total the analysis constituted 2265 function evaluations. Individually, most operations required approximately 15. Each of the seven cases (as described in Table 1) occurred the following number of times respectively: 0, 16, 0, 38, 57, 20, 59. The exact distribution of cases to operations is shown in Fig. 8. In Fig. 8 the labelled operations refer to positions where the case changes. Operations that are not labelled are assumed to be the same as previous ones. Ultimately the best way to visualise this information would be to colour each operation according to the case, however this is not possible here due to the black and white policy. In summary most of the cases occurred consistently. Case 1 and 3 however did not occur at all. Case 4, 5 and 7 occurred most and constituted 154 occurrences in total. This is noteworthy because the objective does not increase immediately after a deviation in those cases. Without looking at the duration of the period whereby the objective is static, this is an indication of a level of robustness. The remaining 36 operations (i.e. spread between case 2 and 6) result in an immediate change to the objective and hence

indicate a level of insensitivity to delay. The duration of the static periods can be summarised by a histogram. For a bin array that increases in units of 5 min the frequencies are as follows: (89, 46, 17, 0, 1, 0, 0, 0, 0, 0, 0, 0). A decreasing pattern is shown with a high proportion of values (i.e. 89) in the zero to five minutes range. This shows a relatively high sensitivity to delays because deterioration starts to occur almost immediately. There is some tolerance to delays up to 15 min but very little indication of a general insensitivity to delay. The 190 profiles can be sorted as mentioned in Section 3.3. For positive SRT deviations, Case 1 is best and Case 3 is worst – but no operations exhibited that behaviour. Case 4, 5 and 7 are generally better than Case 2 and 6 because they have a static interval which is a marker for tolerance and insensitivity to delay. Case 2 is generally superior to Case 6 because Case 2 is always feasible, but for all intents and purposes, they can be combined as one group. Case 4 is generally better than Case 7 for the same reason, however if the interval of no change in Case 7 is longer, it may be regarded as better. Case 5 is generally better than Case 4 and 7 because it remains static until the point of infeasibility. However if the interval of no change is larger, Case 4 or 7 may also be regarded as better. Hence the worst to best ordering of cases is roughly (3, 6–2, 7–4-5, 1), and should be used as a starting point to sort the profiles. A secondary sort is then necessary to compare cases in the same group and of the same type. After applying the aforementioned logic the best and worst operations are summarised in Table 7. All of these results are verifiable visually by inspecting Fig. 7. Train J7 features prominently in the list of the worst operations as it occurs first in the schedule and could affect many other trains. Similarly train J9 also features because it is last. The other poor operations occur around the middle of the schedule on ‘‘critical’’ sections of rail. The best operations generally consist of the early or latter operations in train itineraries, and again occur in the middle of the schedule. This information informs us that the first 20 operations in Table 7 should be supervised closely and should not be delayed. It also tells us that those operations should be buffered if for instance corrective actions cannot be taken or enforced and they are delayed frequently. It should be noted that buffering will of course change the schedules timings and cause deteriorations in the performance of the overall schedule. To demonstrate cost calculations the following examples are provided. Table 7 Sorted Profiles. Worst operation at #1 and best at #190. #

Job

Opn

Case

Duration-1st interval

#

Job

Opn

Case

Duration-1st interval

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

J9 J0 J0 J1 J2 J4 J0 J4 J0 J9 J3 J7 J7 J7 J7 J7 J7 J7 J5 J9

O0 O9 O10 O8 O8 O8 O8 O9 O11 O1 O10 O3 O5 O2 O6 O7 O4 O8 O8 O2

2 6 6 6 6 6 6 6 6 2 6 6 6 6 6 6 6 6 6 2

0 0.794 0.794 0.794 0.794 0.794 0.795 0.882 0.883 2.375 2.559 2.786 2.786 2.787 2.787 2.787 2.788 2.788 2.79 2.875

171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190

J9 J8 J8 J8 J1 J2 J1 J0 J0 J0 J0 J0 J0 J0 J0 J2 J2 J1 J4 J2

O13 O3 O2 O1 O18 O15 O17 O0 O1 O2 O3 O4 O5 O6 O7 O17 O18 O15 O18 O9

4 7 7 7 4 7 4 7 7 7 7 7 7 7 7 4 4 4 5 7

9.805 10.02 10.02 10.02 10.418 10.623 11.682 12.038 12.038 12.038 12.038 12.038 12.038 12.038 12.038 12.419 12.599 13.041 20.626 52.461

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

Example 1. The profile of J0 (O0) is Case 7. The intervals are (0, 12.038), [12.038, 12.836), [12.836, 66.506] and the gradient is 1. The cost of delay is zero for the first interval as the objective is static. At the boundary of the second interval, the makespan is 89.62 min because the gradient is one. Hence for the second interval the cost of delay ranges from $0 to U1(89.62 88.82) = U1(0.8). If U1(Dz) = Dz  XT/T, X = $10,000 per hr and T = 1 h, then the cost of delay is C1i;k ¼ 0:8  $166:67 ¼ $133:33. In the third interval the schedule is infeasible, so corrective action is required. The cost of delay would therefore be based upon the quality of the corrected solution. An alternative approach is to define a ‘‘flat’’ cost of infeasibility.

Example 2. The profile of J6 (O0) is Case 4. The intervals are (0, 0.652), [0.652, 3.599] and the gradient is 1.19. The cost of delay is zero for the first interval as the objective is static. The makespan is 88.83 when the delay is 0.652 and 92.33 when the delay is 3.599. Hence for the second interval the cost of delay ranges from $0 to U1(92.33 - 88.82) = U1(3.51). The cost of delay is C1i;k ¼ 3:51  $166:67 ¼ $585:01. It should also be noted that 3.51  (3.599 - 0.652)  1.19 and hence the objective values are not required for this calculation.

7. Conclusions In this article methods have been developed that quantify the performance of a train schedule when it is subjected to delays of two types. The resulting ‘‘sensitivity analysis’’ can validate the claim that an existing schedule is robust and performs acceptably under various types and magnitudes of delay. The theory in this paper is generic and can be altered relatively easily for other scheduling problems, like job-shops and flow-shops. The main difference between this environment and others is only in the details of the disjunctive graph model that has been used as a delay propagation and evaluation mechanism. A non-delay scheduling policy was applied and the robustness of schedules created from its application, can be analysed and judged. The proposed sensitivity analysis can also be used as part of a proactive scheduling approach. It is based upon the idea, that if you know what will go wrong, and how bad it will go wrong in advance, then you can plan how to deal with it in advance. More specifically for single ‘‘isolated’’ delays, a profile of performance was investigated and quantified for three different objective criteria, namely makespan, total train delay and total time window violation. The results of the sensitivity analysis can also be used to partition operations from best to worst (see Case Study). This ordering describes which parts of the schedule are at most risk, and where rescheduling efforts are required. The results of the analysis can also be used to quantify the cost of delay provided that functions for earliness, lateness and individual train delay have been given. In this paper cost equations have been proposed for each of the objective criteria. In addition a more general measure which includes all three has been proposed. The approach presented in this article and in Burdett and Kozan (2014) is quite comprehensive however the field is still open to further developments. For example an analysis that takes into account both types of delay is an avenue for future research. The effect of multiple delays was not considered and is also a really important avenue. There appears to be a real need for an effective approach to evaluate the interaction between changes in the parameters of two or more operations, in train scheduling applications, and also

111

in other scheduling environments. How to efficiently use the proposed techniques to build a robust train schedule from scratch is also an avenue for future research. It is possible that the proposed sensitivity analysis could be used to explicitly define robustness. At present the entire schedule is analysed but there is no single metric of robustness and no single output profile. The sensitivity analysis is computationally demanding and whether it could be streamlined and used for example in a meta-heuristic approach, for evaluating a schedule and constructing a robust schedule, is an open question. Determining the points of infeasibility require the most computations. The function is most ‘‘uncooperative’’ as there is no gradient information. It is very difficult to predict where the point will occur, i.e. close to left or close to right boundary. Further research effort is perhaps required here. Appendix A. Nomenclature i, k oi,k, mi,k

pi,k, lagi,k, rlti,k, dwelli,k, delayi,k

vi, vi,k xi,k bi,k, ei,k cmax, ttd, twv ; Lsect ; Li;j Ltrain i j

Dt ; t tþ lim lim

Dtþ ; Dt  i;k i;k þðÞ

þðÞ

tInc: ; t inf ; tsta

fl, fu zl, zu, Z

kl ; ku d d;  yobj i;k Cobj i;k T f,

X

b C; C

train and train operation (i.e. stage) index kth operation of train i and the index of the associated section that is traversed. operation attributes: processing time, time lag, release time, dwell and delay time general train speed and specific speed on the kth section numerical value of the node in the disjunctive graph associated with operation oi,k schedule timings for operation oi,k, i.e. begin and end time objective functions train length and section length, ratio of lengths value of delay or deviation upper and lower limits for delay user specified limits for delay special values of delay within the profile, i.e. point of increase and infeasibility, point where profile becomes static solution feasibility parameters objective function value, i.e. lower, upper, original value of the schedule gradient values within the profile, i.e. lower, upper deviation in objective function and deviation tolerance profile of performance for operation oi,k, i.e. a function cost profile for operation oi,k, i.e. a function interval of time fixed charge and overhead cost values of cmax (continued on next page)

112

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

values of ttd

b D; D b V; V

U

values of twv

n

functions that describe the charge incurred when the schedules performance is altered. expected impact of delay

E½Iobj  i;k

functions of (oi,k, t)

I; Z

Appendix B. Algorithm pseudo codes

Algorithm 1. FindMinPLF ðoi;k ; tl ; t u ; zl ; zu ; kl ; ku ; t min ; zmin ; eÞ // The Projection Method Begin found:= false; while (found) begin ðcl ; cu Þ :¼ ðzl  kl tl ; zu  ku t u Þ; // Compute y axis intercept for left and right segments t m :¼ ðcu  cl Þ=ðkl  ku Þ; // Calculate the point of intersection zproj :¼ kl tm þ cl ; // Calculate the projected y axis value at point of intersection EvalDelay (oi,k, tm, zm, fm); // Determine actual y axis value and feasibility tn: = tm + e; // Define nearby point EvalDelay (oi,k, tn, zn, fn); // Determine y axis value and feasibility km :¼ ðzn  zm Þ=ðtn  tm Þ; // Calculate gradient at point if ðzproj ¼ zm _ km ¼ 0Þ begin // Minimum has been found (zmin, found): = (zm, true); // Store the minimum value and set flag if ðkm ¼ 0Þ then // Static interval found GetStaticInterval ðoi;k ; t l ; t u ; t m ; zl ; zu ; zm ; kl ; ku ; tmin ; zmin Þ; else t min :¼ tm ; // Minimum occurs at one point only end else begin // Continue search if ðkm > 0Þ then ðtu ; zu ; ku Þ :¼ ðtm ; zm ; km Þ; // Define right line segment else if ðkm < 0Þ then ðt l ; zl ; kl Þ :¼ ðt m ; zm ; km Þ; // Define left line segment end end End

Algorithm 2. GetStaticInterval ðoi;k ; t l ; tu ; tm ; zl ; zu ; zm ; kl ; ku ; tmin ; zmin Þ Begin if ðkl ¼ 0Þ then ta: = tl; else begin t 0l :¼ ðzm  cl Þ=kl ; // Calculate point of intersection

t a :¼ t 0a ; // Set lower limit of static interval end end if ðku ¼ 0Þ then tb: = tu; else begin t0u :¼ ðzm  cu Þ=ku ; // Calculate point of intersection 0 EvalDelay ðoi;k ; t0u ; z0u ; f u Þ; // Determine actual y axis value 0 if ðzu ¼ zm Þ then tb: = tu; else begin ðt0u ; z0u ; k0u Þ :¼ ðt u ; zu ; ku Þ; // Initialise ðt0l ; z0l ; k0l Þ :¼ ðtm ; zm ; 0Þ; // Fix left boundary FindMinPLF ðoi;k ; t0l ; t 0u ; z0l ; z0u ; k0l ; k0u ; t 0min ¼ ðt 0a ; t0b Þ; z0min ; eÞ; t b t 0b ; // Set upper limit of static interval end end tmin = (ta, tb); End

Algorithm 3. FindPtInc (oi,k, tl, tu, zl, zu, e) // A Bi Section algorithm Begin (tm, zm): = (0, 0); // Initialise mid pt and objective value for mid pt tn: = (tl + tu)/2; // Define new mid pt EvalDelay (oi,k, tn, zn, fn); // Subject the schedule to a delay while (|tn - tm|/tn > e) _ (|zn - zm|/zn > e)begin (tm, zm): = (tn, zn); // Set mid pt and objective value if (zm > zl + e) then (tu, zu): = (tm, zm); // Decrease upper boundary else if (zm 6 zl + e) then (tl, zl)(tm, zm); // Increase lower boundary tn: = (tl + tu)/2; // Define new midpoint EvalDelay (oi,k, tn, zn, fn); // Subject the schedule to a delay end if (zn > zl + e) return tn; else if (tu > tn + e) return tu; // Return point of increase End

Algorithm 4. FindPtInf (oi,k, tl, tu, fl, fu, e) // A Bi Section algorithm Begin (tm, fm): = (tu, fu); // Initialise mid pt and flag tn: = (tl + tu)/2; // Define new mid pt EvalDela y(oi,k, tn, zn, fn); // Subject the schedule to a delay while (|tn - tm|/tn > e) begin (tm, fm): = (tn, fn); // Set mid pt and feasibility status if (fm) then (tl, fl): = (tm, fm); // Increase lower boundary else (tu, fu): = (tm, fm); // Decrease upper boundary tn: = (tl + tu)/2; // Define new midpoint EvalDelay (oi,k, tn, zn, fn); // Subject the schedule to a delay end if (fn) return tn; else return tu; End

0

EvalDelay ðoi;k ; t0l ; z0l ; f l Þ; // Determine actual y axis value and feasibility if ðz0l ¼ zm Þ then ta: = tl; else begin ðt 0l ; z0l ; k0l Þ :¼ ðtl ; zl ; kl Þ; // Initialise ðt 0u ; z0u ; k0u Þ :¼ ðt m ; zm ; 0Þ; // Fix right boundary FindMinPLF   oi;k ; t 0l ; t0u ; z0l ; z0u ; k0l ; k0u ; t0min ¼ ðt 0a ; t 0b Þ; z0min ; e ;

Appendix C. Partial results of sensitivity analysis for SRT increases Below the interval times are in minutes, the slope of the profile line is denoted by grad, and the number of disjunctive graph evaluations is denoted by eval.

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

job=0,op=0,case=7,intl=(0,12.038),[12.038,12.836),[12.836,66.506],grad=1,eval=16 job=0,op=1,case=7,intl= (0,12.038),[12.038,12.838),[12.838,65.756],grad=1,eval=16 job=0,op=2,case=7, intl=(0,12.038),[12.038,12.838),[12.838,64.909],grad=1,eval=16 job=0,op=3,case=7, intl=(0,12.038),[12.038,12.839),[12.839,64.159],grad=1,eval=16 job=0,op=4,case=7, intl=(0,12.038),[12.038,12.834),[12.834,61.489],grad=1,eval=16 job=0,op=5,case=7, intl=(0,12.038),[12.038,12.834),[12.834,60.739],grad=1,eval=16 job=0,op=6,case=7, intl=(0,12.038),[12.038,12.831),[12.831,57.199],grad=1,eval=16 job=0,op=7,case=7, intl=(0,12.038),[12.038,12.840),[12.840,56.449],grad=1,eval=16 job=0,op=8,case=6, intl=(0,0.795),[0.795,36.965],grad=1,eval=17 job=0,op=9,case=6, intl=(0,0.794),[0.794,36.216],grad=1,eval=17 job=0,op=10,case=6,intl=(0,0.794),[0.794,34.041],grad=1,eval=17 job=0,op=11,case=6,intl=(0,0.883),[0.883,33.291],grad=0.9,eval=17 job=0,op=12,case=4,intl=(0,1.853),[1.853,23.421],grad=1,eval=3 job=0,op=13,case=4,intl=(0,4.120),[4.120,22.671],grad=1,eval=3 job=0,op=14,case=4,intl=(0,8.453),[8.453,20.811],grad=1,eval=3 job=0,op=15,case=4,intl=(0,8.453),[8.453,20.061],grad=1,eval=3 job=0,op=16,case=4,intl=(0,8.453),[8.453,13.363],grad=1,eval=3 job=0,op=17,case=4,intl=(0,8.453), [8.453,12.613],grad=1,eval=3 job=0,op=18,case=4,intl=(0,7.273), [7.273,8.451],grad=1.16,eval=3 job=7,op=0,case=2,intl=(0,85.999],grad=1,eval=2 job=7,op=1,case=6,intl=(0,27.857), [27.857,85.499],grad=1,eval=13 job=7,op=2,case=6,intl=(0,2.787), [2.787,84.934],grad=1,eval=17 job =7, op=3, case=6,intl=(0,2.786), [2.786,84.434],grad=1,eval=17 job=7,op=4,case=6,intl=(0,2.788), [2.788,82.654],grad=1,eval=17 job=7,op=5,case=6,intl=(0,2.786), [2.786,82.154],grad=1,eval=17 job=7,op=6,case=6,intl=(0,2.787), [2.787,79.794],grad=1,eval=17 job=7,op=7,case=6,intl=(0,2.787), [2.787,79.794],grad=1,eval=17 job=7,op=8,case=6,intl=(0,2.788),[2.788,74.329],grad=1,eval=16 job=7,op=9,case=6,intl=(0,3.098),[3.098,73.829],grad=0.91,eval=16 job=7,op=10,case=5,intl=(0,6.045), [6.045,72.379],grad=0,eval=14 job=7,op=11,case=5,intl=(0,6.726), [6.726,71.879],grad=0,eval=14 job=7,op=12,case=5,intl=(0,6.072), [6.072,69.254],grad=0,eval=14 job=7,op=13,case=5,intl=(0,6.075), [6.075,68.754],grad=0,eval=15 job=7,op=14,case=5,intl=(0,5.517), [5.517,66.962],grad=0,eval=15 job=7,op=15,case=5,intl=(0,5.517), [5.517,66.462],grad=0,eval=15 job=7,op=16,case=5,intl=(0,3.969), [3.969,60.449],grad=0,eval=15 job=7,op=17,case=5,intl=(0,4.413), [4.413,59.949],grad=0,eval=15 job=7,op=18,case=5,intl=(0,5.726), [5.726,57.174],grad=0,eval=15 job=9,op=18,case=5,intl=(0,3.117), [3.117,44.540],grad=0,eval=15 job=9,op=17,case=5,intl=(0,3.117), [3.117,44.540],grad=0,eval=15 job=9,op=16,case=5,intl=(0,2.139), [2.139,38.600],grad=0,eval=16 job=9,op=15,case=5,intl=(0,2.379), [2.379,38.100],grad=0,eval=16 job=9,op=14,case=4,intl=(0,9.805), [9.805,29.424],grad=1,eval=3 job=9,op=13,case=4,intl=(0,9.805), [9.805,28.924],grad=1,eval=3 job=9,op=12,case=4,intl=(0,1.853), [1.853,18.346],grad=1,eval=3 job=9,op=11,case=4,intl=(0,1.853),[1.853,17.846],grad=1,eval=3 job=9,op=10,case=4,intl=(0,1.853), [1.853,16.396],grad=1,eval=3 job=9,op=9,case=4,intl=(0,1.853), [1.853,15.896],grad=1,eval=3 job=9,op=8,case=2,intl=(0,9.080],grad=1,eval=2 job=9,op=7,case=2,intl=(0,8.580],grad=1,eval=2 job=9,op=6,case=2,intl=(0,6.220],grad=1,eval=2 job=9,op=5,case=2,intl=(0,5.720],grad=1,eval=2 job=9,op=4,case=2,intl=(0,3.940],grad=1,eval=2 job=9,op=3,case=2,intl=(0,3.440],grad=1,eval=2 job=9,op=2,case=2,intl=(0,2.875],grad=1,eval=2 job=9,op=1,case=2,intl=(0,2.375],grad=1,eval=2 job=9,op=0,case=2,intl=(0,0),grad=1.19,eval=2

113

114

R.L. Burdett, E. Kozan / Journal of Rail Transport Planning & Management 4 (2014) 98–114

References Burdett, R., Kozan, E., 2006. Techniques for absolute capacity determination in railways. Transp. Res. B 40 (8), 616–632. Burdett, R., Kozan, E., 2010a. Development of a disjunctive graph model and framework for constructing new train schedules. Eur. J. Oper. Res. 200 (1), 85–98. Burdett, R., Kozan, E., 2010b. A sequencing approach for train timetabling. OR Spectrum 32 (1), 163–193. Burdett, R., Kozan, E., 2014. Determining operations affected by delay in predictive train schedules. Comput. Oper. Res. 41, 150–166. Burdett, R., Kozan, E., Strickland, C. 2012. A practical approach for identifying expected solution performance and robustness in Operation Research applications. ASOR Bull.. June. Cacchiani, V., Toth, P., 2012. Nominal and robust train timetabling problems. Eur. J. Oper. Res. 219, 727–737. D’Ariano, A., Pacciarelli, D., Pranzo, M., 2007. A branch and bound algorithm for scheduling trains in a railway network. Eur. J. Oper. Res. 183, 643–657. D’Ariano, A., Pacciarelli, D., Pranzo, M., 2008. Assessment of flexible schedules in real time traffic management of a railway bottleneck. Transp. Res. C 16, 232–245. Elmaghraby, S.E., 2000. On criticality and sensitivity in activity networks. Eur. J. Oper. Res. 127, 220–238. Fischetti, M., Monaci, M., 2009. Light robustness. In: Ahuja, R.K., et al. (Eds.), Robust and Online Large-Scale Optimization, LNCS 5868. Springer-Verlag, Berlin, Heidelberg, pp. 61–84, 2009. Fischetti, M., Savagnin, D., Zanette, A., 2009. Fast approaches to improve the robustness of a railway schedule. Transp. Sci. 43 (3), 321–335.

Goverde, R.M.P., 2007. Railway timetable stability analysis using max-plus system theory. Transp. Res. B 41, 179–201. Goverde, R.M.P., 2010. A delay propagation algorithm for large-scale railway traffic networks. Transp. Res. B 18, 269–287. Kozan, E., Burdett, R., 2005. A railway capacity determination model and rail access charging methodologies. Transp. Plan. Technol. 28 (1), 27–45. Leon, V.J., Wu, S.D., Storer, R.H., 1994. Robustness measures and robust scheduling for job shops. IIE Trans. 26 (5), 32–43. Liebchen C., Lubbecke M., Mohring R., Stiller, S., 2009. The concept of recoverable robustness, linear programming recovery, and railway applications. In: Ahuja, R.K., et al. (Eds.), Robust and Online Large-Scale Optimization, LNCS 5868, Springer-Verlag, Berlin, Heidelberg, pp. 1–27, 2009. Liebchen, C., Schachtebeck, M., Schobel, A., Stiller, S., Prigge, A., 2010. Computing delay resistant railway schedules. Comput. Oper. Res. 37, 857–868. Liu, S.Q., Kozan, E., 2009. Scheduling trains as a blocking parallel-machine job shop scheduling problem. Comput. Oper. Res. 36 (10), 2840–2852. Liu, S.Q., Kozan, E., 2011. Scheduling trains with priorities: a no-wait blocking parallel machine job shop scheduling model. Transp. Sci. 45 (2), 175–198. Penz, B., Rapine, C., Tystram, D., 2001. Sensitivity analysis of scheduling algorithms. Eur. J. Oper. Res. 134, 606–615. Salido M.A., Barber F., Ingolotti L., 2008. Robustness in railway transportation scheduling. In: Proceedings of the 7th World Congress on Intelligent Control and Automation, June 25–27, Chongqing, China. Sotskov, Y., Sotskova, N.Y., Werner, F., 1997. Stability of an optimal schedule in a job shop. Omega 25 (4), 397–414. Takeuchi, Y., Tomii, N., Hirai, C., 2007. Evaluation method of robustness for train schedules. Q. Rep. RTRI 48 (4), 197–201.