informs
Vol. 42, No. 4, November 2008, pp. 450–465
issn 0041-1655  eissn 1526-5447  08  4204  0450®

doi 10.1287/trsc.1080.0246© 2008 INFORMS

Optimizing the Cargo Express Service of
Swiss Federal Railways
Alberto Ceselli

Dipartimento di Tecnologie dell’Informazione, Università degli Studi di Milano, 26013 Crema, Italy,
ceselli@dti.unimi.it

Michael Gatto

Department of Computer Science, ETH Zürich, 8092 Zürich, Switzerland, michael.gatto@inf.ethz.ch

Marco E. Lübbecke

Technische Universität Berlin, Institut für Mathematik, 10623 Berlin, Germany, m.luebbecke@math.tu-berlin.deMarc Nunkesser

Department of Computer Science, ETH Zürich, 8092 Zürich, Switzerland, marc.nunkesser@inf.ethz.ch

Heiko SchillingTechnische Universität Berlin, Institut für Mathematik, 10623 Berlin, Germany, heiko.schilling@tu-berlin.deT

he Cargo Express service of Swiss Federal Railways (SBB Cargo) offers fast overnight transportation of goods
between selected train stations in Switzerland and is operated as a hub-and-spoke system with two hubs.
We present three different models for planning the operation of this service as a whole. All models capture theunderlying optimization problem with a high level of detail: Trafﬁc routing, train routing, makeup, scheduling,
and locomotive assignment are all addressed. At the same time we respect hard constraints like tight service
time windows and train capacities, and we avoid hub overloading. We describe our approaches for obtaining
provably good quality solutions. Our algorithmic techniques involve branch-and-cut, branch-and-price, andproblem-speciﬁc exact and heuristic acceleration methods. We conclude our study with computational results
on realistic data.
Key words: railway optimization; freight transportation; hub-and-spoke systems; column generation
History: Received November 2006; revision received September 2007; accepted May 2008.

1.

Introduction

and is operated as a hub-and-spoke system (as airlines do) with two hubs. In parallel, SBB Cargo offers
a “classical” freight train system, Cargo Rail, withlooser time constraints and a larger service network.
We brieﬂy describe the Cargo Express service. Its
prominent features are the focus on fast transport,
guaranteed pickup and delivery times for customers,
and a dense railroad network spanning a comparatively small area. While the blocking mode, in whichcars are grouped and routed as a single block throughthe network (and thus are not reclassiﬁed at intermediate classiﬁcation yards), is popular in the UnitedStates, SBB Cargo Express operates in nonblocking
mode, as many European companies do. In fact, cars
are reclassiﬁed in two central hump yards, the hubs.
Each train can be composed by a limited number ofcars. Finally, hub overloading must be avoided bylimiting the number of cars in the yard at the same
time.
These features give a distinguishing ﬂavor to theplanning steps and set the Cargo Express service off

The planning of freight train operations comprises
various difﬁcult intertwined decisions on the strategic (long-term), tactical (midterm), and operational(short-term) level. Figure 1 depicts the tactical level
and the operational level according to the survey article by Cordeau, Toth, and Vigo (1998). Considerable
research has been devoted both to each single aspectof this process and toward integrating as many ofthese planning steps as possible. In general, the particular freight train system at hand determines how
difﬁcult and how important each single aspect is. For
example, there are signiﬁcant differences between theAmerican and the European systems, and even withina given country different systems with different focus
are operated.
The Cargo Express service of Swiss Federal Railways (SBB Cargo) (SBB, Cargo Express) considered inthis paper serves as an example of a freight system:
Cargo Express offers fast overnight transportation of
goods between selected train stations in Switzerland
450

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

451

Strategic

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS“Planning direction”
Freight routing

Blocking

TacticalTrain routing

Route trains
Route cars throughShipment →networkblocks (at each yard) through network
Distribute classification
work among yards

Figure 1

MakeupSchedule
Locomotive
assignment
Coarse: Frequency
Block to train
blockless: Car to of trains
Engines →Fine: Arrival and trains
train
departure of trains

Operational
Empty cars
dispatching,classification
schedules,
fitting of trains
into timetable

Mid- and Short-Term Planning Tasks in Freight Railroadingfrom most of the systems that have been consideredin the literature and the results that have been obtained for these by Assad (1980), Cordeau, Toth, andVigo (1998), and Crainic and Laporte (1997). A classical paper is that of Crainic, Ferland, and Rousseau(1984), who study the interactions between blocking, routing, and makeup for a Canadian freight system. Ahuja, Jha, and Liu (2007) give a particularly
detailed model for the blocking problem arising in
the United States, Campetella et al. (2006) consider an
Italian freight service of a size comparable to that of
SBB Cargo Express, for which they do trafﬁc routing,and planning of service frequency and empty cars,
ignoring train load limits and yard capacities. Noneof these models is intended to produce a detailed
schedule. Therefore, such models are inappropriate
for planning the SBB Cargo Express service. Kwon,
Marthland, and Sussman (1998) describe how an
existing solution consisting of blocking, routing, andmakeup plans and a target schedule can be adapted
to meet the train load constraints. Gorman (1998)
extends a model of Keaton (1992) to incorporatetime constraints on a coarse time scale. He addresses
the blocking, trafﬁc routing, makeup, and schedulingproblem for a given set of candidate routes. He proposes a tabu search approach that he tests on a U.S.
instance. In this instance he ﬁnds cost savings withrespect to a solution used in practice. However, themethod produces operating plans that can violate
constraints on the load of the trains and on time windows. In contrast, in the SBB Cargo Express setting,both time windows and load constraints are hard, andthe time windows are tighter.
In this paper, we present three different ways tomodel the SBB Cargo Express service, all of which
capture the optimization problem at hand with ahigh level of detail. Our models span the whole tactical planning process: Trafﬁc routing, train routing,makeup, scheduling, and the basic engine (locomotive) assignment are all addressed. On the other hand,
we do not model the operational level; i.e., our models do not capture empty cars movements, the precise
shunting operations in the classiﬁcation yards, minimum buffer times between arriving trains at the yard,or the ﬁtting of trains into the timetable. As we modelan overnight service, the latter aspect is not as crucial
as for other freight systems.

In the ﬁrst approach, we provide a compact integer
linear programming (ILP) formulation of the problem
and apply a state-of-the-art general purpose solver toit. In the second approach we hierarchically decompose the problem, as suggested by Figure 1, and provide separate models for distributing the classiﬁcation
work; for the combination of routing, makeup, andlocomotive assignment; and, ﬁnally, for scheduling in
our setting. We develop a branch-and-cut approach
for its hardest subproblem. Our third approach is
to formulate the entire planning problem as a sideconstrained set partitioning problem and solve itusing column generation. The integration of the aboveplanning steps leads to a complex master problem.
To the best of our knowledge, this model is the ﬁrst
to integrate all planning steps of the tactical level.
Although we tailored the solution approach to our
speciﬁc problem, the model itself is applicable withminor modiﬁcations to any nonblocking freight system of larger scale.
Our solution methods span different algorithmic
techniques generally used for ILPs. Moreover, thethree approaches have an increasing level of modelingdetail and require algorithmic techniques of increasing
complexity. We carry out a comparison of the threemethods to identify the best trade-off between modeling detail, computational tractability, and implementation effort.
All experiments are performed on data derived
from a typical day of operation of the SBB Cargo
Express service, provided by SBB Cargo.
The paper is organized as follows. First, we givea detailed problem description. Then, we introduce
our three models and next describe the corresponding solution methods. Finally, we report on computational experiments, comparing the running time of
the algorithms and the quality of the solutions found.
Problem Description. The SBB Cargo Express system is designed for customers who regularly needovernight transportation of freight by train: A typical customer orders a regular transport of shipmentsfrom an origin to a destination station. The customer
announces her ﬁxed daily demand for the lifetimeof the schedule in advance and negotiates an earliest departure time and a latest delivery time for the
shipments with SBB Cargo. Once a year, SBB CargoCeselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

452

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

designs a new operating plan and a schedule that
accommodates all customer demands.
The transport itself works as follows. In the evening, the customers deposit their cars at the departure
stations before the earliest departure times. During
the night, a ﬂeet of trains collects all deposited cars.
In general, each train of the ﬂeet collects several shipments at different stations. The cars of the same
shipment are always transported together, different
shipments at the same station may be transported
by different trains. In addition to the shunting timeneeded at the hubs, a substantial amount of timeis spent in brake testing at each stop that involves
pickup or delivery.
The SBB Cargo Express network is operated as a
hub-and-spoke system: The ﬂeet transports all the
goods to central classiﬁcation yards in Däniken and
Zürich-Mülligen, where inbound trains are reclassiﬁed to form outbound trains. Outbound trains either
go directly to the other hub nearby without en route
pickup or delivery, or they deliver their cars at
the respective delivery stations. Finally, direct trains
between origin and destination stations are also
possible.
The service is operated such as to guarantee thenegotiated earliest pickup and latest delivery times.
Because trains are composed of so few cars, eachengine of the ﬂeet operates in the same way; that is,
the planners consider the ﬂeet to be homogeneous.
Each engine of the ﬂeet can perform only one of thetasks described above in one night, i.e., either goingto and from a hub once, including potentially a few
rides between the nearby hubs, or transporting a shipment directly to its destination. An important constraint of this system is the yard capacity: Only a limited number of cars can stay in each hub concurrently,
and, most of all, hub overloading affects reclassiﬁcation time.
The main costs of an operating plan consist in operating the engines, that is, employing the drivers, andservicing the equipment. In addition to this ﬁxed cost
component, there is a cost per kilometer. The overall
goal is to minimize the total ﬁxed and variable costsof an operating plan.
A schedule is determined and operated with minor
daily changes for a whole year. Currently, the planners do not use any decision support system: They
construct the schedule by hand and slightly adapt it
in a trial period after its implementation. If demands
change over a year or new customers want to beserved, the SBB Cargo team generally succeeds in
adapting (again, by hand) the existing schedule to the
new situation.

