Timetabling for Strategic Passenger Railway Planning

In research and practice, public transportation planning is executed in a series of steps, which are often divided into the strategic, the tactical, and the operational planning phase. Timetables are normally designed in the tactical phase, taking into account a given line plan, safety restrictions arising from infrastructural constraints, as well as regularity requirements and bounds on transfer times.<br><br>In this paper, however, we propose a timetabling approach that is aimed at decision making in the strategic phase of public transportation planning and to determine an outline of a timetable that is good from the passengers’ perspective. Instead of including explicit synchronization constraints between train runs (as most timetabling models do), we include the adaption time (waiting time at the origin station) in the objective function to ensure regular connections between passengers’ origins and destinations. We model the problem as a mixed integer quadratic program and linearise it. Furthermore we propose a heuristic to generate starting solutions. We illustrate the type of solutions found by our approach on two case studies based on the Dutch railway network and analyse trade-offs that are made to balance dwell times and regularity of trains.


Introduction
The public transportation planning process is traditionally subdivided into a number of steps which are assigned to either the strategic, the tactical, or the operational planning phase.
According to Huisman et al. [2005], the strategic phase encompasses a time horizon of two to ten years before implementation and includes infrastructure decisions and line planning.
The timetabling problem is often allocated to the tactical phase (approximately one year before implementation). Timetabling in the tactical phase takes a given line plan, safety restrictions arising from infrastructural constraints, as well as regularity requirements and bounds on transfer times as input.
This paper, however, focuses on strategic timetabling, i.e., the generating of a (preliminary) timetable already in the strategic planning phase. Strategic timetabling can be used to make strategic decisions with respect to timetables, like "What should the headway times be between consecutive trains at a station?" and "Where should good transfer connections between trains be made?". Due to the location of strategic timetabling early in the planning horizon, it can also be used to evaluate line plans and to point to bottlenecks in the infrastructure (and thus to promising infrastructure investments).
The value of strategic timetabling has been recognized in the practice of transportation planning. Following the example of Switzerland, the initiative Deutschland-Takt aims at establishing a so-called 'integraler Taktfahrplan' in Germany from the year 2030 on (Deutschland-TAKT [2019]). The transportation system should be redeveloped in such a way that connections between cities are served every 30 or 60 minutes, and that better transfer connections are provided. Reversing the current planning practice, the creation of a so-called 'target timetable' should precede and guide infrastructure investment decisions (e.g. in additional tracks between stations, or additional platforms at stations). In the Netherlands, a similar approach is used to evaluate infrastructure investments (Beter & Meer [2014]).
However, to the extent of our knowledge, models from academic research on timetabling as well as software tools for timetabling decision support are aimed at timetabling in the tactical (and operational) planning phase. Therefore, they focus on operational feasibility on a given infrastructure, and are rather restrictive in modeling of quality requirements. While these features are suitable for the more restrictive setting of tactical and operational planning (where changes in infrastructure and major changes in passenger behaviour are not desirable or possible), they are not appropriate to find new and innovative timetabling solutions as is desirable in the strategic planning phase.
In this paper we aim to close this gap by presenting an optimisation approach to strategic timetabling. As common in railway timetabling, we aim at finding a periodic timetable, i.e., we require that the timetable follows a repeating pattern and hence the timetable of a base period is repeated throughout the day. Our objective is to find a periodic timetable that minimizes average perceived travel time for a given line plan. Different from most other timetabling models, we include adaption time in the perceived travel time. Adaption time describes the time difference between the desired departure time of a passenger and his actual departure time. This allows us to omit regularity constraints between runs of the same line (or runs of different lines that run in parallel for part of their route), which are otherwise often used to ensure low adaption times in an indirect way, simply by enforcing regular departures.
Note that in case of dense networks and high frequencies, where OD-pairs are served by more than one line, the question which trains should be synchronized with each other becomes far from trivial to answer. In such situations, imposing synchronization constraints may lead to sub-optimal solutions or even infeasibility of the timetabling problem. This is illustrated in the following example. Consider three stations S 1 , S 2 , S 3 , and travel demand between all pairs of stations. Assume that the line plan prescribes a line from S 1 via S 2 to S 3 with a frequency of two trains per hour, and a line from S 2 to S 3 with a frequency of one train per hour. If we synchronize all trains between S 2 and S 3 , the headway between the trains on this part of the route will be 20 minutes, but from S 1 to S 2 the headways are 20 minutes and 40 minutes. This is depicted in the time-space diagram in Figure 1a, where time and distance are shown on the horizontal and vertical axis, respectively. On the other hand, if we synchronize between S 1 and S 2 , we have one headway of 30 minutes and two shorter headways between S 2 and S 3 (Figure 1b). Perfect synchronization on both parts of the network is possible, but only if one of the trains from S 1 to S 2 waits an additional 10 minutes at S 2 ( Figure 1c). Which of this solutions is best with respect to average perceived travel times depends on the size of the travel demand between the stations and the perceived value of adaption time compared to in-train time.
Our timetabling model allows us to find the best trade-off in such situations by explicitly including the adaption time into the perceived travel time, instead of deciding on where to impose regularity constraints heuristically before the optimisation. We define the Strategic Passenger-Oriented Timetabling (SPOT) problem as follows: Given a railway network consisting of stations and links connecting them, and a line plan specifying lines routes and frequencies on the network: find a timetable that minimizes average perceived passenger travel time under the assumption that passengers will choose the route with least perceived travel time. Hereby, perceived travel time is measured from the desired departure time on, that is, it includes adaption time.
Since we consider timetabling in the strategic planning phase, we cannot expect to have accurate demand information. In particular, the exact timing of travel requests is unknown.
Therefore, we think that in this time frame it is appropriate to model passengers' desired departure times as evenly spread over the period and explicitly use this assumption in our mathematical program for the SPOT problem. Both the assumption that passengers indeed arrive randomly at the station and the assumption that they adapt to the communicated timetable to a large extent can be modeled by a parameter in our objective function which relates the perceived duration of waiting at the origin to in-train time.
We model the SPOT problem as a quadratic mixed integer program that extends the traditional PESP model for periodic timetabling (Serafini and Ukovich [1989]). We linearise the model and develop a heuristic to find a starting solution. We test our approach in two case studies based on the Dutch railway network.
Note that infrastructure constraints can be included in PESP (and thus also in our SPOT model) as headway constraints in a natural way (Liebchen and Möhring [2007]). However, due to the strategic perspective we take, we do not include them in our approach for two reasons: to be able to identify promising infrastructure investments, and to keep the model tractable. In later planning phases (tactical and operational planning), when the timetabling focus shifts towards macroscopic and later microscopic feasibility, such constraints can (and should) be added. Furthermore, in our model we omit upper bounds on transfer times and regularity constraints, since this would artificially restrict the solution space and long transfers and irregular departure patterns will already be penalized in the objective function.
Our contribution in this paper is fourfold. First, we formulate the Strategic Passenger-Oriented Timetabling (SPOT) problem for timetabling in the strategic planning phase. Second, we model this problem as a quadratic integer program that integrates timetabling with passenger routing (on perceived-travel-time-minimal routes) and linearise this formulation.
Third, we propose a heuristic to construct a starting solution, in order to find good solutions even for complex large instances. Fourth, we test our model on two case studies based on the Dutch railway network, illustrating the trade-offs between the duration of dwell times and regularity of train service.
The remainder of this paper is organized as follows. In Section 2, we describe literature that is related to our study. We state our problem definition and formulate a quadratic integer programming model for SPOT in Section 3. In Section 4 we linearise this model and describe how we solve it. In Section 5 we evaluate our solution approach and perform a case study on two practical instances from Netherlands Railways. Finally, we conclude the study in Section 6.

