% Template article for preprint document class `elsart'
% SP 2001/01/05

% Modified 2007/01/09 for CTW 2010
% Please find the original file at http://www.elsevier.com/

\documentclass{ctwart}

% Use the option doublespacing or reviewcopy to obtain double line spacing
% \documentclass[doublespacing]{elsart}

% if you use PostScript figures in your article
% use the graphics package for simple commands
% \usepackage{graphics}
% or use the graphicx package for more complicated commands
% \usepackage{graphicx}
% or use the epsfig package if you prefer to use the old commands
% \usepackage{epsfig}

% The amssymb package provides various useful mathematical symbols
\usepackage{amssymb}

\begin{document}

\begin{frontmatter}

% Title, authors and addresses

% use the thanksref command within \title, \author or \address for footnotes;
% use the corauthref command within \author for corresponding author footnotes;
% use the ead command for the email address,
% and the form \ead[url] for the home page:
% \title{Title\thanksref{label1}}
% \thanks[label1]{}
% \author{Name\corauthref{cor1}\thanksref{label2}}
% \ead{email address}
% \ead[url]{home page}
% \thanks[label2]{}
% \corauth[cor1]{}
% \address{Address\thanksref{label3}}
% \thanks[label3]{}

\title{A Heuristic Algorithm for the Train-Unit Assignment Problem}

% use optional labels to link authors explicitly to addresses:
% \author[label1,label2]{}
% \address[label1]{}
% \address[label2]{}

%\author{}
\author[]{V. Cacchiani\corauthref{cor1}}
\ead{valentina.cacchiani@unibo.it}
\author[]{A. Caprara}
\ead{alberto.caprara@unibo.it}
\author[]{P. Toth}
\ead{paolo.toth@unibo.it}

\corauth[cor1]{V. Cacchiani}

%\address{}
\address{DEIS, University of Bologna, Viale Risorgimento 2, 40136 Bologna, Italy}

\begin{keyword}
% keywords here, in the form: keyword \sep keyword
Train Unit Assignment \sep Lagrangian relaxation \sep Heuristic
\sep Assignment Problem
\end{keyword}

\end{frontmatter}

% main text
\section{The Train-Unit Assignment Problem} \label{probldescr}

\vspace*{-.5cm}

We study the Train-Unit Assignment Problem (TUAP) (see e.g.
\cite{Cetal05}), calling for an optimal assignment of train units
(which are self-contained trains with an engine and passenger
seats) to a given set of timetabled train trips. More precisely,
each trip has a departure station, an arrival station, a departure
time and an arrival time, and requires a number of passenger
seats. Train units can be classified in different types: each
train unit has a number of available seats and can be combined
with other units in order to fulfill the seat requests. For each
train trip a maximum number of train units that can be combined is
given. In addition, sequencing constraints between trips must be
satisfied: a pair of trips can be in a sequence for a train unit
if the time elapsing between the arrival of the first one and the
departure of the second one is large enough to allow the
corresponding train unit to travel from the arrival station of the
first one to the departure station of the second one. Finally,
each train unit has to undergo a maintenance operation every fixed
number of days, which requires a certain amount of time, as well
as the transfer to and from the maintenance station. The goal is
to minimize the number of train units globally used, while
satisfying the seat requests, the sequencing constraints and the
maintenance constraints.

\section{Solution method} \label{sol}

\vspace*{-.5cm}

%\subsection{Arc formulation and Lagrangian relaxation} \label{arcILP}
Let $n$ be the number of trips and $p$ the number of train unit
types. Each trip $j\in\{1,\ldots,n\}$ is defined by a request
$r_j$, given by the required number of passenger seats, a maximum
number $u_j$ of train units that can be assigned to the trip, and
its timetable. Each train unit type $k\in\{1,\ldots,p\}$ is
defined by a number $d^k$ of available train units and an
associated capacity $s^k$, given by the number of available seats.
We consider the natural and well-known Integer Linear Programming
(ILP) with arc variables, based on a canonical graph
representation of the problem (see e.g. \cite{CCT10}). Let
$G=(V,A)$ be a directed acyclic multigraph, in which nodes
correspond to trips, and the arc set $A$ is partitioned into $p$
subsets $A^1, \ldots, A^p$, where $A^k$ is associated with train
units of type $k$. The sequencing constraints are implicitly
represented by the graph. In particular, arc $(i,j)^k$ exists
whenever the time between the arrival of trip $i$ and the
departure of trip $j$ allows a train unit of type $k$ to travel
from the arrival station of trip $i$ to the departure station of
trip $j$. Let us introduce an integer variable $x_{ij}^k$
($k\in\{1,\ldots,p\}$, $i,j\in\{1,\ldots,n\}$), that indicates the
number of times that arc $(i,j)^k$ is selected in the solution,
i.e., the number of train units of type $k$ that execute trip $i$
before trip $j$ in the associated sequence. Moreover, let
$c_{ij}^k$ denote the cost of arc $(i,j)^k$, corresponding to the
time in minutes elapsing between the departure of the starting
node and the departure of the ending node. The ILP reads as
follows:
\begin{eqnarray}
\min \sum_{k=1}^p \sum_{i=1}^n \sum_{j=1}^n c_{ij}^k x_{ij}^k, &&
\label{eq:a1obj}\\
\sum_{i=1}^n x_{ij}^k = \sum_{i=1}^n x_{ji}^k, && ~~ k\in\{1,\ldots,p\},\ j\in\{1,\ldots,n\}, \label{eq:a1cons}\\
\sum_{i=1}^n \sum_{j=1}^n c_{ij}^k x_{ij}^k \leq 1440 ~ d^k, && ~~
k\in\{1,\ldots,p\},
\label{eq:a1disp}\\
\sum_{k=1}^p \sum_{i=1}^n s^k x_{ij}^k \geq r_j, && ~~
j\in\{1,\ldots,n\}, \label{eq:a1cap}\\
\sum_{k=1}^p \sum_{i=1}^n x_{ij}^k \leq u_j, &&
~~ j\in\{1,\ldots,n\}, \label{eq:a1bound}\\
x_{ij}^k \geq 0, \mbox{ integer }, && ~~ k\in\{1,\ldots,p\}, i,
j\in\{1,\ldots,n\}. \label{eq:a1pos}
\end{eqnarray}

The objective is to minimize the overall cost of the selected
arcs, equal to 1440 (the number of minutes in a day) times the
number of the train units globally used. Constraints
(\ref{eq:a1cons}) express the flow conservation. Constraints
(\ref{eq:a1disp}) forbid to use more than the available train
units for each type. Constraints (\ref{eq:a1cap}) impose to cover
each trip with the corresponding seat request. Constraints
(\ref{eq:a1bound}) impose a bound on the number of train units
that can be used to cover each trip. Note that the maintenance
constraints are not introduced in the formulation, since they
would be complex to model; we take them into account directly in
the heuristic algorithm.
%\subsection{Lagrangian Relaxation} \label{lagr}

The approach in \cite{CCT10}, based on solving LP relaxations by
general-purpose LP software, considers an alternative equivalent
ILP formulation with path variables, whose LP relaxation is much
faster to solve (modulo using column generation). In this work, we
consider a Lagrangian-relaxation based approach, which nicely
combines with the ILP above. Specifically, we first replace
constraints (\ref{eq:a1bound}) by the weaker \iffalse
$\sum_{i=1}^n x_{ij}^k \leq u_j, (k\in\{1,\ldots,p\},
j\in\{1,\ldots,n\})$ (5'),\fi
\begin{equation}
\sum_{i=1}^n x_{ij}^k \leq u_j,
~~ k\in\{1,\ldots,p\}, j\in\{1,\ldots,n\}, \label{eq:aboundnew}\\
\end{equation}
which impose that each train unit of type $k$ can cover each trip
$j$ at most $u_j$ times. Then, we relax in a Lagrangian way
constraints (\ref{eq:a1disp}) and (\ref{eq:a1cap}), by using
nonnegative Lagrangian multipliers $\lambda_j,
(j\in\{1,\ldots,n\}$) and $\sigma_k, (k\in\{1,\ldots,p\}$),
respectively.
%The relaxed model reads:
%\begin{eqnarray}
%\sum_{j=1}^n \lambda_j r_j - \sum_{k=1}^p \sigma_k ~ 1440 ~ d^k + \min \sum_{k=1}^p \sum_{i=1}^n \sum_{j=1}^n (c_{ij}^k -\lambda_j s^k + \sigma_k c_{ij}^k)x_{ij}^k, && \label{eq:objlagr}\\
%\sum_{i=1}^n x_{ij}^k = \sum_{i=1}^n x_{ji}^k, ~~ k=1,...,p, j=1,...,n, \label{eq:lcons}\\
%\sum_{i=1}^n x_{ij}^k \leq u_j^k , ~~ k=1,...,p, j=1,...,n, \label{eq:abound2}\\
%x_{ij}^k \geq 0, ~~ \mbox{ integer }, ~~ k=1,...,p, i=1,...,n,
%j=1,...,n. \label{eq:lpos}
%\end{eqnarray}
%For each $\lambda, \sigma \geq 0$, the optimal solution value of
%(\ref{eq:objlagr}) - (\ref{eq:lpos}) is a valid lower bound on the
%optimal solution value of (\ref{eq:a1obj}) - (\ref{eq:a1pos}) and
%therefore on the optimal solution value of the original problem.

The resulting Lagrangian relaxed problem decomposes onto $p$
independent subproblems, one for each train unit type. Let
$\tilde{c}_{ij}^k := (c_{ij}^k-\lambda_j s^k+ \sigma_k c_{ij}^k)$
be the Lagrangian cost for each arc $(i,j)^k \in A$. We impose a
null Lagrangian cost ($\tilde{c}_{ii}^k:=0$) for loops $(i,i)^k$
for each vertex $i \in V$ (that had originally infinity cost).
Thus, the subproblem associated with a train unit type $k$ calls
for the minimization of $\sum_{i=1}^n \sum_{j=1}^n
\tilde{c}_{ij}^k x_{ij}^k$ subject to (\ref{eq:a1cons}) and
(\ref{eq:aboundnew}). Moreover, the presence of zero-cost loops
allows us to replace inequality by equality in
(\ref{eq:aboundnew}). \iffalse reads:
\begin{eqnarray}
\min \sum_{i=1}^n \sum_{j=1}^n \tilde{c}_{ij}^k x_{ij}^k, \label{fobjr}\\
\sum_{i=1}^n x_{ij}^k & = & u_j, ~~ j\in\{1,...,n\}, \label{eq:fireq}\\
\sum_{i=1}^n x_{ji}^k & = & u_j, ~~ j\in\{1,...,n\}, \label{eq:seceq}\\
x_{ij}^k \geq 0, & \mbox{ integer } &, ~~i\in\{1,...,n\}, j\in\{1,...,n\}. \label{eq:varsign}
\end{eqnarray}
There is an immediate one-to-one correspondence between solutions
$\bar x^k$ of (\ref{eq:a1cons}) together with (5') and solutions
$\tilde x^k$ of (\ref{eq:fireq})-(\ref{eq:seceq}), given by
$\tilde x_{jj}^k:= u_j-\sum_{i \neq j} \bar x_{ij}^k$ for $j \in
V$. \fi It is well known that all vertices of the feasible region
of this subproblem are integer for $u_j$ integer and the
associated constraint matrix is totally unimodular. Moreover, in
the particular case, arising in our case study, in which $u_j$
does not depend on $j$, we can replace $u_j$ with $u$ for
$k\in\{1,\ldots,p\}$ and $j\in\{1,\ldots,n\}$. This makes the
subproblem equivalent to an Assignment Problem (AP), obtained by
replacing $u$ by $1$ in (\ref{eq:aboundnew}),
%\begin{eqnarray}
%\min \sum_{i=1}^n \sum_{j=1}^n \tilde{c}_{ij}^k y_{ij}, \label{objas}\\
%\sum_{i=1}^n y_{ij} & = & 1, ~~ j=1,...,n, \label{fireqas}\\
%\sum_{l=1}^n y_{jl} & = & 1, ~~ j=1,...,n, \label{seceqas}\\
%y_{ij} \geq 0, & \mbox{ integer } &, ~~i=1,...,n, j=1,...,n.
%\label{varsignas}
%\end{eqnarray}
the correspondence between solutions $\bar x$ of the subproblem
and $\bar y$ of AP being given by $\bar x=u \bar y$. As is well
known, the assignment problem can be solved in $O(n^3)$ time.
%, using, e.g., the Hungarian Algorithm.
%\subsection{A Lagrangian-based Heuristic} \label{heur}
In order to find good Lagrangian multipliers we apply a standard
iterative subgradient procedure.

%% \cite{HK72}.
%At each iteration of the subgradient procedure we solve assignment
%subproblems as described above, finding the optimal solution $\bar
%x$ of (\ref{eq:objlagr}) - (\ref{eq:lpos}) corresponding to the
%current Lagrangian multipliers. We check if $\bar x$ violates some
%of the relaxed constraints and update the multipliers by using the
%well-known formula: \vspace*{0.2cm}
%
% $ \lambda_j =
%\max\left\{\lambda_j+\delta \frac{UB-LB}{||\mu||^2}\\\mu_j,
%0\right\}$, $ \sigma_k = \max\left\{\sigma_k+\delta
%\frac{UB-LB}{||\tau||^2}\\\tau_k, 0\right\}$,\vspace*{0.2cm}
%
% where $\delta$ is a step parameter, LB is the best lower bound value found so far, UB
%is the best upper bound, i.e., value of a feasible solution of
%(\ref{eq:a1obj})-(\ref{eq:a1pos}), determined by a heuristic
%procedure to be illustrated later and $\mu, \tau$ are the
%subgradient vectors associated with $\bar x$, namely $\mu_j := r_j
%- \sum_{k=1}^p \sum_{i=1}^n s^k \bar x_{ij}^k$, for $j=1,\ldots,
%n,$ $\tau_k := \sum_{i=1}^n \sum_{j=1}^n c_{ij}^k \bar x_{ij}^k
%-1440 d^k$ for $k=1,\ldots, p$.
Besides yielding a valid lower bound on the optimal solution
value, the best Lagrangian multipliers $\lambda^*, \sigma^*$ found
throughout the iterations and the corresponding reduced costs
$\bar c_{ij}^k := c_{ij}^k - \lambda_j^* s^k + \sigma_k^* c_{ij}^k
- \bar w_i^k -\bar v_j^k$, where $\bar w^k$ and $\bar v^k$ are the
optimal AP dual variables associated with the assignment
constraints, are used to drive the following constructive
Lagrangian heuristic algorithm. We order the train unit types for
decreasing capacity values, and construct the workload for each
train unit, until either all the trips have been covered or all
the train units have been used (in this case, we apply a local
search procedure to obtain a feasible solution).
%Finally, we apply a
%local search procedure in order to improve the solution found.
The construction starts with the selection of a trip whose
departure is in the beginning of the day. Then, we choose the
following trips by assigning to each one a score that takes into
account the reduced costs, the original costs and how ``well" the
current train unit can cover the trip (taking into account which
other train units are still available). In addition, we give a
prize to the arcs that allow to perform maintenance, until we
reach the number of necessary maintenance operations. Every time a
trip is selected in the solution, we update its number of seats
and the reduced cost of each arc entering the corresponding node.
In order to decide how to end the workload of the current train
unit type, for each selected trip in the workload we solve the
relaxed problem on the residual train units and trips. The
solution value of the reduced relaxed problem gives a lower bound
on the global number of train units that are needed to cover all
the residual trips. When no more train units are available of the
current type, we end the workload with the trip giving the
smallest solution value.

\vspace*{-.5cm}

\section{Computational Experiments} \label{res}

\vspace*{-.5cm}

We present some preliminary computational experiments on a set of
real-world instances for an operator running trains in a regional
area, and compare the results with the approach presented in
\cite{CCT10}, with and without imposing the maintenance
constraints. The tests were performed on a PC Pentium 4, 3.2 GHz,
2 GB RAM, and using Cplex 9.0 as an LP-solver in \cite{CCT10}. The
results are presented in Table \ref{ta:res}, showing that the
proposed approach produces comparable results within much shorter
computing time (expressed in seconds).
%, where \emph{value}
%represents the number of train units, and \emph{time} is expressed
%in seconds.

\begin{table} [ht!]
\small
%\baselineskip=2pt
\begin{center}
\begin{tabular}{|c|c|c|c|c|c|c|c|c|c|c|}
\hline & && \multicolumn{2}{|c|}{Lagr. heur.} &
\multicolumn{2}{|c|}{\cite{CCT10} heur.} &
\multicolumn{2}{|c|}{Lagr. heur. maint.} &
\multicolumn{2}{|c|}{\cite{CCT10} heur. maint.}\\
inst.\ & $n$ & $p$ & value & time & value & time & value & time & value & time\\
\hline
1 & 85 & 1 & 2 & 0 & 2 & 0 & 2 & 0 & 2 & 0\\
2 & 120 & 1 & 4 & 0 & 4 & 0  & 4 & 1 & 4 & 5\\
3 & 302 & 1 & 17 & 4 & 17 & 288  & 17 & 5 & 17 & 544\\
4 & 208 & 2 & 26 & 4 & 25 & 17 & 27 & 3 & 25 & 19 \\
5 & 364 & 2 & 20 & 18 & 20 & 1912  & 21 & 20 & 20 & 3899\\
\hline
\end{tabular}
\end{center}
\caption{\label{ta:res} Comparison on a set of real-world
instances in the case without or with the maintenance
constraints.}
\end{table}

% The Appendices part is started with the command \appendix;
% appendix sections are then done as normal sections
% \appendix

% \section{}
% \label{}

\vspace*{-.5cm}

\begin{thebibliography}{00}

\vspace*{-.5cm}

% \bibitem{label}
% Text of bibliographic item

% notes:
% \bibitem{label} \note

% subbibitems:
% \begin{subbibitems}{label}
% \bibitem{label1}
% \bibitem{label2}
% If there is a note, it should come last:
% \bibitem{label3} \note
% \end{subbibitems}

\bibitem{CCT10} Cacchiani, V., Caprara, A., and Toth, P., ``Solving a
Real-World Train Unit Assignment Problem", {\em Mathematical
Programming Series B}, to appear (2010).

\bibitem{Cetal05} Caprara, A., Kroon, L., Monaci, M., Peeters, M., and Toth, P.,
``Passenger Railway Optimization", in C. Barnhart and G. Laporte
(eds.), {\em Handbooks in OR \& MS}, Vol. 12, Elsevier Science
(2006).

\end{thebibliography}

\end{document}