2.

ModelsCapturing all of the characteristics of the SBB Cargo
Express system in a single model is a challenging

task. In this section we ﬁrst lay the groundwork fora precise mathematical description of the problem byintroducing the necessary common notation. Then we
present a suite of three models for our application.
The ﬁrst one, Model 0, is an ILP formulation for the
problem, neglecting hub overloading issues and forcing exactly one reclassiﬁcation step for each shipment.
In Model 1 we decompose the problem into threeconsecutive planning steps that we treat separately.
Model 2 is a set partitioning integer linear programthat exploits all of the optimization degrees of freedom by considering all of the decision levels at once.
2.1. Common Notation
We are given a railroad network N = V  E . The
node set V represents stations, hubs, and junctions;the edge set E represents the tracks connecting those; E → + is a length function on the edges. In ourmodels we consider the problem in which an arbitrary number of nodes of the network represent hubs,
and we denote the set of hubs by H ⊆ V . In the following we give a list of parameters and features that
we consider.
—The set S of shipments contains an element foreach transportation request. Each shipment s ∈ S has
the following properties:• origs: the origin station,• dests: the destination station,• departs: the earliest possible pickup time at station origs,• arrives: the latest possible delivery time at station dests,• vols: the number of cars (volume) composing s.
—The maximum train load Lmax bounds the total
number of cars that any engine can pull.
h—The shunting time at a hub Tshunt
is the additional
shunting time an outbound train has to wait before
departing from hub h ∈ H after its last shipment has
arrived. This time is assumed to be independent of
the number and volume of the shipments. Similarly,han engine needs Tshunttime units to be decoupled froman inbound train and coupled to an outbound train.
v—The couple time at the stations Tcouple
is the additional time incurred by picking up or delivering any
set of shipments at a station v ∈ V . It is independentof the number and volume of the shipments.
—The hub capacity caph speciﬁes, for each hub
h ∈ H , the maximum number of cars that can concurrently stay at h without overloading the hub.
—The engine cost Cengine represents the operating
cost of one engine (ﬁxed cost).
—The average speed v is used to calculate the traveling times of trains on the tracks.
—The kilometer cost c̄ represents the travel cost of
an engine per kilometer (variable cost).

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways
Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

The fundamental part of a solution is a route ofa train through the network. By route we mean agraph theoretic walk, which can in particular containrepeated edges and nodes, and a set of shipments
served at the visited nodes. We distinguish between
pickup and delivery routes, respectively, bringing
shipments from the stations to the hubs or from thehubs to the stations, and hub-connecting routes, whichmove cars between two hubs. Given v ∈ V , we abusenotation slightly and write v ∈ r to indicate that route r
visits node v; similarly, given s ∈ S we write s ∈ r
to indicate that shipment s is served by route r. The
length lr of a route r is the sum of the lengths of the(possibly repeated) edges used in r, and the cost cr of
the route is computed as c̄ · lr . The volume volr of aroute r is the sum of volumes of the shipments served
by r. We call r volume-admissible if it satisﬁes the load
constraint: s∈r vols = volr ≤ Lmax .
A complete solution to our problem consists ﬁrst ofsets of pickup routes Rx , delivery routes Ry , and hubconnecting routes Rh . Second, it speciﬁes arrival anddeparture times arriver v, departr v at each nodev ∈ r for these routes. We call a route together withsuch timing information a scheduled route. If a scheduled route respects all time windows and the couv
, we call it time-consistent. A scheduled
ple time Tcouple
pickup route r x to hub h is compatible with a scheduled delivery route r y from h if either r x and r y haveno common shipment, s ∈ S s ∈ r x ∧ s ∈ r y , or theyrefer to the same hub and respect the precedence,
including shunting time at the hub: departr y  h ≥h. The same notion of compatibilarriver x  h + Tshunt
ity analogously applies to connections involving hubconnecting routes. Formally, our railroad problem isdeﬁned as follows:
Deﬁnition 1 (Multihub and Spoke Optimization
Problem (MHSOP)). Given a railroad network N =V  E , a set of hubs H ⊆ V , a set of shipments S,hv, Tcouple, caph , Cengine ,and the parameters Lmax , Tshuntv, and c̄ as deﬁned above, ﬁnd a feasible solutionof minimum cost. A feasible solution consists of thesize k of the necessary engine ﬂeet and a set of pairwise compatible, scheduled, time-consistent, volume
admissible routes, such that all shipments are transported from their respective origin to their respective destination, the hub capacity limits are respected,
and the routes can be driven by k engines. The costof a solution is given by the sum of costs of theselected routes plus the cost of operating k identicalengines.
MHSOP is strongly NP-hard because it contains
problems like the traveling salesman, bin-packing,and diverse scheduling problems (Garey and Johnson
1979). This justiﬁes computationally expensive approaches like mixed integer programming.

453

2.2.

Model 0: Edge-Based Classical Vehicle
Routing ModelILP formulations have proven to be a powerful andversatile tool for modeling the whole host of NP-hardproblems. In addition, as long as standard solvers canbe applied, they represent an attractive choice froman implementation cost point of view. Therefore, we
ﬁrst evaluated the potential of such approaches. Tothand Vigo (2002) discuss various formulations for simpler but related vehicle routing problems (VRP), suchas the capacitated VRP and the VRP with time windows. Next, we adapt and extend one classical formulation therein to model MHSOP.
Let K be an upper bound on the number of engines
used in an optimal solution. Clearly, K ≤ S. We use athree-index vehicle ﬂow formulation (Toth and Vigo
2002) for MHSOP. This formulation has OE · Kinteger variables uet , each counting the number oftimes track e ∈ E is traversed by a train with engine t,
and OS · K · H  binary variables vsth , each taking
value 1 if shipment s is served over hub h by a train
with engine t. Moreover, to handle time windows,
we introduce OK · V  nonnegative variables wtv ,each representing the arrival time of engine t at station v. These sets of variables allow the formulation
with a polynomial number of constraints.
However, the use of multiple hubs introduces several complicating issues. First, pickup and delivery
routes require different sets of variables; the hubconnecting routes have to be modeled, needing OS ·K · H 2  binary variables hs t startt endt . Furthermore,delivery routes and pickup routes must be compatible, which makes further OK · S · H  binary variables necessary to model the dependency between
in- and outbound trains. Direct connections further
complicate the model. Finally, in this formulation we
do not model hub overloading issues. The complete
model is available in Ceselli et al. (2006).
The resulting formulation is “compact” in the sensethat it involves only a polynomial number of variables and constraints. The approach proved to bevery useful to obtain feedback from our partners
quickly, to gain solid understanding of the problem,
and to evaluate the optimization potential. However,
the approach is completely impractical already onvery modestly sized instances, as we will see in §4.1.
2.3.

Model 1—Hierarchical Decomposition:
Cluster-First, Route-Second, Schedule-Third
A natural approach is to follow the planning process
depicted in Figure 1 and translate it to the particular
problem at hand. For SBB Cargo Express, this means:First, partition the shipments and the engines among
the hubs. Second, route trains separately for the hubs.
Third, schedule the routes. In this section, we give one
model for each step. As with most such hierarchical approaches, it is important to robustly design the454

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

objectives and constraints of the models for the earlysteps such that they anticipate feasibility problems
that might occur in later steps.
Partition. The partitioning is performed with threeobjectives: keeping the average travel time of each
shipment short, reclassifying shipments from or to thesame station in the same hub, and balancing the classiﬁcation work among the hubs (as in the classicalblocking problem).
Because too much trafﬁc at a hub may yield infeasible scheduling problems, we treat load balancing as a
hard constraint. Then, we trade off the two remaining
objectives by computing a kilometer equivalent costfor connecting each station to a hub. The resulting
problem is a variation of a facility location problem
having a load balancing constraint.
Routing. The partitioning step creates H  singlehub subproblems that we route separately; i.e., we do
not plan dedicated hub-connecting routes. As solely
the routing, not the schedule, is decided in this step,
one cannot fully foresee the compatibility of generated routes and the feasibility for the size of theengine ﬂeet. Still, it is necessary to take these aspects
into account. To this end, we a priori ﬁx the numberof engines to a value K to guarantee that the same
number of engines enter and leave the hub. Second,
we impose a maximum distance constraint on both
pickup and delivery routes to avoid construction oflong routes. The reason for this is that routes that are
too long cannot be scheduled in a compatible way.
An advantage of this approach is that the routing
now further decomposes into routing of pickup andof delivery trains. The two routing problems are symmetric; hence, we only give a deﬁnition of the pickup
problem.
Deﬁnition 2 (Pickup Train Routing Problem
with Fixed Train Fleet (TRP)). An instance consists of a network N = V  E , a speciﬁed hub nodeh ∈ H ⊆ V , a set of shipments Sh , a maximum trainload Lmax , a maximum trip distance Dmax , the averagev
speed v, the couple time Tcouple, and the ﬂeet size K.
A feasible solution to the TRP consists of a set of K
volume admissible routes Rx , with each ending in hsuch that the following properties hold:1. For each shipment s ∈ Sh , there is a route r ∈ Rx
serving s: ∀ s ∈ Sh ∃r ∈ Rx  s ∈ r.
2. No route is longer than Dmax . The length of a
pickup route r ∈ Rx is deﬁned as the length of theroute plus the following term accounting for the couv
pling time: v ∈ r  ∃s ∈ r origs = v  · Tcouple· v.
The cost of a solution is the sum of the lengths of theroutes. The TRP asks for a solution of minimum cost.
In §3.1 we discuss a branch-and-cut approach forthis TRP.

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