Related Work
In this section, we give an overview on related research. In Section 2.1 we describe other attempts to timetabling in the strategic planning phase. Section 2.2 gives a brief overview on periodic timetabling. Section 2.3 describes how passenger routing can be combined with timetabling and how this is done in existing literature.

Strategic timetabling
While in the practice of public transportation planning, strategic attempts on timetabling are not uncommon (see, e.g., Deutschland-TAKT [2019]), research on railway planning normally allocates timetabling to the tactical phase in the public transportation planning process. For this reason, many timetabling models take into account a given line plan, safety restrictions arising from infrastructural constraints, as well as regularity requirements and bounds on transfer times. However, there are also approaches (see, e.g., Robenek et al. [2017Robenek et al. [ , 2016, Pätzold et al. [2017], , Borndörfer et al. [2017]) which neglect, e.g., infrastructural constraints in a railway timetabling setting, which may make them more suitable for strategic than for tactical planning.

Periodic Timetabling
The Periodic Event Scheduling Problem (PESP) model that is commonly used for periodic railway timetabling was introduced in Serafini and Ukovich [1989]. Overviews on how to model timetabling constraints and extensions in a PESP framework can be found, e.g., in Odijk [1996], Peeters [2003]. Liebchen and Möhring [2007] provide a discussion on what can be included in the PESP framework, and what cannot, like symmetry of timetables and maximum headway times between consecutive trains.
The integer program for PESP can be extended with an objective function, to find good timetables, as several of the aforementioned approaches do. For example, a weighted sum of the activity durations can be minimized. See, e.g., Peeters [2003], Nachtigall [1999], Liebchen [2008], Liebchen and Peeters [2009], Caimi et al. [2017]. More details about the PESP model are provided in Section 3.

Timetabling and passenger routing
Timetables are evaluated according to different criteria. Following Cacchiani and Toth [2012], the most prominent are that the timetable should be (1) 'efficient' and (2) 'robust'. Overviews on approaches to deal with robustness can be found in Cacchiani and Toth [2012], Lusby et al. [2018]. Efficiency can be aimed both at costs and travel time aspects (or a trade-off of both).
In the following we give an overview on how the literature addresses one aspect of 'efficiency', namely minimizing the passenger travel times, since this is also the objective we use in our model.
Early approaches to find timetables minimize passenger waiting times by assigning weights, modeling passenger numbers, to activities in the objective function (see the aforementioned references). This approach, however, neglects that passengers choose their routes based on the timetable, which makes it difficult to assign a-priori weights to activities.
Thus, better results can be obtained when timetable and passenger routing are determined simultaneously. Several approaches have been published regarding an integrated approach, both in periodic and aperiodic settings.  integrate passenger routing in aperiodic timetabling. Passenger demand is a priori assigned to a departure event, and passengers are routed from that point onwards. For periodic timetabling, a similar approach is taken by Borndörfer et al. [2017]. In this approach, train capacities are used to determine frequencies of train lines. Furthermore, many performance criteria are introduced regarding several routing methods. A more recent approach where a viewpoint on applicability in practice is taken is by . The authors study the effect of including only a subset of the OD-pairs, in order to reduce the computational effort and to obtain good timetables in a short time. Hartleb et al. [2019b] also integrate timetabling with passenger routing, but here passengers are not routed along shortest paths, but a logit distribution is used. An alternative to integrating timetabling and passenger routing in one integer programming model, is to iterate timetabling and passenger routing. Kinder [2008], Lübbe [2009], Siebert [2011], Siebert and Goerigk [2013] determine passenger flows by routing passengers through the network on, for example, shortest paths with respect to the travel time. After this, the timetable is optimised and passenger are rerouted, until a stopping criterion is reached.
The division of the public transportation planning into several sub problems (like line planning, timetabling, vehicle scheduling, etc.) is likely to lead to globally suboptimal solutions. Therefore, there are attempts to also integrate line planning and vehicle scheduling into timetabling with passenger routing (see, e.g., Schöbel [2017], Lübbecke et al. [2018]).
However, for real-life instances this leads to models that are hard to solve, and in these cases each sub problem is solved separately.
None of the aforementioned approaches considers adaption time, although a few ap-proaches exist who explicitly consider this. Some of them focus on a single corridor and schedule the trains such that the average adaption time is minimized (Barrena et al. [2014b,a]).
In these papers, a mathematical programming model and a large-neighbourhood search algorithm to find good solutions is provided. A similar situation is considered in Zhu et al. [2017], where the authors consider a bi-level model. In the upper layer, a timetable is found based on passenger demand. In a lower layer, passenger arrival times are updated such that passengers arrive shortly before their train departs, in order to minimize adaption time. Another approach to solve timetabling with integrated passenger routing including adaption times is given in Gattermann et al. [2016], where the problem is modelled as a Satisfiability problem and solved with a dedicated solver. Here passengers are assigned to a time slice and a penalty is given if a passenger cannot depart in that time slice and has to adapt to a different slice. Borndörfer et al. [2017] and Hartleb et al. [2019a] discuss and compare alternative evaluation functions for passenger-oriented timetabling. Wang et al. [2015] propose an approach to reduce the operation costs of train which models demand as time-dependent and includes route choice. Finally, Yin et al. [2017] include passenger demand and adaption time minimization into an approach to optimise energy efficiency.
In this paper, we integrate timetabling and passenger routing in a mathematical model and include the adaption time of passengers. By discarding the current infrastructural situation, solutions to our model can be used to support decision making in the strategic planning phase of railway planning.

Problem Definition and Integer Programming Model
In this section, we model the SPOT problem as a mathematical program. We introduce the necessary notation and terminology in Section 3.1. After that, we present the quadratic integer programming formulation in Section 3.2, and its linearisation.

Notation and Terminology
The approach we take in solving the SPOT problem is by formulating a quadratic mixed integer programming model integrating periodic timetabling with passenger routing. Passengers are routed along perceived-travel-time-shortest paths in the timetable. The perceived travel time is composed of adaption time, in-train time, transfer time and a possible penalty for a transfer. The penalty is indicated by a parameter γ t . Whether it is assumed that passengers have a strong preference for the requested departure time, or will adapt to the communicated timetable to a large extent, can be modelled by the parameter γ w which relates the perceived duration of adaption time to in-train time.

Periodic timetabling
We use the well-known PESP model (Serafini and Ukovich [1989]) to model the timetabling part of the problem. Following the PESP approach, we model the timetabling constraints of the problem as arcs in an event-activity network G = (V, A), with (periodic) events V and activities A. Each activity (i, j) ∈ A is a relation between events i and j, stating that the time difference between these two events should be in a given (periodic) time interval, bounded by a lower bound ij and an upper bound u ij . Examples of activities are drive, dwell and transfer activities. They respectively state a time difference between a departure and the next arrival, the time a train may dwell at a station, and a restriction on how long a transfer may take for a passenger. It is also possible to include headway activities, ensuring a certain time distance between trains. Overviews on how to model timetabling constraints in a PESP framework can be found in Odijk [1996], Peeters [2003].
The base period of a periodic timetable is denoted by T . The aim in PESP is to assign event times to all events, satisfying all these activities, i.e., to find an assignment i.e., the time difference between events i and j should be within the T -periodic interval [ ij , u ij ]. An alternative way of formulating (3.1) is by introducing a term T p ij , where p ij is an integer variable, representing the modulo operator. The result can then be written as The additionally introduced variable y ij represents the activity duration for activity (i, j) ∈ A.
Because we consider the strategic timetabling problem, we consider driving times, minimum and maximum dwell times, and minimum transfer times as input to the problem, but do not consider headway activities. We do not set upper bounds to the dwell times and transfer times, in order not to not limit the search space.

Passenger routing
Next to timetabling, we have a set of variables and constraints dealing with the passenger routing. Demand is specified in an OD-matrix OD, providing for each k ∈ OD an estimate of the number of passengers d k that want to travel from the origin to the corresponding destination per time period. We assume that the demand is assumed to be uniformly distributed over the time period.
For each OD-pair k ∈ OD, we precompute a set of possible routes. A route r ∈ R is a directed path through the Event-Activity Network. It consists of a sequence of trip and dwell activities, possibly with transfer activities, so r is an (ordered) subset of A. Note that the perceived duration of a route is timetable dependent. To precompute the possible routes, we use the method described in Warmerdam [2004]. This method first determines all direct travel options. Next, this is extended by all options with 1 transfer, then with 2 transfers, and so on. After each step, a check is done whether some travel options are dominated by others (based on expected travel time and number of transfers), and if so, they are removed.
Note that a different method would be possible as well, as long as the paths can be used as input to the model.
The set of routes for OD-pair k ∈ OD is denoted by R k . The set of all routes is denoted by R and is determined as The total (timetable-dependent) perceived duration Y r of a route is determined as the sum of all activity durations it uses, possibly plus a transfer penalty, i.e.: The function 1 t (a) is an indicator function, denoting whether activity a ∈ A is a transfer activity or not.
For each OD-pair k, we denote the set of relevant departure events, that is, all first departure events j(r) of the routes in r ∈ R k by We assume that each passenger will choose the route that minimizes his perceived passenger travel time. Note that passengers are not a priori assigned to depart at a given departure event, this is chosen together with the timetable. The route passengers choose depends on what are the next trains departing, and this can only be known once a timetable is available.
For all passengers, arriving between two potential departure events in V k , the same route will be optimal, and therefore they can all be grouped together. In order to do so, we determine the time differences between all departure events. For OD-pair k, we denote by A k v the time span between the time of a relevant departure event v ∈ V k , π v , and the departure time of the previous departure event. I.e., where α v ,v is a binary variable modeling the modulo operator. To impose an order between events, even if two departures happen at the same time, we require Then, using the assumption of uniformly distributed demand over the period, the number of passengers of OD-pair k, for who the next relevant departure event is v, is given as There may be several routes to the destination of OD-pair k starting with departure event v. We denote by Y k v := min r∈R k v Y r the time duration of the shortest route from departure event v to the destination of OD-pair k. However, starting with the next train that leaves from the origin station may not be the best option for a passenger, if the train is slow, or if the corresponding route contains a long transfer. To correctly model that each passenger takes the route with the lowest perceived travel time, we thus allow a passenger to wait at the station beyond the next relevant departure event v ∈ V k . We denote by the perceived time from π v to arrival at the destination. Note that this includes the perceived duration of the chosen route r, Y r , as well as the adaption time between event v and departure event v of route r.
Then, the perceived total travel time of a passenger whose next relevant departure event is v is computed as the sum of total adaption time until π v , and perceived travel time from π v until arrival:

Mathematical Program
Using the notation and constraints introduced above, our mathematical model for SPOT can now be formulated as follows: Constraints (3.6i)-(3.6n) state the domains of the variables.
The model stated in (3.6) is non-linear. The objective is quadratic as it contains the Next to that, (3.6e) and (3.6h) contain one or two minimums. For our computations, we linearise the objective and reformulate (3.6e) and (3.6h). For the linearisation, we define Q v,v as the periodic time difference between events v and v , i.e., (3.7) Then we can replace (3.6e) by the restrictions The first restriction represents the minimum and the second ensures that all time differences add up to T . Note that the latter is already a valid restriction in (3.6).
In order to linearise (3.6h), we introduce new binary variables z k v,v ,r , denoting if passengers for OD-pair k ∈ OD, arriving before event v use route r, starting with event v . We refer to Appendix A for details. For the linearisation of the objective, we introduce new variables For the details, we refer again to Appendix A. The model stated in this section determines a timetable that minimizes the total perceived travel time of all passengers. No synchronisation constraints are added to the model, instead, the objective is designed such that the optimal spread of trains over time is determined.

Solution approach
For real-life instances, even for networks of small size, the size and nature of the models easily exceed the capabilities of commercial solvers to find good solutions. Also due to the complex nature of the models, we do not expect to solve the models to optimality in a reasonable amount of time.
In this section we present the approach that we use in our experiments in order to find good solutions in reasonable time. First of all, we set a time limit T L. Secondly, we can simplify our SPOT model in various ways. These simplifications are described in Section 4.1.
Third, we use a heuristic method to generate a feasible starting solution, which is described in more detail in Section 4.2.

Reduced versions of SPOT
In this section we discuss some simplifications to the SPOT problem, which lead to a reduced model size and therefore possibly speed up the solution process. These can be used as heuristic approaches towards solving the full SPOT model. The first two simplifications use a subset of the OD-pairs instead of the full set. In the first simplification, we consider only passengers with a direct travel option. Note that passengers from this set do not need to travel directly if a better connection is available for them.
The second simplification is motivated by the observation that in practice, the distribution of OD-pair sizes is very skewed: only a few OD-pairs are very large and cover a large part of the passengers. We expect that the timetable largely depends on the large OD-pairs, and that including the remaining OD-pairs would only lead to some minor changes to the timetable. To choose a subset of OD-pairs, we introduce a parameter λ ∈ [0, 100]. We then include the OD-pairs which are largest in passenger size such that in total at least λ% of the passengers is included. If λ is small enough, only a few OD-pairs are included, while a large part of passengers is taken into account.
The third simplification is to require in the model that passengers always take the first relevant train leaving from the station: in that case they are not allowed to wait for a later departure. Note that also in this simplification, the order of trains departing from a station is not fixed. The intuition is that for the majority of the passengers, waiting for a later train is in general not beneficial. Therefore, we expect this to be a simplification that does not sacrifice much in terms of quality of the solution, while still reducing the complexity of the model significantly. To implement this simplification, the first minimum in (3.6h) is taken . As a fourth possible simplification, we choose to include only direct routes and do not allow for transfers. To achieve this, one could set the penalty γ t to a very large value, thus allowing transfers, but making them very expensive. We chose to reduce the sets R k v , such that it includes only direct routes. This implies that OD-pairs for who no direct travel option exists cannot be included.

Heuristic generation of a starting solution
In this section we describe a heuristic approach to solve the integer program for the SPOT problem. When trying to solve SPOT to optimality, the heuristically generated solution can be used as starting solution for the IP solver and in this way, help to speed up the solution procedure.
In order to state our approach, we first group the variables of the integer programming model for SPOT into the following sets: (4.1e) The first two sets contain variables that relate to the timetable itself. The variables in the set A are used to determine time differences between two events correctly. Finally, the sets X and Z are introduced in the linearisation of our model, and are related to the passenger routing.
To generate a good starting solution, we consecutively solve partial relaxations of the SPOT model, as outlined below. Since for all steps we require that the variables in the sets Π and P are integer, the timetabling constraints (3.6b), (3.6c), (3.6l), and (3.6m) are fulfilled.
Thus, in each step, we find a feasible timetable.
Note that as soon as a timetable is fixed, it is possible to evaluate it according to objective function (3.6a). To this end, for each OD-pair we compute the lengths Y k v of perceived-traveltime-minimal routes from each relevant departure event v ∈ V k to the destination by solving a shortest path problem. Furthermore, we order the relevant departure events, and thus compute the time difference between π v and the departure time of the previous departure, A k v , as well as the average waiting time for these passenger W k v = A k v /2. This allows us to compute the objective value of the timetable as specified in (3.6a).
We evaluate each timetable generated in the solution as described in (3.6a). The best solution according to this evaluation is stored as the incumbent and only replaced when a better solution is found.
The steps of the heuristic are detailed below. To give a quick overview, Table 1 displays for each step what type the variables are in that step, i.e., whether they are continuous (R) or integer (Z) or mixed.
Each step is solved with a time limit, that is based on the overall time limit T L. Furthermore, the solution for each step is used as a warm start for the next step. The heuristic is a variant on the 'Relax-and-Fix' heuristic, as explained in Belvaux et al. [1998], Wolsey [1998].
1. In the first step, a solution is found that is feasible with respect to all timetabling constraints. Therefore, we relax all variables to continuous variables, except for the timetabling related variables, i.e., those in Π and P. This model is solved to optimality, Step Π Fixed Fixed Fixed Z Z 0.0 Table 1: Overview of the integrality of variables in the heuristic or until a time limit of T L/10 is reached.
2. In this step, we improve the time differences between trains to get a better passenger routing, by changing a subset A of the variables in A into integers, which we initialize as A = ∅.
In order to determine this set, we check for each pair of trains t 1 and t 2 whether their geographical routes overlap. If so, let v 1 , v 2 ∈ Π be the departure events of trains t 1 and t 2 , respectively, at their first shared station.
We change all variables in A into integers and set A = A \ A and A = ∅. Then we re-optimise with a time limit of T L/10 or until an optimality gap of 1.0% is reached.
3. In the previous step, a part of the α variables is changed into integers, but the majority is still continuous. In this step we iteratively change the remaining variables in A into integers, using a 'most-fractional' rule, i.e., we start by changing these variables of which the value in the incumbent solution is closest to 0.5.
In each iteration of this step, we define the set of variables that are to be changed to integers as Here, val(α) denotes the value this variable α attains in the incumbent solution.
Again, we change the nature of all variables in A to integers, we set A = A \ A , A = ∅ and re-optimise with a time limit of T L/10 or until an optimality gap of 1.0% is reached. This is continued until |A | ≤ 50, in which case we set A = A in order to limit the number of iterations. Furthermore, this ensures that after these loops all α-variables have integer values.
4. In this step, we fix all variables in Π, P and A to the value they attain in the incumbent solution (according to the evaluation with all OD-pairs). Next, we change all variables in X and Z to integers and reoptimise this model to optimality.
In order to better understand the heuristic, we highlight the rationale behind it. As headway constraints are not considered, the timetabling part is relatively easy in our model.
Therefore we first find a timetable that is feasible with respect to the timetabling constraints, and include the passenger routing part only as a continuous relaxation. As this is a relatively easy task, we try to find an optimal solution here. This can however lead to a bad timetable, as the time differences between events can be determined incorrectly, due to the continuous nature of the variables in A. Therefore, in the next steps we try to improve this.
First, we determine where train lines meet for the first time. By making the corresponding variables integer, we aim at better spreading different train lines over time. The expectation here is that by changing only a few variables to integers, a good gain in terms of quality can be obtained, without making it too difficult. The places where train lines meet are these places where frequencies on the tracks can change and therefore the expectation is that these are crucial decisions. The next step turns the remaining variables into integers.
By experiments we found that the majority of the variables in A is close to integer, and that the remaining variables are generally rather close to 0.5. Therefore we select these variables that have 0.05 ≤ val(α) ≤ 0.95. Iteratively these variables are changed to integers. When there are not many variables left, we change the remaining variables into integers in order to limit the number of iterations needed. Finally, for the best found timetable, we determine the best routing and the heuristic finishes.

Computational Results
In this section, we apply our approach to two instances based on the network operated by Netherlands Railways which we describe in Section 5.1. We use these instances to computationally evaluate the use of a heuristic starting solution in Section 5.2.1 and to investigate the effect of solving restricted versions of our model in Section 5.2.2. In the case studies in Section 5.3 we look more in detail into the solutions created by the SPOT model, discuss our findings, giving some insights on how our approach can be used in strategic railway planning.
In all experiments we discretise time to minutes and use a period length of one hour, i.e., T = 60.
For the perceived travel time, values for adaption time and transfer penalty have to be set (γ w and γ t ). According to De Keizer et al. [2015], the resistance for a transfer depends on many factors, but on average a penalty of 23 minutes (including 2 minutes of transfer time) is appropriate. We use a minimum of 3 minutes for a transfer, so we chose to use a value of γ t = 20 in our models. We want to put emphasis on the regularity of trains to reduce adaption time, but not over-stress it because it already appears in the objective as a quadratic term, so we use γ w = 3.
In our implementation and when reporting objective values in this section, we only report the 'additional time'. That means, we subtract constant terms from the objective function to improve numerical stability. For trip time, that implies that we subtract the minimal duration of the shortest possible trip for that OD-pair. If some OD-pair needs at least one transfer, we subtract a penalty, and only penalize additional transfers. Finally, for the constant for the adaption time, we assume that departure events are spread evenly over time, which leads to the lowest possible adaption time, and subtract the corresponding value for the adaption time. This also explains why we do not report relative gaps. If we do not subtract the constant terms, all gaps would be relatively small. In our experiments, the lower bound is often close to zero and hence relative gaps are very large.
Our computations are carried out on a machine with an Intel Xeon Silver 4110 2.10Ghz processor and with 96 GB of RAM installed. The integer programs are solved by Cplex 12.9.0 under default settings, using up to 15 parallel threads.

Instances
The instances that we use in this study are real-life instances of Netherlands Railways (NS), the largest operator of passenger trains in the Netherlands. The first instance is the so-called 'A2-corridor', a network that contains 5 Intercity lines, that all share part of their route. The second instance is the 2019 Intercity network of Netherlands Railways (NS). In the remainder of this section, we describe the two instances in more detail.

A2 corridor
The first instance we consider in our study contains the so called 'A2-corridor', which is the part of the Dutch railway network between Eindhoven (Ehv) and Amsterdam Centraal (Ehv).
The line plan for this instance is shown in Figure 2a. The used abbreviations for the stations are mentioned mentioned in Table 2.   The reason to study this instance is that the 'A2-corridor' has very high passenger numbers and it has been subject to intensive study in practice recently, since Intercityfrequencies increased from four to six trains per hour here. In Asd and Ehv, four of the six trains continue to Amr and Std, respectively. This raises the question what the headway times should be between consecutive trains, both on the corridor itself and on the remainder of the network. As an example, if the headway times between all consecutive trains on the corridor is 10 minutes upon arrival in Ehv, and trains do not get additional dwell time there, the pattern between Ehv and Std will be irregular, headway times alternate between 10 and 20 minutes. In order to get a more regular pattern, trains would have to get a longer dwell time in Ehv. We study these kind of situations to find out what is the best solution from a perspective of total perceived passenger travel time.
In this instance, we consider only OD-pairs that travel either in the southbound or the northbound direction, and not in both directions. For example, OD-pairs Nm to Ehv and vice versa are not considered, as they would have to travel via Ut. In total the network has 34 relevant stations and we consider 891 OD-pairs. The average number of routes per OD-pair is 6.3, with a maximum of 24 routes. The Event-Activity Network contains 1344 events and 1700 activities, of which 376 are transfer activities.

Intercity network
As second instance, we consider the 2019 Intercity Network of NS. In this network, there are many OD-pairs without a direct connection. We thus expect that arrival and departure times at important transfer stations will be influenced by the need to make good transfer connections for these passengers. The network includes 95 stations and 76 trains in total.
The geographical network is depicted in Figure 2b. There are 8870 OD-pairs.

Evaluation of the solution approach
In this section, we evaluate our solution approach. In Section 5.2.1 we evaluate computationally the benefit of generating a heuristic starting solution instead of a cold start. In Section 5.2.2 we explore the effects of solving several restricted versions of the SPOT model on the quality of the timetable.

The benefit of generating a starting solution
To evaluate the benefits of using a starting solution, we compare running times of the linearised SPOT model, with and without starting solution. We do this on three different cases: the A2-corridor instances, the intercity network instance with OD-pairs which have a direct travel option, and the intercity network instance with all OD-pairs.
For generating a starting solution, we set λ = 30, i.e., at least 30% of the passengers in the network are included. Given the distribution of the OD-pair sizes (only a few OD-pairs account for a large portion of the passengers) and after performing several tests, this turned out to be a reasonable number to use for this purpose. This leads to including only a small subset of the OD-pairs in the model, while still ensuring that a large portion of the passengers is covered. We then employ the heuristic procedure described in Section 4.2.
To guide the search when no heuristic starting solution is generated, we first spend 20% of the allotted time on a model where all dwell times are set to their lower bound. The remaining 80% of the time is spent on solving the full model.
For the A2-corridor, we set a total time limit of 2 hour for the computations. For the Intercity case, we set a total time limit of 10 hours.
The results of our computations are shown by means of convergence plots in Figure 3.
The horizontal axes display time in seconds on a logarithmic scale. Note that the heuristic solves a strongly restricted problem with a subset of the passengers, and therefore the objective values of the heuristic and the objective value of the full model cannot be compared.
Therefore, every time a new timetable is found in the solving process, its objective value (3.6a) is evaluated based on the full set of OD-pairs, in the way described in Section 4.2. To indicate how much time is consumed for the different steps of the heuristic, we indicate the time taken by the different steps in Figure 3 by means of shaded bars. Each bar displays the step of the heuristic as well. As the third step iteratively changes a subset of the α variables, we also display the iteration number, i.e., 3-2 denotes step 3, iteration 2. Only a few iterations are needed to turn all α variables into integers.
Step 4 is not displayed because this interval is too short to be visible on the logarithmic time scale we use. More details on how much time is spent on each step can be found in Table 3. This table also displays the number of OD-pairs used for the heuristic, and its percentage of the total number of OD-pairs in the instance.  We see in Figure 3 and Table 4 that in all three cases, better solutions are found when generating a heuristic starting solution first. Even stronger, the heuristic finds a good solution, before the full model finds the first feasible solution.

Solving reduced versions of SPOT
In the previous section, we motivated the use of a two-stage approach: generate a heuristic start solution, then use this to warm-start the linearised SPOT model. In this section, we follow the two-stage approach, and experiment with solving different reduced versions of SPOT in the second stage. The rationale behind this is that on the one hand, within the same time limit, we may be able to get closer to optimality when working on a reduced   Table 4: Comparison between using a heuristic starting solution or not problem version. On the other hand, we hope that when reducing in the 'right' way, little relevant information is lost, such that the timetable we find is good when evaluated for the full problem.
To reduce the model, we experiment with four different parameters: (1) We take only passenger who have a direct travel option (only Direct), or all passengers (all OD); (2) We take λ ∈ {95, 99, 100}, i.e., we restrict the number of passengers that we take into account; (3) We either allow passengers to transfer (trans) or not (noTrans); (4) We either force passengers to take the first departing train (noWait), or let them wait for a later train (wait). See also Section 4.1 for detailed explanations of the simplifications.
Note that we do take these parameters already into account when constructing the heuristic starting solution. E.g., when we do not consider OD-pairs with transfer options (only Direct), only these are considered in the heuristic.
In the following we compare the evaluation values for our approach under different pa-rameter setting on the Intercity network instance. We chose this instance because this is the most difficult instance to find good solutions for. In order to properly compare the resulting timetables, we have evaluated each of them considering all passengers with full route choice.
The corresponding evaluation values are displayed in Figure 4.  Note that when including all OD pairs, but not allowing transfers, this leads to many OD-pairs not having a transfer option in the model. Therefore, we leave out these situations.
As can be seen in this figure, no single parameter choice seems to lead to clearly superior results. As a tendency, in this instance it appears that the combination of 'noWait' and 'Trans', i.e., forcing passengers to take the first train that is leaving, and allowing transfers, leads to lower objective values. A possible explanation for this is that the cases in which waiting at the origin station is beneficial will be very limited, especially since we are looking at an Intercity network in which trains run at the same speed. The option 'noWait' hence leads to a smaller model without sacrificing much in terms of quality. Next to this, transfer options can significantly improve route choices and providing them seems to be relevant. As for the passenger sets to include in the model, it seems that excluding the smallest OD-pairs can indeed lead to better results, since they will have minor influence on the overall quality of the timetable. However, the difference between the solutions found with the different parameter settings is very small, and we have not tested the different parameter settings on different instances, therefore we have to be very careful in our conclusion.

Case studies
In this section, we look into the solutions that our approach finds in more detail. In Section 5.3.1 we focus on the A2-corridor and demonstrate the trade-offs that our model is able to make, between additional dwell times and regularity of train service. In Section 5.3.2 we demonstrate how our method makes choices regarding transfer connections between trains on the intercity network.

Balancing regularity and dwell times
We illustrate the outcome of our approach on the A2-corridor instance in more detail. Specifically, we focus on what happens at the two locations Asd and Ehv where frequencies change from six to four trains per hour. Remember that, when passing from the corridor with frequency six to the part pf the network with frequency four, if the trains arrive with 10 minutes headway time between each pair of consecutive trains and no additional dwell time at the station is allowed, the headway times will be irregular outside the corridor, alternating between 10 and 20 minutes. If we want the patterns to be regular both on and off the corridor, additional dwell times are required at Asd and Ehv, in order to make the transition between the different frequencies. In order to shed light on the trade-offs at different values of adaption time, we visualize two timetables for the A2-corridor. To find the first one, we ran our solution approach for different values of parameter γ w , which relates the perceived duration of adaption time to in-train time. In the first situation γ w = 3, thus adaption time is considered to be less pleasant than being in the train itself. This will put a higher emphasis on the regularity of trains. In the second situation γ w = 1, thus adaption time is valued equal to in-train time.
Time space diagrams of the timetables we found for both situations are shown in Figure 5.
Time is shown on the horizontal axis, distance on the vertical axis, where also the relevant stations are shown. Each train line is plotted with a different color. More detailed timetables for the stations Asd and Ehv are shown in Tables 5 and 6. Interesting to note is that the arrival pattern in Asd with γ w = 1 is perfectly regular, all headway times between consecutive trains are exactly 10 minutes. However, when continuing towards Amr, the pattern is becomes very irregular, as there is no additional dwell time added, and the headway times now alternate between 10 and 20 minutes. Also in Figure 5b, this irregularity is clearly visible. With γ w = 3, the irregularity north of Asd is reduced to headway times of 13 and 17 minutes. In this case the headway times on the corridor are no longer equal and vary between 8 and 13 minutes.
The arrival headways at Ehv are fairly regular, they vary between 9 and 11 minutes, for both values of adaption time. Also here, departure headways are not perfectly regular.
Instead, for γ w = 3, they are even as large as 21 minutes, which is larger than the departure γ w = 3   Interesting to note in Figure

Insights on the Intercity network
The Intercity network of Netherlands Railways contains many train lines that are linked to each other, because they share part of their route and because they can offer good connections.
We have analysed the timetable that was evaluated best in Section 5.2.2 (that is, the one found for parameter setting 'only Direct', λ = 95%, 'trans', 'noWait') to demonstrate how our method can be used to generate insights on the desired timetable. For the purpose of this paper, we analyse the timetable at two different stations, namely Leiden (Ledn) and Zwolle (Zl).    otherwise. In the timetable at Ledn, we have not seen these prolonged dwell times, as this would be less beneficial there due to the higher frequencies.

Leiden (Ledn)
The above discussion indicates how SPOT can be used to design timetables, and give valuable insights in the strategic timetabling phase. We observe that in Ledn and Zl, the solution is often close to regular, but that exceptions from these patterns can improve the timetables in some cases. This implies that we may overlook good timetables, when imposing regularity constraints. Especially in Ledn, we see that these irregular patterns allow the alternating order of trains, and in turn the alternating connections between trains.
Secondly, we see that for some stations longer dwell times are good to ensure transfers, especially if the alternative for missing the transfer would be a long waiting time, when frequencies are low.

Conclusions and further research
In this paper we introduce the Strategic Passenger Oriented Timetabling (SPOT) problem.
This problem aims at finding a timetable pattern which is optimal for passengers, explicitly including adaption time into the perceived passenger travel time. In our approach to solve the SPOT problem, we formulate a quadratic integer program. We linearise it and we propose and test an approach for solving it. We have shown in our case studies, how the solutions generated by the SPOT model can be used to learn about desirable patterns at key points of the network.
Due to the strategic nature of the problem at hand, we formulate the SPOT problem without including headway constraints, so that the underlying timetabling problem is relatively simple. However, the inclusion of adaption time in the model formulation leads to a quadratic objective, making the model harder to solve again. We achieve improvement with respect to the solution time by warm-starting the model with a heuristically achieved solution. Still, in none of our instances we were able to prove optimality of the solution found, with a lower bound far off the best solution found. It may be promising to investigate further solution methods, possibly working directly on the quadratic formulation of the program.
Most timetabling models for the tactical planning phase do not include adaption time into the perceived travel time, thus implicitly assuming that passengers will fully adapt to the timetable and not suffer from inconveniently placed departure times. To overcome this questionable assumption, the SPOT model could also be applied to a setting that includes headway constraints in the PESP model. However, we expect that these will make the model even harder to solve. In addition, our current model formulation uses the assumption that passengers arrive uniformly distributed over the period for the definition of passenger groups and average waiting times. To include more detailed passenger information that may become available when entering the tactical or operational planning phase, different modeling approaches may be needed.
A different idea on how to move towards the tactical planning phase is to consider the timetable obtained with SPOT as an ideal timetable, and adjust it, where needed, to make it 'fit' with infrastructure requirements. A method to do this is currently under research.

A Linearisation
The SPOT model in (3.6) contains a quadratic objective and has several minimums in the constraints. In Section 3.2, constraints (3.6e) is linearised. In this section, we linearise the remainder of the model.

A.1 Objective
The objective function (3.6a) is a quadratic function. Using that W k v = A k v /2, the objective can be written as We linearise this expression by writing A k v as a sum of binary variables, defined as For a stronger formulation, we can impose the additional restrictions that Using these new variables, we can write Substituting this in (A.1) results in a multiplication of binary variables x k v,d by bounded variables Y k v . This can be resolved by introducing new variables R k v,d = Y k v · x k v,d . The objective then becomes to minimize and additional restrictions have to be added to correctly determine the value for R k v,d : where l k v and u k v are the lowest and highest possible values for Y k v respectively.

A.2 Minimums
Constraints (3.6e) and (3.6h) both contain a minimum, which we can linearise. (3.6e) is already linearised in Section 3.2. Constraints (3.6h) are replaced by the following set of restrictions for every k ∈ OD and every v ∈ V k : We introduced new binary variables z k v,v ,r , which correspond to the route that is chosen. That means, if z k v,v ,r = 1, passengers wait from event v to v (which can be the same), and take route r, starting at v . For computational stability, the newly introduced constants M k v have to be chosen as small as possible, but still large enough to make the second set of constraints redundant if z v,v ,r = 0, i.e., we can take where Y r , Y r denote the highest and lowest possible value for the variable Y r respectively.
As we are minimizing the perceived passenger travel time, we can exclude the newly introduced constraints (A.6a), (A.6c) and (A.7a) in order to reduce the model size.
To summarize, in the linearisation we take several steps. First of all, the objective (3.6a) is replaced by (A.5). Here, additional variables x k v,d and R k v,d are introduced, with additional restrictions (A.3), (A.6b) and (A.6d). Secondly, the minima are replaced by linear restrictions.

B List of symbols
This appendix summarizes all notation used in the paper.

V
The set of events (indexed by v, i or j) A The set of activities (indexed by a or (i, j)) OD The set of all OD-pairs (indexed by k) R The set of all routes (indexed by r) R k The set of routes for OD-pair k (indexed by r) V k The set of departure events for OD-pair K (indexed by v or v ) Constants T The cycle period d k The number of passengers for OD-pair k ∈ OD γ w The objective coefficient for adaption time γ t The penalty for using a transfer Variables π v The timing of event v ∈ V p ij A modulo parameter used for the shift from one cycle period to another, for activity (i, j) ∈ A y ij The duration of activity (i, j) ∈ A Y r The duration of route r ∈ R A k v The number of minutes before event v, in which no other departure event for OD-pair k takes place W k v The expected waiting time for passenger for OD-pair k, for who event v is the first departure event Y k v The perceived travel time for passengers of OD-pair k, from the timing of event v onwards α v,v An integer variable ensuring the correct determination of the time difference between event v and v Q v,v The time difference between event v and v