Scheduling. A complete schedule for given solutions Rx for the pickup TRP and Ry for the delivery
TRP speciﬁes the departure and arrival times of eachroute at each station. However, it is not necessary to
specify a schedule in such detail: The arrival time anddeparture time windows are one-sided, in the sensethat there is a priori no latest pickup time or an earliest delivery time for the shipments. For this reason,
it is never convenient for a train to slow down on
the tracks or to wait outside a station until it is possible to enter it. Therefore, we can completely specify
a schedule by giving the arrival and departure timesof the routes at the hub and assume w.l.o.g. that thetrains travel in the fastest possible way according tos.
v and Tcouple
The objective of the scheduling step is to minimize
the maximum hub load.
Deﬁnition 3 (Train Shunting and Scheduling
Problem (TSSP)). Let the solutions Rx and Ry to thecorresponding pickup and delivery TRPs for a set S
of shipments via a hub h be given. A feasible solution
to TSSP deﬁnes an arrival time arriver x  h for eachr x ∈ Rx and a departure time departr y  h for eachr y ∈ Ry such that:
1. The (inferred) arrival and departure times ofeach route at the stations are time consistent; i.e., they
respect the time windows of the shipments.
2. All routes are compatible w.r.t. the chosen arrivaland departure times at the hub.
3. There are always enough outbound engines
available: xr ∈ Rx  arriver x  ≤ t h ∀ t ∈ + ! (1)≥ r y ∈ Ry  departr y  ≤ t + Tshunt
The cost of a solution equals the maximum numberof cars that are in the hub at the same time:




max
volr x  −volr y  ! (2)t∈+r x  arriver x ≤tr y  departr y ≤t

An optimal solution to TSSP is one for which this costis minimum.
We remark that the number of relevant constraintsin (1) is ﬁnite and that the maximum in (2) can becomputed over a ﬁnite set of relevant points in time,for example, the arrival and departure times of the
trains at the hub plus the times derived by adding or
subtracting the shunting time.
In this model we treat the capacity limit as a soft
constraint by declaring it as the objective. In fact, assketched in §1 it can be argued that the capacity of aclassiﬁcation yard is not determined by a single number of cars that, if exceeded, renders the yard inoperable but rather by a range in which it increasingly
gets more and more difﬁcult to fulﬁll the requested
classiﬁcation.

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

455

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

Model 1 has some limitations: First, as in all suchapproaches, optimization potential is lost in the hierarchical planning process. Second, on instances withtight time windows and capacity constraints, the separation between routing and scheduling may lead tofeasibility issues in the scheduling step, because therouting is oblivious to the compatibility and timeconsistency problems its routing decisions may create. A further restriction of Model 1 is that it neither
considers hub-connecting nor direct trains. These limitations are discussed in some more detail in §§4.2
and 4.4.
2.4. Model 2: Path-Based Set Partitioning Model
When tight time windows are imposed, it may be
impossible to ensure an even load at the hubs using
Model 1. To overcome this problem and exploit theadditional optimization potential of direct connections and hub connections, we present a model that
describes the MHSOP as a whole.
A scheduled route is the fundamental notion in our
problem, and it is only natural to base a model on this.
A variable in this path-based approach represents
a time-consistent, load-admissible, scheduled route,of which there are exponentially many. We select compatible pairs in such a way that a usage limit on eachhub is respected. The resulting ILP for MHSOP isgiven in Figure 2. To deal with a ﬁnite number of constraints and variables in the ﬁrst place, we consider
a set T = 0 ! ! !  " − 1 of " points in time, and wedeﬁne as time slots each interval between subsequent
points. We assume trains to always arrive or leave thehubs at the ends of a time slot. This is a common andmild restriction, as we do not have full control overthe precise times anyway.
In this model, startr and endr denote the startand end node of a scheduled route r. Each binary
variable xrt indicates whether the scheduled pickup
route r, arriving at hub endr during time slot t, isselected or not; similarly, the binary variable yrt indicates whether the scheduled delivery route r, leaving hub startr during time slot t is selected or not.
Moreover, each binary variable hatr indicates whether
the hub-connecting scheduled route r, departing fromhub startr and arriving at the destination hub endr
at time a:t, is selected or not. We may refer to one
and the same variable also by its departure time d:tfrom its origin hub. Following the policy of SBB Cargo
we do not consider pickups or deliveries on hubconnecting routes. Finally, for each shipment s ∈ S we
introduce a binary variable ds modeling the possibility of transporting s with a dedicated engine directlyfrom its origin to its destination. Such a direct path is
not associated with any time because it is always possible to deliver a shipment on the direct path, respecting the time windows.

minr∈R t∈T$s  ds +c̃r xrt + c̃r yrt + c̃r hdtr +t∈T  r∈R s∈r$̌s  ds +t∈T  r∈R s∈r

%sh 

t∈T  endr=h s∈r&hst 

s∈Sxrt = 1∀ s∈S

yrt = 1∀s ∈ S

xrt +hatr −tc̃s dst∈T  startr=h s∈r

at

xr1 + hr 1 −r s∈r endr=ht1 ≥tyrt +hdtr =0t∀ s ∈ S h ∈ H

d:t

yr2 + hr 2 ≤ 0r s∈r startr=ht2 ≥t+T h
shunt∀ t ∈ T  h ∈ H s ∈ S&ˇ hst t

d: t

yr1 + hr 1r s∈r startr=h
t1 ≤t



−

t2
r

x +h


a:t2
r

≤0

∀t ∈T  h∈H s ∈Sr s∈r endr=ht2 ≤t−T hshunt'th 

t

t1 ≤t
endr=h+ volrhr 2 ≤ caph
xrt + ha:t−rr∈R endr=ht  ≤txrt  yrt  h·tr  ds ∈ 0 1
Figure 2t

volryr2

t2 ≤tstartr=hdt

(ht 

at

volrxr1 + volrhr 1 −∀t ∈ T  h ∈ H
yrt + hd:t≥0r∀t ∈ T  h ∈ Hr∈R startr=ht  ≤t+T hshunt∀ r ∈ R t ∈ T  s ∈ S

ILP Model 2; See Comment in the Text for Interpreting theSummations

Abusing notation slightly, we denote the set of allscheduled routes by R, irrespective of their type; furthermore, all summations in Figure 2 are meant tobe over feasible routes separately for the summands.

That is, a sum of type r∈R t∈T xrt + hatr is to be readas the sum of all variables corresponding to feasible
pickup routes with arbitrary arrival time at any hub
plus the sum of all variables corresponding to feasible hub-connecting routes with arbitrary arrival timeat any hub. All sums of this type in the constraints
should be seen as a shorthand notation for two separate sums. For later use for the corresponding dual
variables, we refer to the constraints using the Greek
letters, indicated at the left of the model.
Constraints labeled $ and $̌ are set partitioning constraints stating that each shipment has to bepicked up and delivered, respectively.
The %-constraints are global ﬂow conservation constraints for the shipments at the hubs. Together with
the &-constraints they ensure the time-consistent inﬂow-outﬂow of shipments at hubs.
The & -constraints enforce that a shipment that
arrives after time t has a corresponding outboundh. The &-constraintsˇ
represent
train after time t + Tshuntthe symmetric statement for outbound trains.

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

456

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

The (-constraints play a role for the engines similar
to that the &-constraints play, except that we allow
engines to end their duty at the hub. Together, the& and ( constraints enforce the compatibility of thechosen scheduled routes.
The '-constraints limit the usage for each hub h ∈ H
to caph cars in each time slot. We remark that introducing these constraints provides an a priori guarantee to avoid hub overloading. With this model thedecision maker can evaluate different scenarios allowing for different levels of maximum hub load.
Of the three types of constraints, &,ˇ & , and %, every
pair of types implies the third type. Therefore, wechose to discard constraints &.
ˇ This implication was
expected, because & and &ˇ constraints are nothingother than a variation of the classical generalized ﬂow
conservation constraints for networks with intermediatestorage for ﬂows over time problems; see, for example,
Hall, Hippler, and Skutella (2003).
We charge the engine cost to the pickup route;
together with the ( constraints, this allows using anengine for pickup operations only.
Finally, we remark that shipments with a hub as an
origin or destination need special treatment. Although
this aspect has a minor impact on the model, it introduces particular subtleties in the solution methods,
which we do not address here.
The advantages of this formulation are that it considers both shipment and engine ﬂows in one modeland that the ﬂow conservation carries over to thefractional variables; this means that all fractionallyvalued routes in a solution have compatible counterparts, which simpliﬁes the extraction of integral solutions. Furthermore, using this decomposition, we can
efﬁciently solve the issue of time consistency andvolume admissibility in independent subproblems,as described in §3.2.

3.

Solution Methods

In this section we investigate the various algorithmic
challenges inherent in the three models: What are the(sub-)problems that arise, how do we address themboth algorithmically and implementation-wise, andwhat further techniques do we use to achieve a good
performance. Because Model 0 is meant to be optimized using a (commercial) standard solver, we focus
on Models 1 and 2.
3.1. Model 1: A Branch-and-Cut Approach
The result of the decomposition approach in Model 1
is that the MHSOP is reduced to a partitioning problem, two independent TRPs, and a TSSP for each hub.
We consider the problems in this order.
Partitioning. The partitioning problem of §2.3 canbe formulated as follows.

We introduce binary variables )sh , indicating whether shipment s is served by hub h, and binary variables
pvh and dvh , indicating whether hub h is connected tostation v by a pickup or delivery route, respectively.
connect
be the cost of connecting hub h to staLet Cvh
tion v and qsh the cost of the cheapest trip for transporting the shipment s over hub h.
The problem of partitioning the set of shipments,
such that no hub serves more than a fraction f of thetotal volume, can be formulated as follows.
minv∈V  h∈H

s.t.



h∈H

connect
Cvh·pvh +dvh +)sh ≥ 1)sh ≤ pvhs∈Sh∈H

qsh ·)sh

∀s ∈S

(3a)

∀ h ∈ H  ∀ s ∈ Sv ∈ V  origs = v-

(3b))sh ≤ dvh ∀ h ∈ H  ∀ s ∈ S v ∈ V  dests = v- (3c)vols·)sh  ≤ f · vols ∀ h ∈ H !
(3d)
s∈S

s∈S

Constraints (3a) ensure that all shipments arerouted via a hub. Constraints (3b) and (3c) enforce
consistency of hub assignments and routes: If a shipment s is routed via hub h, the origin and the destination station must be connected to h. Constraint (3d)enforces the volume balance for each hub. Here, we
connect
as the combination of a ﬁxed offset
choose Cvh
for connecting a hub to a station, plus the connection
costs in terms of distance from the hub.
This problem is NP-complete, as it contains problems like three-partition and facility location (Ausielloet al. 1998; Garey and Johnson 1979). However, thisformulation can be effectively optimized by a general
purpose solver (see §4).
Routing. Our solution approach to the TRP consists of a transformation to the following well-known
distance constrained capacitated vehicle routing problem (DCVRP), which allows us to use existing reliable
software (Ralphs et al. 2003).
Deﬁnition 4 (Distance Constrained Capacitated Vehicle Routing Problem (DCVRP)). Given acomplete network N DCVRP = V  E c , where c is acost function on edges and  is a length function onedges, a speciﬁed hub node ĥ ∈ V , demands di  i ∈ V, and a maxion the nodes, a maximum load LDCVRP
max
DCVRP
mum distance Dmax .
Find K elementary circuits starting in ĥ with minimum total cost, such that each customer node is visited by exactly one circuit, the sum of the demands, andon each circuit does not exceed the load LDCVRP
max
DCVRP
.
no circuit exceeds the length Dmax
In the following, we give a transformation / of
TRP into DCVRP, such that an optimal solution of
any TRP instance ITRP can be derived from an optimal

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

457

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

solution of the corresponding transformed DCVRP
instance / ITRP .
An optimal solution to the routing problem consists of two parts, the pickup routes and the delivery
routes. For simplicity, we again describe the transformation for the pickup case only.
Roughly speaking, the task of the transformation is
to do the following: Translate a problem deﬁned on asparse graph for which the solution consists of a setof circuits covering the network to a problem on thecomplete graph for which the solution consists of aset of paths covering the network. Moreover, we haveto correctly translate the length and load constraints.
The transformation / applies the following typesof modiﬁcations to instance ITRP to achieve the above
goals:
1. Add all missing edges to N TRP . The length ofsuch a new edge e = u v is set to the length of theshortest u v path in N TRP .
2. Increase the length of each edge of the networksby Tcouple· v.
3. Partly merge shipments with identical origin thatwill deﬁnitely be transported by the same train.
4. Replace each station with j shipments, j > 1, bya j-clique with edges of cost and length zero. Identify
each shipment of such a node with one of the nodesof the new clique by assigning the shipment’s volume
to the demand of the clique-node.
5. Add K nodes with demand M  to the graph. The
nodes are connected to the hub by edges of weight−M and to the rest of the network by the full bipartite
graph with zero-weight edges.
Steps 1 to 4 aim at making the network complete,have the couple time included in the distances, andhave one shipment for each node. The idea behind
the K extra nodes introduced in Step 5 is to force
all circuits to start or end with such a node followed
by a “free jump” to the starting node of the corresponding TRP path. Therefore, each of the K vehicles
must visit exactly one of these nodes. This is achievedby setting the weight of the hub-extra node edgesto −M M  0, their demand to M  > s∈S vols,
DCVRP
LDCVRP= M  + Lmax , and Dmax
= Dmax − M.
max
The correctness of transformation / is establishedin the following lemma, whose proof can be found in
Gatto, Jacob, and Nunkesser (2005).
Lemma 1. Let a TRP instance ITRP be given. Let &DCVRP
be an optimal solution to the DCVRP instance / ITRP 
of cost c. Then, an optimal TRP solution for the pickup
has cost c + K · M and can be reconstructed from &DCVRP
in linear time. The same statement holds for an optimal
solution for the delivery TRP.
Given an optimal solution of an instance of DCVRP,
an optimal solution of the TRP is constructed by considering the paths restricted to the nodes of the original network.

The construction can be slimmed down a bit: First,
it is sufﬁcient that the set of K nodes introduced can
induce a matching on any subset Q ⊂ V of size K.
Second, for the original graph nodes, all edges can beremoved, for which picking up both shipments of theincident nodes causes a violation of the load or length
constraints.
Solving the DCVRP. The above construction allows
us to focus on solution approaches for DCVRP.
Many exact algorithms for the DCVRP were proposed in the 1980s (Christoﬁdes, Mingozzi, and Toth
1981; Laporte, Desrochers, and Nobert 1984; Laporte,
Nobert, and Desrochers 1985). There are several software packages for the general vehicle routing problem, commercial as well as free ones; see Blasum andHochstättler (2000), Hall (2004), and Ralphs (2006)
for a survey. A reliable and successful open source
package for the VRP is the branch-and-cut code by
Ralphs et al. (2003). We base our implementation on
this package and thus describe the “two index formulation” of the DCVRP on complete undirected graphs(Toth and Vigo 2002), which we use to extend the VRP
package.
minc e xe
e∈E

s.t.



xe = 2e=i j ∈E∀ i ∈ V \ĥ

xe = 2K(4a)(4b)

e=ĥ j ∈E

e=i j ∈E i∈Q jQ

xe ≥ 2rQ∀ Q ⊂ V \ĥ  Q = 
xe ∈ 0 1

∀ e ∈ E!

(4c)
(4d)

Binary variables xe indicate whether a given edge e ∈ E
is chosen. Equations (4a) and (4b) enforce the correct
degree at the nodes and at the hub, respectively.
Constraints (4c), the capacity cut constraints, play asimilar role for the VRP as the subtour elimination
constraints do for the TSP (Lawler et al. 1985). The
left-hand side, evaluated at a solution vector, gives thenumber of edges in that solution that cross the graphtheoretic cut 5Q V \Q6. Note that every vehicle serving customers in Q contributes with two edges to thesize of the cut. The right-hand side should therefore
represent the minimum number of necessary crossings of vehicles because of the connectivity requirement, capacity reasons, and the distance constraints.
The value rQ can be understood as the maximumof two values: dQ, which accounts for the maximum distance constraints, and 7Q, which accounts
for the capacity constraints (and also for the connectivity constraints).

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

458

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

There are several valid but not equivalent choices
for a deﬁnition of dQ and 7Q. In fact, there is awhole hierarchy of possible values for 7Q that leadto different families of valid inequalities—the mostcommon being 7Q =  v∈Q dv /Lmax . This gives therounded capacity inequalities; see Naddef and Rinaldi(2002) for a more detailed discussion. The value dQis the minimum value k ∈  such that the objectivek
DCVRP
value vTSP
of a k-TSP problem on Q divided by Dmax
and rounded up equals k; see Toth and Vigo (2002):





k-TSPQ
dQ = min k ∈   k ≥
!
(5)
DCVRP
Dmax

(and, thus, circuits violating the length bound), andall subpaths resulting from integral components that
violate the length bound. As the solutions arise in abranch-and-cut setting, some of the branching decisions might not be optimal. Therefore, these cuts have
local validity only in the branching tree. To overcome
this problem, we introduce globally valid cuts as follows. We consider the graph induced by the nodes ofan integral component Qi and compute lower bounds
for the distance of the tour needed to serve thesenodes. If the lower bound exceeds the distance thatcan be served by the number of vehicles serving thecomponent, we can introduce a valid cut. We use
two easy methods to compute lower bounds for the
tours: the relaxation of TSP to 1-trees and its classical relaxation to the assignment problem; see Lawler
et al. (1985).
The transformation / and the extension of theexisting software to handle distance constraints
allowed us to effectively handle the TRP. The results
are presented in the experimental §4.

Separation Heuristics. The separation problem forrounded capacity inequalities and inequalities ofType (4c) is NP-complete. For this reason, we focus oneffective separation heuristics that try to ﬁnd violated
inequalities of Type (4c). By their heuristic naturethese methods cannot guarantee that a cut is foundeven if one exists. As the cutting plane generation isembedded into a branch-and-bound framework, thismight yield weaker bounds but does not compromise
the correctness of the algorithm.
Because the capacity constraints are effectively handled by the existing software package we use, we
focus on the distance constraints. We devised the following separation heuristics. Given a solution, we
consider the support graph induced by all nonzeroedge variables. Then, we remove the hub node andconsider each connected component Qi separately.
Because we are only interested in instances of
DCVRP arising from our transformation / , we tuned
our cuts to handle these speciﬁc instances, thus discarding the nodes introduced in Step 5 and the
involved edges of size −M in the length computation. Nevertheless, the cuts can be applied to the general DCVRP. For the ﬁrst set of cuts, we compute thelength of each component by weighing the length ofeach edge by the value of its associated variable xe .
If the resulting length exceeds what can be traveledby the vehicles serving the component, we know thatat least an additional vehicle is needed for the component, and we can enforce this with an inequalityof type (4c). We introduce these cuts for three different cases: fractional components, integral components

1

2

3

4

5

Solving the TSSP Problem. The TSSP has an interesting connection to minimum cut linear arrangement
problems. To see the connection we consider an easier variant of TSSP in which we drop Condition 3,h= 0, and only ask for the sequence of inset Tshunt
and outbound trains at the hub that minimizes thenecessary hub capacity. The resulting sequencing task
can be depicted by the bipartite in-out graph Gio =U ∪ V  E in Figure 3. The inbound trains (in Rx ) correspond to nodes of the top partition, the outboundtrains (in Ry ) to nodes of the bottom partition. Eachyedge e = rix  rj  has a volume vole corresponding toythe number of cars that train rj receives from train rix .
We call Gio a uniformly directed bipartite graph, because
all edges are directed from U to V .
The problem is equivalent to ﬁnding a linear
arrangement of the graph G, i.e., an embedding of thegraph onto the horizontal line, such that all edges are
directed from left to right. For such an arrangement,
the maximal number of edges crossing any verticalline is the cut- width, and it corresponds to the maximal number of cars residing in the shunting yard. The
width of a graph G is given by the minimal width ofa linear arrangement of G.

6
3

4

7

5

3

4
5

1

2

3

4

5

6

1

2

2

3

Figure 3
An In-Out Graph Together with a Linear Arrangement of Cutwidth 15
Note. All nonlabeled edges have unit weight.

5

7
3

1

4

5

4

6

6

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

459

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

Without the directions and the restriction to bipartite graphs, this problem is known as the minimum
cut linear arrangement, a well-studied NP-complete
problem (Garey and Johnson 1979, problem GT44)
that was shown to remain NP-hard for graphs ofdegree 3 (Makedon, Papadimitriou, and Sudborough1985) and even for planar graphs of degree 3 (Monien
and Sudborough 1988). In Gatto, Jacob, and Nunkesser(2005) the authors extend these results as follows.
Theorem 1. It is NP-hard to decide if a uniformly
directed bipartite planar graph of out-degree 3 and indegree 2 admits a linear arrangement of width .
Despite this result, we can solve the instances of the
TSSP problem that arise in our setting by using a simple ILP formulation, because these instances are not
too large. For this formulation, we discretize the time
horizon into " points in time T = 0 ! ! !  " − 1 andmake use of the previously deﬁned in-out graph Gio .
We introduce binary variables atr and drt  that model
arrival (and departure) of the trains r ∈ Rx (r  ∈ Ry ,respectively) at times t ∈ T . Here, we assume that theHshunting time Tshunt
is given in time slots. We refer to
E as the edge set of the in-out graph Gio .
min Cs.t. atr ≤ at+1r

∀r ∈ Rx  t ∈ T(6a)drt ≤ drt+1∀r ∈ Ry  t ∈ T

(6b)

t+T H
atr ≥ dr  shuntH∀t ∈ 0!!!" −Tshunt−1  ∀rr   ∈ E(6c)

H
drt  = 0 ∀r  ∈ Ry  ∀t ∈ 0!!!Tshunt−1tt
volear −dr   ≤ C ∀t ∈ T(6d)atri = 0 ∀r ∈ Rx  arriveH r > ti

(6f)

(6e)

e∈E
e=rr  tdri = 1 ∀r  ∈ Ry  depr   < tiH t  t+Tshunt
H
ar ≥
dr
∀t ∈ 0!!!" −Tshunt
−1

r∈Rx

(6g)
(6h)

r∈Ry

a0r = 0

a"−1
=1
r

∀r ∈ Rx

dr0 = 0

dr"−1
=1


∀r  ∈ Ry

all a· d· ∈ 01 !

(6i)
(6j)

Equations (6a), (6b), and (6i) impose that, for everyedge e, the variables a·e and de· form a monotonesequence starting with 0 and ending with 1.
The arrivals and departures of trains are scheduled
at the 0–1 transition of the respective variables; i.e,
for a pickup route r x ∈ Rx we set arriveH r x  = t  if
atr x+1 − atr x = 1 and symmetrically for delivery routes.
Constraints (6c) and (6d) ensure that an outboundtrain can only depart if all its cars have arrived andHthat Tshunttime units are available for shunting thosecars. Constraints (6e) represent the capacity constraint

over all time slots, which is the objective value. Constraints (6f) and (6g) introduce time constraints forthe earliest arrival/latest departure of trains; i.e., fromthe time windows we infer a constraint of typearriveH r > t  on the arrival (departure) times at the
hub and express it in the form of Constraints (6f)and (6g). Constraints (6h) ensure that a train can only
depart from the hub if there is an engine available.
Our experiments show that for the TSSP instancesthat arise from the solutions to the TRPs on the SBB
Cargo Express instance, we can calculate a schedulethat minimizes the maximum hub load and respects
the time windows in, at most, a few minutes.
3.2.

Model 2: A Column Generation-BasedApproach
Because Model 2 involves an exponential numberof variables, we adopt column generation techniques(Desaulniers, Desrosiers, and Solomon 2005): we startwith a restricted problem, containing only the direct
shipment variables ds , then solve its linear relaxation (the so-called restricted master problem, RMP),
obtaining a vector of (optimal) dual variables. We use
this dual information to identify new variables of negative reduced cost. If any such variables are found,they are included in the RMP and the whole process isiterated. Otherwise, the value of the linear relaxationof the RMP is a valid lower bound. Note that including all direct shipment variables has the advantage
of making the RMP feasible from the beginning. In
contrast to the decomposition-based Model 1, the linear relaxation of Model 2 delivers a lower bound on
the value of an optimal solution. Furthermore, good
integer solutions can be found by combining a subsetof the dynamically generated columns that satisfy theconstraints of Model 2.
Because an optimal solution in which each shipment is served by one route always exists, we cansafely relax the set partitioning constraints to set covering constraints. Then, the corresponding dual variables are restricted in sign, which leads to an easiercost structure in the generation of columns. Furthermore, it is also easier to obtain feasible RMPs.
The Pricing Problems. The problem of ﬁnding
columns of (most) negative reduced costs is called thepricing problem. For the sake of clarity, we denotethe dual variable corresponding to each constraint byadding indices to the Greek letter indicating the constraint itself. Because the pricing problems involving
pickup and delivery routes are completely symmetric,
we discuss the pickup pricing only.
Here, the reduced cost of a (column encoding a)
pickup route r to hub h scheduled at time t  is
Cengine + le · c̄ − $s + %sh  −&h ste∈r− volrt≥t s∈r t≤t s∈r'th −

t≥t (h t (7)

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

460

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

where $s  (ht ≥ 0, 'th  &hst ≤ 0, and %sh ∈  represent the dual variables associated to the constraintsof Model 2. This reduced cost has three components:a cost, which depends only on the arcs used in theroute, a prize,rcs =

s∈r


s∈r

$
s + %sh + vols ·


t≥t 

'th +t≤t &
h st which depends only on the shipments
 picked up,
and a constant contribution Cengine − t≥t (h t , which
depends only on the destination hub and arrival time
slot.
Deﬁnition 5 (Pickup Pricing Problem). Given a
network N , a hub h ∈ H , an arrival time t  , parameters Lmax , Tcouple , v, c̄, a set of shipments S, and a prizercs for each shipment s ∈ S, ﬁnd a time-consistent,
volume-admissible, scheduled route in N of minimum reduced cost.
This problem is a variation of the well-known
NP-hard resource constrained shortest path problem
with time windows (RCSPP) (Feillet et al. 2004). We
devised a particular dynamic programming shortest
path algorithm, similar to those presented in Boland,
Detheridge, and Dumitrescu (2006) and Feillet et al.
(2004), for solving the pickup pricing problem. Eachlabel l represents a partial route and encodes a 5-valued state =l >l "l ?l &l, where =l ≥ 0 is
the cost of the partial route, >l ≥ 0 the collectedprize, "l ≥ 0 the elapsed time, ?l ≥ 0 the used volume, and &l the set of visited nodes in which pickup
operations occur. Each label refers to a particular nodeof the network. Because it is not surprising for a railroad problem, our routes cannot be always elementary,that is, without node and edge repetitions. However,
we only construct routes that we call “shipment elementary,” which means that no shipment is picked up
twice. Forcing the routes to be shipment elementary
is enough to preserve the quality of the lower boundwe obtain from the RMP.
Initialization. Let S1 v and S2 v be two subsets ofpickup shipments at node v; then, S1 v dominates
S2 v if
rcs ≥
rcs ∧
vols ≤vols
s∈S1 v

s∈S2 v

s∈S1 v

s∈S2 v

and all the shipments in S1 v have no later pickup
time with respect to each shipment in S2 v. This is nota restriction in our case, because all of the shipmentsin the same station share the same earliest pickup
time. That is, picking up all shipments in S1 v insteadof all shipments in S2 v requires no more resource
consumption and gives no smaller prize. We begin
by enumerating all nondominated subsets of shipments
for each pickup station in the network. Furthermore,
we create a set L of N  labels corresponding to theinitial state 0 0 0 0 , one for each node of the
network.
Extension Procedure. Let L be the set of labels createdso far and l∗ ∈ arg minl∈L =l − >l ! In each roundthe label l∗ is pushed to every neighboring node. Letl∗ refer to node i∗ , let j be one of its adjacent nodes,
and let e be the arc connecting them: a set of labels iscreated and added to L. This set contains a label foreach nonempty, nondominated subset of shipments
Sk j at node j, encoding a state=l + le · c̄ >l +rcs "l + le · vs∈Sk j+ Tshunt  ?l +vols &l ∪ js∈Sk j

and one additional label encoding a state =l +le · c̄ >l "l + le · v ?l &l corresponding
to no pickup operation. Each label with ?l > Lmaxor "l > t  is discarded. Label l∗ is ﬁnally removed
from L and stored in a separate list L̄.
Label Pruning. As in Boland, Detheridge, andDumitrescu (2006), we delete any label l at node iwith "l + t̄ih > t  . Furthermore, let l be the value ofthe incumbent RCSPP solution. Following Lübbecke’s
ﬁndings (2005), during the creation of each label l we
compute an upper bound > on the best prize thatcan still be collected by ﬁlling the remaining volume;
this requires solving a fractional knapsack problem(Pisinger 1997). If the value =l − >l − > is still
higher than l, label l can be discarded, because it cannot yield improvements on the incumbent solution.
Dominance Rule. In our algorithm, a label l1 dominatesa label l2 , if they refer to the same node of the network,
if =l1  ≤ =l2  and =l1  − >l1  ≤ =ls  − >l2  and if"l1  ≤ "l2 , and ?l1  ≤ ?l2 . Moreover, in an optimal
MHSOP solution no route performs pickup operations
more than once at the same station, although a particular setting of the prizes rcs may yield the generation of routes containing cycles. As mentioned above,
forbidding such cycles (i.e., restricting the search toelementary routes) leads to substantially better lower
bounds. Unfortunately, this comes at the price of making the RCSPP computation much harder; in fact,
we must enforce &l1  ⊆ &l2  as a further condition forlabel l1 to dominate label l2 . We stress that our technique requires elementariness only on a small subsetof the nodes of the network: This makes the problem tractable from a computational point of view. Following Beasley and Christoﬁdes (1989), the set &lof each label is represented as a vector of binary
resources, one for each pickup station of the network;
each of them is set as consumed as soon as pickup
operations are performed at the corresponding station. Finally, as was done in the work of Irnich andCeselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways
Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

Desaulniers (2005), we tighten the dominance rule byincluding in the set &l all stations that cannot bereached anymore because of resource limitations.
Termination. Let L̄h  be the subset of labels in L̄
referring to node h . When no new label is created,a label l∗ ∈ arg minl∈L̄h  =l − >l encodes an optimal solution to the pricing problem.
Dominance of Time Slots. We also check the following simple dominance rule for entire time slots. Lethub h be ﬁxed. If t2 > t
1 and for all s ∈ S we have
(
=
0
and
&

≥
h st , that is,that
t1 ≤t<t2 h tt≤t2 h stt≤t1 &

&

=
0,
an
optimalsolutionfort
2 cannot be
t1 <t≤t2 h st
worse than an optimal solution for t1 . Therefore, we
may discard time-slot t1 in the search for the mostnegative reduced cost column.
Pricing of the Hub-Connection Routes. Because no
pickup or delivery operations occur in hub connections, these routes always follow the shortest path
between the hubs. The pricing problem simpliﬁes
to a knapsack problem, which we solve using the
MINKNAP algorithm (Pisinger 1997), adapted to handle fractional prizes as described in Ceselli and Righini(2005). Furthermore, the time-slot reduction procedures described above can be readily adapted to hubconnection pricing.
Finally, we remark that because the direct shipment
variables are included in the initial RMP we do notneed to consider them during pricing.
Acceleration Techniques. The literature on column generation abounds with acceleration techniques(Desaulniers, Desrosiers, and Solomon 2001). We
brieﬂy discuss some of our methods:
Preprocessing Strategies. We condensed the original
SBB Cargo Express network to a smaller network
without loosing optimality guarantees. The details ofsuch reductions are discussed in §4.
Heuristic Pricing. As most of the time is spent inpricing we devised heuristic speedup techniques. At
each pricing step, we stop the pricing algorithm assoon as a ﬁxed number of negative reduced costcolumns has been found. Moreover, we limit the number of generated labels, as commonly done in routing problems (Toth and Vigo 2002). We reduce thenumber of demands by aggregating shipments during pickup and delivery pricing. The heuristics are
discarded if no negative cost label is found.
Perturbation. For the calculation of each RMP we
perturb the right-hand side of the &-constraints bysmall random values. Consistent with the literature(Desaulniers et al. 1999), this yields a signiﬁcant
speedup of the LP-solving steps.
Columns Management. To keep the RMP small, we
subject the columns to aging. If a column keeps being
nonbasic for a given number of pricing iterations, it isremoved from the RMP and added to a column pool.
Before pricing, we scan this pool: If any previously

461

generated column is found with a negative reduced
cost, we insert it in the RMP and we skip pricing.
Stabilization. A common problem in column generation is that the dual variables tend to oscillate andto assume extreme values. Using an interior point
method to solve the RMP is a possible remedy, but weencountered severe numerical problems when using
the barrier algorithm in CPLEX. Instead, we adapted
the interior point stabilization approach described in
Rousseau, Gendreau, and Feillet (2007) to our problem, obtaining better convergence. However, balanced
dual values required to solve harder pricing problems
and this technique did not pay off with respect to theoverall performance.
Primal Heuristics. As solving only the linear relaxation of Model 2 to optimality takes a long time,
we experimented with heuristics based on the diveand-ﬁx paradigm (Wolsey 1998). We augmented thestandard rounding procedure by performing column
generation steps and including problem-speciﬁc rules,
obtaining a dive-and-generate heuristic. Given a fractional solution for the RMP, the heuristic iterates thefollowing steps until an integral solution has beenfound.
Step 1. Round up the column c with the highest
fractional value, given that it is consistent with thecolumns rounded up so far in the RMP constraints.
Step 2. Round down all columns in the RMP thatare not time-consistent with c.
Step 3. Solve the remaining RMP using the dual
simplex algorithm.
Step 4. If the solution is integral, stop.
Step 5. Otherwise, remove the shipments included
in c from the pickup, delivery, or hub-connection
pricing problem, according to the type of the routeencoded by c.
Step 6. Perform a ﬁxed number of pricing iterations
to include new columns into the RMP and go back to
Step 1.
Note that by rounding up and down fractional
valued columns, it can happen that the RMP may
become infeasible. If this is the case, we apply a technique that tries to restore the feasibility of the RMP
by considering optimal dual rays of the (unbounded)
dual problem. In particular, we set all the objectivefunction coefﬁcients to 0, and we simulate Phase I of
the simplex algorithm.
Because the number of pricing steps is limited,
dive-and-generate may still fail in ﬁnding a feasible
solution. Therefore, after each rounding step we execute a fast local improvement procedure that tries
to complete the current partial solution by further
rounding. To this aim, the procedure ﬁrst greedily rounds up more columns, maintaining the time
consistency between arrival and departure times ofthe shipments at the hubs. Uncovered shipments areCeselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

462

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

transported by direct paths. Inconsistencies of shipments arriving at one hub and leaving from anotherhub, as well as inconsistencies on the engines serving
each hub, are solved by using hub-connecting routes.
Hub capacity problems are addressed by shifting theroutes in time as far as possible. If capacity problems
remain, the heuristic introduces direct paths to lower
the necessary hub capacity.
We remark that our dive-and-generate approach
can be applied without problem-speciﬁc knowledge
to any ILP solved by column generation techniques.
Such heuristic is unconventional if compared to otherheuristics in column generation settings that rely
more on the compact formulation or metaheuristics
that are initialized and guided by the column generation process, see Danna and Pape (2005). Indeed, we
try to exploit, ﬁrst, the skill of the pricing algorithm ingenerating high-quality routes; second, the involved
structure of the RMP to combine the routes respecting time consistency and hub-capacity; and, third, theproblem-speciﬁc structure for local improvement.

4.

Experiments

In this section we report on the experimental resultswith the three models. All computations were carried out on a standard Linux PC with a 3-GHz processor and 4-GB memory. CPLEX 9 was used as ILP
solver for Model 0 and as LP solver for Model 1 andModel 2.
We start by describing the SBB Cargo Express
instance.
The planners of SBB Cargo Express service provided us with real data, i.e., the actual railroad network, their timetable, and the demand matrix of a
speciﬁc month. The railroad network for the SBB
Cargo Express service has 651 nodes and 1,488 edges,
and they currently operate with two hubs, locatedin Däniken and Zürich-Müllingen. In a preprocessing phase, we ﬁrst calculate the all-pairs shortest
paths among the nodes with shipments and the hubs.
Edges that do not occur on any such shortest pathcan be safely ignored. In the resulting graph we contract degree-two nodes if they are neither a hub nor
an origin nor a destination of any shipment. The
preprocessing condensed the network to 121 nodes
and 332 edges.
In Figures 4 and 5 we show, respectively, the original network and a detail of it, together with thecondensed network that we extracted. Dark nodes
represent stations with shipments in the SBB Cargo
Express service and stations without shipments thatwere retained in the condensed network. The two bigger square nodes correspond to hubs. The light colored nodes and edges are stations and connections
that are not retained in the condensed network.

Figure 4

The Original Railroad Network Together with the Condensed
Network That We Extracted

Because the given demand matrix only comprises amonth total, we divided the supply in a daily average and rounded fractional numbers. This resulted
in a 200-shipment instance. Because the actual time
windows were not available, we deﬁned the earliest
pickup and latest delivery times by relaxing the
pickup and delivery times on the currently implemented plan by one hour.
To simulate the real setting, we ﬁxed the parameters to the following realistic values: The maximum
train load is 25 cars, and the average train speed is set
to v = 60 km/hour. As a coupling time at stations wevchose Tcouple= 27 minutes. The length of each time slot
is set to 15 minutes. The shunting time at the hubs ishset to Tshunt= 27 minutes, equal at both hubs. We considered a hub to be overloaded if more than caph = 80
277:BUE
4434:EMB
5435:PFU
4 4
610:NW
5
270:W
34415:OWT 414:EL
17 3436:WWUE
6 287:STMR
3437:WTOE
2 407:WGR 93413:AD10
286:DI
1 6276:NG
6
6
24 5
7 5
4 302:OF3 4 301:BUD 275:RUEM
272:KL18 269:KE 612:KOL
7
4303:WUER
2
3300:REG 3
4271:BSD
4
613:TU
164:KLW
4
6
12
274:GLB
5
5268:EF
298:RBL
54
9273:ZSEB
223:WS8
4267:IL5
3 -1:ZMUE
297:DT
2299:SCHL
1
197:ZMUE9
12 13 222:DUE
3
4
408:SD
5 194:ZVB
195:ZIND
49 3
73221:SCWE 266:FET
4
10
10
6
3265:PFZ 409:BMA
186:BIR 645:ZGH
1567646:HOEB
2228:ZWOL
196:ZTB
220:UST
4
4
648:Z
1 3MA
1E 4
219:AAT
3 264:KN
288:BAET
187:BON 3 229:KIL
2616:KUES
218:WZ2
637
4
649:ADW
230:RUES
1
4224:HI7
188:HED3
231:TW
5
6
189:AF
3
650:LNG
5 200:HERL
3
7 232:HG199:ML
3
217:BUB
3216:RUE
7
32
213180:HGO
4
201:UET
202:MAEN
1 2
16
20 179:SBG 233:AU
3
203:STAE
9
7
234:WAE
225:RW13
5 3235:RI
4
4
12
4263:WOL
236:FRBS
4 2
226:PF
6
183:SINS 190:STE
6
262:SAM
3
227:LA
4
3
23
3
177:CHAM
5
178:ZG
261:SCHI
3643:KABO
184:OI 5
260:BIB
16
10
4
5

Figure 5

An Enlarged Detail of the Complete Network of Figure 443

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

463

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

cars are in the yard during the same time slot. Finally,
we set Cengine = 1000 and c̄ = 1. The three orders ofmagnitude difference in cost make the minimizationof distance a secondary objective over minimizationof engines.
Taking into account only the number of shipments,
the SBB instance is about the size of the largest optimally solved (much simpler) VRP instances; however,trying to solve the entire tactical planning process
turns it into a very challenging problem.
4.1. Model 0
We implemented Model 0 using the OPL Modeling language (Hentenryck 1999) and solved it using
CPLEX. On a toy instance with 11 nodes, 23 edges,
and 11 shipments, we did not get any feasible solution in seven hours. However, a feasible (not optimal)
solution for the ﬁrst ﬁve and the last six shipments
separately can be found in about 20 minutes. In fact,
this suggests a decomposition approach, and it was a
motivation to develop and implement more involved
approaches like Models 1 and 2.
4.2. Model 1
We implemented the branch-and-cut approach for
TRP, using SYMPHONY 5, a framework by Ralphs
and Ladányi (1999). We used CPLEX as LP solver for
SYMPHONY, and LEDA 4.5 for computing the minimum spanning trees and the assignment problems
needed for the cuts.
The computation ﬂow through the three steps ofModel 1 is summarized in Table 1. The table is built
by ﬁve blocks, one for each step in the solution procedure, as indicated in the heading row.
In the ﬁrst block we report the number of cars composing shipments managed by the same hub after thepartitioning step. For the partitioning we enforced abalance constraint of 55%; i.e., at most 55% of all shipments can be assigned to one hub. The resulting split
led to four instances for TRP, pickup, and delivery for
Däniken and Zürich-Müllingen.
In the second block we indicate, for each hub andfor both the pickup and delivery phases, the number of nodes N , shipments S, and engines E in
DCVRP
valthe corresponding DCVRP instance. The Dmax
ues in the DCVRP subproblems (column “Dmax ”) wereTable 1chosen by taking the maximum distance from each
hub to a shipment assigned to it and increasing itby the equivalent of time needed for shunting twice.
The rationale behind this choice is to allow serving
two shipments with the engine serving the farthest
one. The shipments for the instances were obtainedby aggregating, in cooperation with the SBB Cargo
planners, some of the demands with the same origin
for the pickup problem and with the same destination
for the delivery problem.
In the third block we report the results obtainedby solving the TRPs in terms of overall distance traveled by the ﬂeet (column “Dist.”), CPU time needed
to obtain such solution (column “CPU time”), andgap between the best known upper and lower bounds
if optimality was not proved (column “Gap”). The
engine ﬂeet size of each instance was incremented
sequentially until the number of infeasibilities withrespect to the time windows did not decrease anymore. The solution process for all instances took about
10 days in total. Because we observed dramaticallyhigher computation times for instances with fewer
engines (and more infeasibilities), the best approach is
to decrement the ﬂeet size from an initially large value
and stop when the number of infeasibilities increases.
The few remaining infeasibilities were solved by avery simple ﬁxing heuristic similar to the one presented in §3.2 for Model 2. In the fourth block of
Table 1 we indicate the overall traveled distance,
the number of used engines and the overall cost on
the whole solution. This approach resulted in using anengine ﬂeet of 11 engines for the shipments assigned
to Zürich-Mülligen (ZMUE) and 24 for the shipments
assigned to Däniken (DK), a ﬂeet of 35 engines intotal.
The resulting TSSP instance takes the time windows and the routes determined by TRP as an input.
As reported in the ﬁfth block of the table, in thetimetable that minimizes the maximum hub load,
at most 80 cars are in ZMUE and 106 cars are in DK
at the same time. This causes a slight overload in DK.
The solution was found in less than a minute.
To sum up, we were able to produce a solutionthat is close to the operational constraints. This comesat the cost of a large computation time; it is also difﬁcult to estimate losses in objective value by the threephase decomposition.

Computation Flow Through All Phases of Model 1

Partition

TRP instances

TRP results

Totals after ﬁxing

TSSP

Hub

Cars

Mode

Dmax (km)N

S

E

Dist. (km)CPU time

Gap (%)

Dist. (km)

E

Cost

Max load

DK

253

Pickupdelivery

378
296

26
31

25
30

21
21

1,660
2,590

<1 min90 min

0
0

7,686

35

42,686

106

ZMUE

208

Pickupdelivery

276
319

15
17

13
14

9
9

936
1,628

<1 min
15 h

0
5.9

80

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways

464

Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

4.3. Model 2
We implemented our column generation algorithm
for Model 2 using the SCIP library by Achterberg(2007), with CPLEX as LP solver. It is the ﬁrst time
SCIP is used to implement a CG algorithm.
Results. First, the full instance, on which we apply
Model 2, is considerably larger than the aggregatedinstances for the TRP, as it considers 200 shipments
and the preprocessed network. With our current
implementation we reach the tailing-off phase in theroot node of the branch-and-price tree after a calculation time of around four days. At this time the valueof the relaxation is 29,217.
The best integral solution found with the diveand-ﬁx heuristic of §3.2 has a cost of 35,276. Thissolution was achieved with eight intermediate pricing steps after each variable ﬁx. Similar results canbe found with more steps, but experimental results
showed that it is ineffective to price more than 12
times. The best solution was found after 85 hours
of computation. However, very similar quality results
can be found within 12 hours. The best solution
uses 27 engines, two of which are direct connections.
Twelve engines travel to and from Däniken, 13 to andfrom Zürich–Mülligen. Furthermore, there are eight
hub-connecting trains in each direction. The total traveled distance is of 8,276 km.
The local improvement is mainly useful when thedive-and-ﬁx heuristic fails or for obtaining good solutions in the early computation phases.
4.4. Comparison of Results
The advantage of Model 2 over Model 1 is twofold.
First, although the approach worked ﬁne for the hub
ZMUE, we could not avoid some overload of theslightly bigger hub DK.
Despite the shorter overall distance, Model 1 thus
requires substantially more engines than Model 2. For
the SBB Cargo Express instance, this results in a lossof 7,410 (or 20%) with respect to the best solution
obtained by Model 2. Second, the need to iterate the
TRP and TSSP solution process with an increasing
number of engines and handling the infeasibilities
leads to an enormous computation time. Thus, theintegrated approach of Model 2 performs signiﬁcantly
better on the given large real-world instance.

5.

Conclusion and Future Work

Our computational experience conﬁrms that capturing all the modeling details of multiple hub-andspoke systems while still providing reliable solution
methods is a challenging task.
However, on the other hand, we also demonstrated
that this task is by far not hopeless.
Model 1 shows both the appealing and limiting
features of distributing tactical decisions to separate
levels, as is commonly done in practice. In fact, thesolution methods of Model 1 rely on consolidated andeffective algorithmic techniques and robust existing
software packages. The optimal solution at each decision level can be found with an affordable computational effort. However, a tough parameter tuning has
to be carried out to make the algorithms of the earlier
decision levels produce feasible instances for the later
decision levels.
On the other hand, with Model 2 we are able to consider all interacting tactical decisions as a whole andprovide a solution algorithm that consistently produces feasible solutions of provably good quality.
A direct comparison between our solutions and thecurrently implemented plan is not possible at thisstage, as our experiments were carried out on average valued instances; moreover, we decided to ignore
some of the seemingly secondary aspects in our models: We do not consider the problem of engine driverassignment. We do not consider the problem of switchbacks and furcations. It is rather straightforward, however, to incorporate the handling of switchbacks andfurcations in the pricing, because these only cause aroute to take longer to be carried out.
Our column generation algorithm exploits and extends very recent optimization techniques: one issue
we did not discuss is the question of branching rules.
A natural branching rule in our setting is the assignment of shipments to hubs, based on the fractional
assignment of the solution at hand, which we implemented. However, we only did few experiments on
small instances with different branching rules andfocused on heuristically obtaining integral solutions
already in the root node for the whole instance.
A near-optimal solution of a practical problem ofthis scale and complexity was entirely out of scopesome years ago. Even though we focused on a single instance of a particular problem, we think that
research like ours will contribute to algorithmic andmodeling improvements, leading to state-of-the-art
solvers that are capable of solving a problem like ours
“out of the box.”
Acknowledgments

The authors thank Kathrin Schulte from SBB Cargo Basel forher great support in providing them with data and answering many questions about freight rail systems. Furthermore,they thank Tobias Achterberg for his support and crashcourse for SCIP and Ted Ralphs for his support with SYMPHONY. This work was partially supported by the Futureand Emerging Technologies Unit of EC (IST priority—
6th FP), under Contract FP6-021235-2 (Project ARRIVAL).
Partially supported by a Marie Curie scholarship while visiting TU Berlin.

References
Achterberg, T. 2007. Constraint integer programming. Ph.D. thesis, Technische Universität Berlin, Berlin. http://opus.kobv.de/tuberlin/volltexte/2007/1611/.

Ceselli et al.: Optimizing the Cargo Express Service of Swiss Federal Railways
Transportation Science 42(4), pp. 450–465, © 2008 INFORMS

Ahuja, R. K., K. C. Jha, J. Liu. 2007. Solving real-life railroad blocking problems. Interfaces 47(5) 404–419.
Assad, A. A. 1980. Models for rail transportation. Transportation Res.
Part A: General 14A 205–220.
Ausiello, G., P. Crescenzi, G. Gambosi, V. Kann, A. MarchettiSpaccamela, M. Protasi. 1998. Complexity and Approximation.
Springer-Verlag, Berlin.
Beasley, J., N. Christoﬁdes. 1989. An algorithm for the resource constrained shortest path problem. Networks 19 379–394.
Blasum, U., W. Hochstattler. 2000. Application of the branch
and cut method to the vehicle routing problem. Technical Report zpr2000-386, Zentrum für angewandte Informatik,
Köln, Germany.
Boland, N., J. Detheridge, I. Dumitrescu. 2006. Accelerated label setting algorithms for the elementary resource constrained shortest path problem. Oper. Res. Lett. 34(1) 58–68.
Campetella, M., G. Lulli, U. Pietropaoli, N. Ricciardi. 2006. Freight
service design for the Italian railways company. R. Jacob,
M. Müller-Hannemann, eds. ATMOS 2006—6th Workshop on
Algorithmic Methods and Models for Optim. Railways. Internationales Begegnungs- und Forschungszentrum fuer Informatik(IBFI), Schloss Dagstuhl, Germany. Retrieved January 1, 2006,
http://drops.dagstuhl.de/opus/volltexte/2006/685.
Ceselli, A., G. Righini. 2005. A branch and price algorithm for thecapacitated p-median problem. Networks 45(3) 125–142.
Ceselli, A., M. Gatto, M. E. Lübbecke, M. Nunkesser, H. Schilling.
2006. Implementation of the edge-based classical vehicle routing model and of the TSSP problem for SBB Cargo Express in
OPL. OPL ﬁles, http://www.inf.ethz.ch/personal/mnunkess/
SBB-Cargo/.
Christoﬁdes, N., A. Mingozzi, P. Toth. 1981. Exact algorithms forvehicle routing. Math. Programming 20 255–282.
Clarke, G., J. W. Wright. 1964. Scheduling of vehicles from a central
depot to a number of delivery points. Oper. Res. 12(4) 568–581.
Cordeau, J.-F., P. Toth, D. Vigo. 1998. A survey of optimization models for train routing and scheduling. Transportation Sci. 32(4)
380–404.
Crainic, T., G. Laporte. 1997. Planning models for freight transportation. Eur. J. Oper. Res. 97 409–438.
Crainic, T., J.-A. Ferland, J.-M. Rousseau. 1984. A tactical planningmodel for rail freight transportation. Transportation Sci. 18(2)165–184.
Danna, E., C. L. Pape. 2005. Branch-and-price heuristics: A casestudy on the vehicle routing problem with time windows. Column Generation. Springer-Verlag, New York, 99–129.
Desaulniers, G., J. Desrosiers, M. M. Solomon. 2001. Accelerating
strategies in column generation methods for vehicle routing
and crew scheduling problems. C. C. Ribeiro, P. Hansen, eds.
Essays and Surveys in Metaheuristics. Kluwer, Boston, 309–324.
Desaulniers, G., J. Desrosiers, M. M. Solomon, eds. 2005. Column
Generation. Springer-Verlag, New York.
Desaulniers, G., J. Desrosiers, A. Lasry, M. M. Solomon. 1999. Crewpairing for a regional carrier. Computer-Aided Transit Scheduling.
Lecture Notes in Economics and Mathematical Systems. SpringerVerlag, Berlin, 19–41.
Feillet, D., P. Dejax, M. Gendreau, C. Gueguen. 2004. An exact algorithm for the elementary shortest path problem with resource
constraints: Application to some vehicle routing problems.
Networks 44(3) 216–229.
Garey, M. R., D. S. Johnson. 1979. Computers and Intractability. W. H.
Freeman and Company, San Francisco.
Gatto, M., R. Jacob, M. Nunkesser. 2005. Optimization of a railway hub-and-spoke system: Routing and shunting. Technical
Report 477, ETH Zürich, Zürich.
Gorman, M. F. 1998. An application of genetic and tabu searches tothe freight railroad operating plan problem. Ann. Oper. Res. 78
51–69.

465

Hall, R. 2004. On the road to recovery. OR/MS Today (June),
http://www.lionhrtpub.com/orms/orms-6-04/frsurvey.html.
Hall, A., S. Hippler, M. Skutella. 2003. Multicommodity ﬂows
over time: Efﬁcient algorithms and complexity. Proc. 30th
Internat. Colloquium on Automata, Languages and Programming ICALP’03. Lecture Notes in Computer Science, Vol. 2719.
Springer-Verlag, Berlin, 397–409.
Hentenryck, P. V. 1999. The OPL Optimization Programming Language. MIT Press, Cambridge, MA.
Irnich, S., G. Desaulniers. 2005. Shortest path problems withresource constraints. Column Generation. Springer-Verlag,
New York, 33–65.
Keaton, M. H. 1992. Designing railroad operating plans: A dual
adjustment mehod for implementing Lagrangian relaxation.
Transportation Sci. 26(4) 263–279.
Kwon, O. K., C. D. Marthland, J. M. Sussman. 1998. Routing andscheduling temporal and heterogeneous freight car trafﬁc onrail networks. Transportation Res. Part E: Logist. TransportationRev. 34(2) 101–115.
Laporte, G., M. Desrochers, Y. Nobert. 1984. Two exact algorithms
for the distance-constrained vehicle routing problem. Networks
14(1) 161–172.
Laporte, G., Y. Nobert, M. Desrochers. 1985. Optimal routing undercapacity and distance restrictions. Oper. Res. 33(5) 1050–1073.
Lawler, E. L., J. K. Lenstra, A. H. G. Rinnooy Kan, D. B. Shmoys.
1985. The Travelling Salesman Problem: A Guided Tour of Combinatorial Optimization. John Wiley & Sons, Chichester, UK.
Lübbecke, M. 2005. Dual variable based fathoming in dynamic
programs for column generation. Eur. J. Oper. Res. 162(1)122–125.
Makedon, F. S., C. H. Papadimitriou, I. H. Sudborough. 1985.
Topological bandwidth. SIAM J. Algebraic Discrete Methods 6(3)418–444.
Monien, B., I. H. Sudborough. 1988. Min cut is NP-complete foredge weighted trees. Theoret. Comput. Sci. 58(1–3) 209–229.
Naddef, D., G. Rinaldi. 2002. Branch-and-cut algorithms for thecapacitated VRP. The Vehicle Routing Problem. SIAM Monographs on Discrete Mathematics and Applications. SIAM,
Philadelphia, 53–84.
Pisinger, D. 1997. A minimal algorithm for the 0–1 knapsack problem. Oper. Res. 45 758–767.
Ralphs, T. 2006. Branch cut and price references. Retrieved
January 1, 2006, http://branchandcut.org/reference.htm#
implementations.
Ralphs, T. K., L. Ladányi. 1999. SYMPHONY: A parallel framework
for branch and cut. White paper, Rice University, Houston.
Ralphs, T., L. Kopman, W. R. Pulleyblank, L. E. Trotter, Jr. 2003.
On the capacitated vehicle routing problem. Math. Programming
94(2–3) 343–359.
Rousseau, L. M., M. Gendreau, D. Feillet. 2007. Interior point stabilization for column generation. Oper. Res. Lett. 35 660–668.
SBB. 2006. Cargo Express. Retrieved January 1, 2006, http://www.
sbbcargo.com/en/index/ang_produkte/ang_produkte_express.
htm.
Toth, P., D. Vigo. 2002. An overview of vehicle routing problems.
The Vehicle Routing Problem. SIAM Monographs on Discrete
Mathematics and Applications. SIAM, Philadelphia, 1–26.
Toth, P., D. Vigo. 2002. Branch-and-bound algorithms for the capacitated VRP. The Vehicle Routing Problem. SIAM Monographs on
Discrete Mathematics and Applications, SIAM, Philadelphia,
29–51.
Toth, P., D. Vigo. 2002. The Vehicle Routing Problem. SIAM Monographs on Discrete Mathematics and Applications, SIAM,
Philadelphia.
Wolsey, L. A. 1998. Integer Programming. John Wiley & Sons,
New York.