# Optimal routing and scheduling in a time-dependent and directed multi-graph:
a multi-agent temporal rural postman problem

Somnath Buriuly IITB-Monash Research Academy, India
Email<EMAIL_ADDRESS>Leena Vachhani Systems and Control
Engineering
IIT Bombay, India
Email<EMAIL_ADDRESS>Arpita Sinha  Systems and Control
Engineering
IIT Bombay, India
Email<EMAIL_ADDRESS>Sunita Chauhan  Department of Mechanical and
Aerospace Engineering
Monash University, Australia
Email<EMAIL_ADDRESS>Sivapragasam Ravitharan  Department of
Mechanical and Aerospace Engineering
Monash University, Australia
Email<EMAIL_ADDRESS>
###### Abstract

Time-dependent (or temporal) properties may arise in many practical static
network search problems, for example routing and scheduling of railway track
inspection where availability of track for inspection depends on the rail-
schedule. In the absence of side constraints this routing and scheduling
problem is a multi-agent rural postman problem on a temporal directed network;
having temporal attributes motivated from railway track unavailabilities. This
work presents a novel Mixed Integer Linear Programming (MILP) formulation
model for a Multi-agent Temporal Rural Postman Problem (MTRPP), and proposes
theory for temporal graphs to support the model. A modified Benders’
decomposition algorithm is also proposed to address the non-convex
unavailability constraints occurring in its sub-problems. Furthermore, a
polynomial time algorithm is designed to address this sub-problem, while
restraining the NP-Hard part to the master problem. A simulation study showing
application in rail-track inspection using track-based-agents is presented,
alongside comparison study showing improved results in smaller and sparser
graphs compared to Time-Dependent Rural Postman Problem (TDRPP) which is
developed over temporal attributes arising due to congestion in roadway
networks.

##  I Introduction

Railways feature a vast infrastructure capable of servicing multiple
transportation objectives. Capacity restrictions of the infrastructure and
increasing demand of services motivates sharing of this resource, hence
introducing diverse challenges. Some of these challenges are reduced to
optimization problems, like efficient utilization of allocated resources and
optimal sharing among different units. Efficient management of this vast
railway network is extremely difficult for any isolated computational tool.
Hence, a distributed approach is a most preferred zeroed-choice, leading to
hierarchical optimization problems i.e. important ones being solved first. It
is in the hands of Infrastructure Manager (IM) to decide the hierarchy,
segregate sub-problems, and assign them to suitable servicing units called
Railway Undertakings (RUs). (Schlechte, 2012 [1]) and (Peeters, 2003 [2])
describe various planning phases in details while introducing the railway
system in their respective dissertations.

Urbanization has increased the demand of railway services therby increasing
wear and tear of the infrastructural resources. To satisfy demand while
maintaining network reliability, frequent monitoring of the infrastructure is
necessary. This article addresses a routing and scheduling challenge faced in
inspection planning of railway tracks in a time-shared environment. The core
problem is a Rural Postman Problem (RPP) with time-dependencies occurring due
to sharing of railway infrastructure for various services. Under simplifying
assumptions that ignores capacity and other practical restrictions, this
problem can be stated as a multi-agent time-dependent rural postman problem.
In brief, given a graph $G$ with arc set $A$, a set of arcs ($A_{*}\subset A$)
to be visited, a set of agents $\mathcal{K}$ and unavailability period of
arcs; visit each arc once using any one of the agents, in minimum total fuel
cost.

Temporal (or time-dependent) attributes may arise in a network-based problem
due to many factors, like congestion in roadways and unavailabilities in
railways. (Cordeau, 2014 [3]) presents Time-Dependent Traveling Salesman
Problem (TDTSP) and (Calogiuri, 2019 [4]) presents Time-Dependent Rural
Postman Problem (TDRPP) that reflects temporal attributes in roadway network.
Similarly, (Lannez, 2015, [5]) presents a routing and scheduling problem for
inspection of railway tracks using ultrasonic inspection vehicles, with
temporal attributes modeled with integer variables using unavailability
calendars. However, block unavailability due to scheduled trains as well as
(maintenance) possessions contribute to the temporal attributes, requires
modeling the constraints with both integer and linear variables.

The major contributions of this work are: (1) a new way of modeling multi-
agent temporal rural postman problem, (2) an algorithm, having structure
similar to Benders’ decomposition algorithm [6], but with modified steps to
improve computation times of the sup-problems, and address non-convex
constraints generated due to unavailabilities, (3) simulation study that shows
the benefits of this approach of addressing temporal graphs, and (4) an
application of proposed work in a realistic problem of scheduling inspection
vehicles for railway network without disturbing the regular train schedule.

### Literature survey

The problem of railway track inspection is modeled as a routing problem where
multiple subsequent assignments of interdependent decisions are obtained.
Routing problems are categorized into two categories: Vehicle Routing Rroblem
(VRP) and Arc Routing Problem (ARP). The main focus in VRPs are point
locations (vertices), while in ARPs are the inter-relations linking these
points (arcs). The primary objective in these routing problems is to optimize
traversal cost of the vehicle(s) such that a subset of either vertices (for
VRP) or arcs (for ARP) of the network (graph) are visited at least once. The
general problem is very ideal and ignores many practical restrictions like
network and vehicle capacity, inspection and maintenance time slots,
heterogeneous vehicle characteristics, etc.

Based on practical considerations, many variants of VRPs and ARPs have been
introduced in the literature, like VRP with Time-Window (TW), Capacitated Arc
Routing Problem (CARP), etc. The two article series by (Eiselt, 1995a [7] and
1995b [8]) serves as a good reference to understand categories of arc routing
problems. Literature exists for both, the general problem (mostly NP-Hard) as
well as the practical problem (with many side-constraints) - proposing exact
solutions as well as heuristics. Since most practical problems have vast
network size, exact solutions are rarely attempted in such work due to
provably large computation times.

The railway track inspection problem may be modeled as a VRP or an ARP (with
our focus on RPP based representations). A general hard problem in roadway
setting has been studied by (Calogiuri, 2019 [4]), which attempts a time-
dependent RPP in a traffic environment. The authors use a branch and bound
algorithm to produce an exact solution for this problem. Unlike TDRPP, where
the temporal variant of the graph is traversable if the underlying graph is
connected, railway based unavailability scenarios can result in disconnected
sub-graphs. For example: a track section, which is unavailable for $T$ time,
can never be traversed by an agent that has traversal time more than $T$.
Hence we propose a stronger connectivity condition for agent-network
compatibility. Furthermore, TDRPP addresses a single-agent case and also
assumes a First-In-First-Out (FIFO) rule i.e. if a vehicle leaves a vertex
(say $u$) for a vertex (say $v$) at a given time, leaving vertex $u$ at a
later time implies arriving later at vertex $v$. Due to such subtle
differences in the problem statements and also for comparison convenience, the
core portion of the time-dependent routing and scheduling problem, tackled in
our work, is termed as Multi-agent Temporal Rural Postman Problem (MTRPP). A
general version of this ARP, for railway setting, can be observed from article
by (Lannez, 2015 [5]), where an exact solution is proposed using Danzig-Wolfe
decomposition and column generation. Lannez’s article models this inspection
routing and scheduling as an Integer Programming (IP) problem by disretizing
time, which is avoided in MTRPP, hence modeled in a Mixed Integer Linear
Programming (MILP) framework.

From a practical perspective, aside from addressing large size network, many
side-constraints are mandatory for implementing inspection schedules; for
example, headway constraints are one of the basic constraints that preserve
signalling aspect of railway protocols (implicitly preventing overtaking among
vehicles). As a consequence, heuristics are proposed to generate a plan in
reasonable computational time span. (Lannez, 2015 [5]) considers inspection
vehicle (with crew) routing and scheduling problem with worker shift hours as
well as capacitated vehicles in the problem setting. Lannez’s article presents
a Danzig-Wolfe decomposition and column generation based matheuristic for
generating the route and schedule for a French railway network with 760 nodes
and 2100 arcs.

In literature, maintenance scheduling problems are modeled in a macro-
perspective i.e. the path generation is ignored in planning phase. As a
consequence, these scheduling problems are modeled as VRP, where the vertices,
reflecting service-requirement, must be visited by the agents in an optimal
manner while satisfying practical constraints. (Pour, 2017 [9]) and (Budai,
2005 [10]) developed mathematical models for preventive maintenance scheduling
problem. (Pour, 2017 [9]) used a MILP model for Danish railway system, and
obtained schedules over an eight-week period using a modified version of
constraint programming technique. Several managerial parameters are also
considered while constructing cost and constraints e.g. time-dependencies on
crew schedule, crew competencies, task splitting over multiple days, etc.
(Budai, 2005 [10]) used an IP framework for scheduling periodic maintenance
tasks, and analyzed performances of various heuristic strategies. (Dao, 2018
[11]) and (Higgins, 1999 [12]) also developed a maintenance scheduling model
while considering track availability. Such track availability based
constraints are modeled to avoid large possession time for maintenance,
causing minimum disruption in regular schedules of the trains; especially
suited to busy network. (Higgins, 1999 [12]) used a tabu search based
heuristic by defining neighbourhood solutions and searching over it. On the
other hand, (Dao, 2018 [11]) used standard commercial solver to optimize the
MILP model.

Some inspection tasks are more suitable in the VRP setting (rather than ARP)
due to large service time and small headway gap in the network. For example,
(Peng, 2014 [13]) modeled a routing and scheduling problem for periodic
inspection of railway track, and proposed a heuristic algorithm for efficient
solution. The claim is backed up by real-world case studies that significantly
outperforms commercial solvers as well as manual methods. In all these
literature, on inspection and maintenance scheduling, sufficient priority is
given to reliability causing disruption in regular train schedules to a
minimum.

Our work on MTRPP is suitable for frequent inspection tasks with small service
times, hence prioritizing regular train schedules over inspection tasks. This
also motivates introduction of track-based inspection robots (agents) to
railway system; for example camera based inspection of missing track-
components developed by (Resendiz, 2013 [14]) and alternating current field
measurement based sensors developed by (Papaelias, 2009 [15]; and Nicholson,
2012 [16]) can be utilized for supporting autonomous inspection tasks.

Readers are familiarized with some preliminaries in the next Section (Sec-II);
followed by the main formulation (Sec-III) and an exact algorithm (Sec-IV).
The methodology is supported by simulation examples, case studies and
comparison studies in the succeeding section (Sec-V). Finally a conclusion
(Sec-VI) summarizes the main findings, alongside its projection into future
research.

##  II Preliminaries

The next sub-section describes some commonly used terms to build our work.
Benders’ decomposition technique is reviewed in the following sub-section to
add clarity in understanding the modifications proposed for addressing
temporal RPPs.

###  II-A Definitions

Directed multi-graph ($G$): A directed multi-graph is a graph which may have
multiple arcs connecting same pair of vertices, in same direction. It is
generally denoted using one vertex set ($V$), and an arc multi-set111A multi-
set may contain multiple entires of same element, e.g. $\\{1,1,2,3\\}$.
consisting of two tuples which is a subset of $V\times V$. In this work,
vertex and arc sets of directed multi-graph are represented using two distinct
sets related by the incidence matrix ($B$). Hence, elements of arc set $A$ are
uniquely represented.

Incidence Matrix (B): A matrix of dimension $|V|\times|A|$, representing the
relationship between vertices (rows) and arcs (columns). Its entries are
composed of three elements $\\{-1,0,1\\}$. Every column representation of an
arc has one negative entry and one positive entry - the negative entry
correspond to the row/vertex at which the arc originates and the positive
entry correspond to the row/vertex where the arc terminates.

Strongly connected: A directed graph is said to be strongly connected if for
any pair of vertices, there exist a path connecting the two, in both
directions.

Components: Subsets of a graph or subgraphs in which any two vertices within
the subgraph are strongly connected while vertices picked from different
subgraphs are disconnected.

Service arcs ($A_{*}$): A task, such as inspection of railway tracks, that
requires an agent to visit certain arc of a graph is called service, and the
corresponding arc is called service arc. Visiting a service arc implies that
the agent is servicing the arc, thus making the task intrinsic to selection of
an arc.

Deadhead arcs: Movement of an agent over an arc without doing any service is
called cruising, and the corresponding arc is called cruise arc or deadhead
arc. Similar to servicing, cruising is also intrinsic to arc selection.

Running time: Time taken by an agent to traverse an arc is called the running
time of the agent for that particular arc. For service arc, the time quantity
is termed as service time, and for deadhead arc, its termed as deadhead or
cruise time.

Head and tail vertex: For any arc $a$, the vertex to which an arc directs is
called head vertex (sat $v_{2}$) while the other is called tail vertex (say
$v_{1}$), see Figure 1. The arc $a$ is incoming at the head vertex $v$ and
outgoing at the tail vertex $u$.

Figure 1: Head vertex ($v_{2}$) and tail vertex ($v_{1}$) of the arc $a$.

Vertex-disjoint path: A path i.e. sequence of alternating vertices and arcs,
that does not repeat same vertex more than once is called a vertex-disjoint
path.

Cone: A cone is a set $S$ of vectors whose conic combination is closed i.e.
for any two vector $u_{1},u_{2}\in S$, all its conic combination
$\alpha_{1}u_{1}+\alpha_{2}u_{2}$ also lies in the cone; where
$\alpha_{1},\alpha_{2}\geq 0$.

Extreme ray: A vector $u\in S$ is called an extreme ray if it cannot be
represented using convex combination of non-collinear vectors of $S$, i.e. no
$u_{1},u_{2}\in S$ exist such that $u_{1}\neq\alpha u_{2}$ and
$u=\alpha_{1}u_{1}+(1-\alpha_{1})u_{2}$, with $\alpha\geq 0$ and
$\alpha_{1}\in(0,1)$.

###  II-B Benders’ decomposition

The key idea in this decomposition technique is to exploit the structural
properties of a given LP/MILP optimization problem. It is an iterative process
where an optimization problem is broken into one reduced master problem and
one or more sub-problems. The feasible region of reduced master problem
decreases at every iteration till the optimal solution of the reduced master
problems is same as that of actual optimization problem (verified by zero
duality gap). Following text explains the steps and advantages in details,
however readers are requested to refer standard text on Benders’ decomposition
by (Conejo, 2006 [6]) for further details.

Let an MILP minimization problem have the following structure:

$\displaystyle\underset{X,Y}{\text{min}}$
$\displaystyle\mathcal{W}^{\top}X+\mathcal{C}^{\top}Y;\ $
$\displaystyle\mathcal{W}\in\mathbb{R}^{m},\mathcal{C}\in\mathbb{R}^{n}\ \&\
\mathcal{C}\geq 0$ s.t. $\displaystyle AX\leq b;\ $ $\displaystyle
A\in\mathbb{R}^{q_{1}\times m},b\in\mathbb{R}^{q_{1}}$ $\displaystyle
A_{10}X+A_{11}Y\leq b_{1};\ $ $\displaystyle A_{10}\in\mathbb{R}^{q_{2}\times
m},A_{11}\in\mathbb{R}^{q_{2}\times n},b_{1}\in\mathbb{R}^{q_{2}}$
$\displaystyle X\in\\{0,1\\}^{m},\ \ Y\geq 0;\ $ $\displaystyle
Y\in\mathbb{R}^{n},$

where $\mathbb{R}$ represents set of real numbers.

Since variable $X$ is common among all constraints, its computed separately
from variable $Y$ (the advantage of which will be clear in the sub-problem
structure). This is called reduced master problem II-B, given as:

$\displaystyle\underset{z,X}{\text{min}}$ $\displaystyle z$ s.t.
$\displaystyle z\geq\mathcal{W}^{\top}X$ $\displaystyle X\in S_{X}:=\\{AX\leq
b,X\in\\{0,1\\}^{m}\\}$

The variable $Y$ is associated with the second set of constraints only, see
II-B. Given any $\bar{X}_{i}$ ($i$ refers iteration number) from the set
$S_{X}$, the remaining problem can be formulated as a separate LP sub-problem:

$\displaystyle\mathcal{W}^{\top}\bar{X}_{i}+$
$\displaystyle\underset{Y}{\text{min}}$ $\displaystyle\mathcal{C}^{\top}Y$
s.t. $\displaystyle A_{11}Y\leq b_{1}-A_{10}\bar{X}_{i}$ $\displaystyle Y\geq
0$

Benders’ decomposition is most suitable if the sub-problem matrix is block
diagonal, while $X$ serves as a linking variable. In such a case, the sub-
problem SP can be decoupled and solved independently for each matrix block of
constraint; this enhances the computation time considerably.

Since the main problem is a minimization problem with non-negative
$\mathcal{C}$ and all non-negative variables, only two cases are possible (by
LP duality): (1) the sub-problem has a finite solution (the corresponding dual
also has finite solution) and (2) the sub-problem has no solution (the
corresponding dual has unbounded solution).

Case I: If dual is finite, the following relations are established using
strong duality:

$\displaystyle\mathcal{C}^{\top}Y\geq\bar{\lambda}_{i}^{\top}(b_{1}-A_{10}\bar{X}_{i})\textit{,
for feasible }Y\textit{ and }\bar{\lambda}_{i}$ (1) and
$\displaystyle\mathcal{C}^{\top}\bar{Y}_{i}=\bar{\lambda}_{i}^{\top}(b_{1}-A_{10}\bar{X}_{i})\textit{,
if }\bar{Y}_{i}\textit{ and }\bar{\lambda}_{i}\textit{ are respective
optimizers }$

Note that, $\bar{\lambda}_{i}$ is the dual optimizer, given as
$argmax\\{\lambda^{\top}(b_{1}-A_{10}\bar{X}_{i})|A_{11}^{\top}\lambda\leq\mathcal{C},\lambda\leq
0\\}$. Since the dual feasible set is independent of $X$, the above relations
(see Equations 1) are true for any arbitrary $\bar{X}_{i}$, irrespective of
the iteration number. Hence the following inequality can be established for
the cost function of II-B:

$\displaystyle c(X,Y)$ $\displaystyle:=$
$\displaystyle\mathcal{W}^{\top}X+\mathcal{C}^{\top}Y$ $\displaystyle\geq$
$\displaystyle\mathcal{W}^{\top}X+\bar{\lambda}_{i}^{\top}(b_{1}-A_{10}X)=:e_{i}(X),$
$\displaystyle\forall X\in S_{X}\textit{ and }Y\in\\{Y\geq 0|A_{11}Y\leq
b_{1}-A_{10}X\\}$

and

$\displaystyle c(\bar{X}_{i},\bar{Y}_{i})$ $\displaystyle=e_{i}(\bar{X}_{i})$

In other words, $e_{i}(X)$ is smaller than $\underset{Y}{{min}}\ c(X,Y)$
except at $\bar{X}_{i}$ where its equal. Hence if the master problem is
evaluated again (or iterated) with this additional constraint (see II-B), all
$X$’s that have $\mathcal{W}^{\top}X<e_{i}(\bar{X}_{i})$, yield a smaller $z$
value compared to $\bar{X}_{i}$. If no such $X$ exist at the $i^{th}$
iteration, then $z$ ($=e_{i}(\bar{X}_{i})$) is optimal, and the optimizer of P
is $(\bar{X}_{i},\bar{Y}_{i})$.

$\displaystyle\underset{z,X}{\text{min}}$ $\displaystyle z$ s.t.
$\displaystyle z\geq\mathcal{W}^{\top}X$ $\displaystyle
z\geq\mathcal{W}^{\top}X+\bar{\lambda}^{\top}(b_{1}-A_{10}X),\
\forall\bar{\lambda}\in\bar{S}_{i}$ $\displaystyle AX\leq b$ $\displaystyle
X\in\\{0,1\\}^{m}$

where, $\bar{S}_{i}$ collects all dual variables achieved till $i^{th}$
iteration.

Case II: If the II-B is infeasible, the dual optimizer of the following
feasibility problem (every feasible solution is optimal) is computed:

$\displaystyle\underset{Y}{\text{min}}$ $\displaystyle 0^{\top}Y$ s.t.
$\displaystyle A_{11}Y\leq b_{1}-A_{10}\bar{X}_{i}$ $\displaystyle Y\geq 0$

In case of infeasible II-B, (by LP duality) there exist a non-trivial dual
optimizer $\bar{\lambda}_{i}:=\alpha\hat{\lambda}_{i}$
$(||\hat{\lambda}_{i}||_{\infty}=1)$, which is feasible for
$\alpha\in[0,\infty)$. $\hat{\lambda}_{i}$ is called unbounded extreme ray as
it is an extreme ray (see Section II-A) that can be scaled to infinity inside
the dual feasible set. Note that $A_{11}^{\top}\bar{\lambda}_{i}\leq
0\leq\mathcal{C}$, hence its feasible for dual of II-B and the dual cost is
unbounded. This also confirms the LP duality argument.

Again, using duality, the inequality in Equation (2) holds for all $X$,
because $\bar{\lambda}_{i}$ is given as
$argmax\\{\lambda^{\top}(b_{1}-A_{10}\bar{X}_{i})|A_{11}^{\top}\lambda\leq
0,\lambda\leq 0\\}$ whose feasible region is independent of $X$.

$\displaystyle 0\geq\hat{\lambda}_{i}^{\top}(b_{1}-A_{10}X),\ \forall X\in
S_{X}$ (2)

All feasible $Y$ produces zero cost, hence corresponding dual variables should
produce zero dual cost. This implies, on flipping the inequality in Equation
(2), the $X$’s that produce infeasible sub-problem are eliminated. Adding this
constraint to the master problem, an iterative Benders’ decomposition
algorithm is obtained, with reduced master problem given as follows:

$\displaystyle\underset{z,X}{\text{min}}$ $\displaystyle z$ s.t.
$\displaystyle z\geq\mathcal{W}^{\top}X$ $\displaystyle
z\geq\mathcal{W}^{\top}X+\bar{\lambda}^{\top}(b_{1}-A_{10}X),\
\forall\bar{\lambda}\in\bar{S}_{i}$
$\displaystyle\hat{\lambda}^{\top}(b_{1}-A_{10}X)\geq 0,\
\forall\hat{\lambda}\in\hat{S}_{i}$ $\displaystyle AX\leq b$ $\displaystyle
X\in\\{0,1\\}^{m}$

where, $\bar{S}_{i}$ collects all dual optimizers $\bar{\lambda}_{i}$ obtained
from II-B, and $\hat{S}_{i}$ collects all unbounded extreme rays
$\hat{\lambda}_{i}$ obtained from II-B; achieved till $i^{th}$ iteration.

The Benders’ decomposition converges when the lower bound $z$ obtained from
II-B meets the upper bound cost, which is obtained by evaluating
$c(\bar{X}_{i},\bar{Y}_{i})$.

##  III Proposed MTRPP Formulation

### Problem formulation

The general problem of MTRPP can be described as follows: given a directed
multi-graph network $G=(V,A,B)$ whose arc set $A$ represents various directed
routes, and the vertex set $V$ represents connection points of these arcs. The
elements of arc set $A$ are denoted by single variable instead of $2$-tuple,
such that every parallel arc is assigned a unique identifier, unlike $2$-tuple
based representation where arcs between same vertex pairs have same
identifier. The relation between $V$ and $A$ is captured in the incidence
matrix $B$. A list of unavailabilities ($Z=\cup_{a\in A}Z_{a}$) is provided
for each arc such that an entry $(\underline{\omega},\overline{\omega})\in
Z_{a}$ describes the time period of duration
$\overline{\omega}-\underline{\omega}$ between which arc $a$ is unavailable. A
heterogeneous set of agents $\mathcal{K}$ is also provided that starts their
respective journeys from one of the locations in their depot set $D_{k}$ and
end at any one of the locations given in its exit set $E_{k}$. The running-
time data is available for each agent-arc pair in column vector $R$ of size
$|\mathcal{K}||A|$. The graph-agent data given by $(G,Z,R)$ is analyzed for
compatibility such that traversality of agents in the network can be
guaranteed. We propose two properties for the graph-agent data $(G,Z,R)$,
namely well-defined and $T_{p}$-periodically connected in Section IV, which
guarantees this compatibility. The objective for this optimization problem is
to find trajectories for each agent that penalizes movement and waiting, while
ensuring that the service arc set $A_{*}(\subset A)$ is traversed once by any
of the agents. In this work, it is assumed that the vertices are always
available and they don’t have any capacity restrictions.

### Methodology

In order to avoid discretization of time, this work presents an ARP
formulation in MILP framework. The spatial (movement) and temporal (time-
stamp) descriptions must be captured in the solution to denote trajectories of
the agents. Spatial description is captured by arcs in all ARP based
formulations. Since the structure of graph naturally encodes the
interrelations among the arcs and vertices, assigning temporal variables at
the vertices can store time-description of movement of an agent. For example,
$\\{^{1}u_{v_{2}},^{1}x_{a_{4}},^{1}u_{v_{3}},^{1}x_{a_{6}},^{2}u_{v_{2}},^{1}x_{a_{3}},^{1}u_{v_{1}}\\}$
indicates one such trajectory in a three-vertex multi-graph, shown in Figure
2, where $u$’s are time-stamps at the vertices and $x$’s are indication of
arcs chosen for traversal. In particular, the agent starts from vertex $v_{2}$
at time ${}^{1}u_{v_{2}}$ (the superscript ‘$1$’ indicates order of visit),
follows arc $a_{4}$ to arrive at vertex $v_{3}$ at time ${}^{1}u_{v_{3}}$,
then follows arc $a_{6}$ to re-visit vertex $v_{2}$ at time ${}^{2}u_{v_{2}}$,
and so on. Thus running-time over any occupied arc can be observed by
evaluating the time difference between the two temporal variables, at the head
and tail vertices of this arc, with surplus time indicating waiting at head
vertex. The trajectories for each agent is described by a sequence of
alternating (and maybe repeating) vertices and arcs, where the arcs represent
the movement in space while the vertices represent the time stamp at the end
of each movement.

Search for such a sequence for each agent $k\in\mathcal{K}$, from depot set
$D_{k}$ to its exit set $E_{k}$, requires multiple spatiotemporal descriptions
for each vertex and arc in the graph $G$. Single variable for describing each
spatial and temporal description endorses an intrinsic vertex-disjoint
attribute in the required solution. Assigning multiple spatial and temporal
decision variables for same arc and defining interlinking constraints among
all such decisions can resolve this intrinsic requirement. These
interrelations can be visualized in a graph based representation similar to
that of TSN, see (Lannez, 2015 [5]; Peng, 2011 [17]; and Michail, 2016 [18]).
One of our major contribution is in constructing the graph based
representation, termed as Replicated graph which forms the base of our MILP
formulation for temporal graph multi-agent routing problem, which is one of
the major contributions in this work. Replicated graph can be viewed as one
among many possible coalitions of network data with agent data.

Figure 2: An example graph picked from article by (Lannez, 2015 [5]); vertices
$v_{1},v_{2},v_{3}\in V$, arcs $a_{1},a_{2},a_{3},a_{4},a_{5},a_{6}\in A$, and
service arcs $a_{2},a_{5}\in A_{*}$.

\-

Figure 3: Replicated graph $\mathcal{G}$ constructed for a simple 3-vertex
6-arc graph $G$ (see Figure 2) and two agents having vertex $v_{1}$ in both
depot and exit set. The deadhead arcs of $G$ connect the intra-layer vertices
in $\mathcal{G}$ while the service arcs connect the inter-layer vertices. The
arcs incident at the source and sink vertices are virtual arcs with zero
running-time (can be modeled to show fixed cost of an agent). The source and
sink vertices are assumed to be in layer zero of each component, such that
collapsing all non-zero layers gives the parent graph $G$. An arc connecting
source to sink vertex reflects decision for non-utilized agent; this action
maybe penalized if needed.

#### Construction of replicated graph

Consider the deadhead graph $G^{\prime}:=(V,A\backslash A_{*},B^{\prime})$
constructed from the directed multi-graph $G$, where service arcs $A_{*}$ are
removed from the arc set $A$, and $B^{\prime}$ is incidence matrix based on
this new arc set. This deadhead graph $G^{\prime}$ is re-drawn over
$|A_{*}|+1$ layers for each agent, and connected downwards using arcs in
$A_{*}$. The resulting graph has $|\mathcal{K}|$ (graph) components, one for
each agent, as shown in Figure 3. This new layered graph is termed as
replicated graph
$\mathcal{G}=(\mathcal{V},\mathcal{A},\mathcal{B},\mathcal{W})$; where
$(v,k,l)\in\mathcal{V}$, $(a,k,l)\in\mathcal{A}$ and $(w,k,l)\in\mathcal{W}$
represent vertices, arcs and arc-weights respectively; for agent
$k\in\mathcal{K}$ and layer $l\in\mathcal{L}$
($\mathcal{L}:=\\{1,\dots,|A_{*}|+1\\}$), with entries $v$ and $a$ picked from
sets $V$ and $A$ respectively. The arc-weights represent the minimum time an
agent $k$ requires to traverse the particular arc of layer $l$; consequently
the cardinality of $\mathcal{A}$ equals that of $\mathcal{W}$, i.e.
$|\mathcal{A}|=|\mathcal{W}|=(|A||\mathcal{L}|-|A_{*}|)|\mathcal{K}|$. This
set is constructed using the matrix $R$ consisting of agent running-time data.
The incidence matrix $\mathcal{B}$ encodes the network structure i.e. the
relationship between set $\mathcal{V}$ and $\mathcal{A}$. Pair of virtual
vertices, called source ($s$) and sink ($t$) vertices, are introduced for each
agent in the zeroth layer as $(s,k,0)$ and $(t,k,0)$. Source and sink vertices
are connected to the depot set and exit set respectively using virtual arcs
(zero running-time) such that trajectory always starts at source, visits a
vertex in depot set, services some arcs, visits a vertex of exit set, and
finally ends at the sink vertex. A source-sink virtual arc is also introduced
for each agent to model no movement scenario. The decision variables are
represented as follows:

  * •

$X_{akl}$ \- (spatial) binary decision variable, where $X_{akl}=1$ implies
that arc $a$ of layer $l$ is occupied by agent $k$

  * •

$Y_{vkl}$ \- (temporal) real variable which captures the time when agent $k$
leaves vertex $v$ of layer $l$

The spatial and temporal decision variables are represented in vector form as
$X$ and $Y$ respectively. Note that throughout this article, the terms vertex
and arc are used for elements of vertex set and arc set of replicated graph
$\mathcal{G}$, while elements from directed multi-graph $G$ are explicitly
specified. We next formulate the optimization problem by designing a cost
function and suitable constraints.

#### Cost function

The running-time data $R$ is available apriori which depends on the length of
the track section represented by the arc. Since fuel consumed by an agent is
proportional to the distance traversed by it, which in turn is proportional to
the time for constant speed agents, the cost of fuel can be represented using
the arc-weights of graph $\mathcal{G}$. Furthermore, waiting at any vertex
without doing inspection requires monitoring, tracking and presence of crew
(for agents that are not fully autonomous). Hence, a penalty term $\beta$ is
introduced to the cost function as shown in Equation (3). Balance between the
penalties for fuel and waiting depends on resources required in both
situations, however, in general, the penalty for fuel is assumed larger than
the penalty for waiting i.e. $\beta\in(0,\infty)$. In this work, $\beta$ is
chosen as $1$ implying double penalty for per unit running-time (both service
and deadhead) compared to waiting time. Since service cost is constant over
all feasible solutions, this contribution to the overall cost is ignored.

total cost $\displaystyle=\beta\mathcal{W}^{\top}X+\mathcal{C}^{\top}Y$ (3)
$\displaystyle(\mathcal{C})_{vkl}$ $\displaystyle=\begin{cases}1&\textit{for
sink vertices }(t,k,0),\forall k\in\mathcal{K}\\\
0&\textit{otherwise}\end{cases}$

The vector $\mathcal{C}$ selects the sink vertices of each agent and sums up
its temporal values, thus evaluating the completion times. The spatial part of
the cost function determines the surplus cost on movement.

#### Constraints

All restrictions of this problem statement are modeled with four types of
constraints, the first two being spatial, while the remaining two being
spatiotemporal. The flow, and service constraints are of spatial nature, i.e.
they only depend on the spatial decision variables. The running-time and
unavailability constraints depend on both the spatial and temporal decision
variables. Although this problem is formulated as an MILP problem, the
unavailability constraints are intrinsically non-convex. In particular, the
unavailability constraints can be viewed as a set of disjunctive (convex)
hyper-cubes in temporal space, due to single (temporal) variable based
constraint expressions.

The flow constraints ensure that the subgraph given by the occupied arcs has
sum of in-degree and out-degree as even for all its vertices, except source
vertices $(s,k,0)$ and sink vertices $(t,k,0)$ of each component. Equation (4)
describes the flow balance constraints for each vertex using the incidence
matrix. This flow conservation constraint ensures balanced selection of
incoming and outgoing arcs, necessary for a valid trajectory. This solution
subgraph might have lot of disconnected (graph) components, each component
satisfying the criteria for an Eulerian path (Eiselt, 1995b [8]), i.e. each
component is traversable by one agent. A feasible solution to this problem
will have at least $|\mathcal{K}|$ components in the solution subgraph.

$\displaystyle(\mathcal{B}X)_{vkl}=\begin{cases}-1,&\textit{if }v=s\ \&\
l=0;\forall k\in\mathcal{K}\\}\\\ 1,&\textit{if }v=t\ \&\ l=0;\forall
k\in\mathcal{K}\\}\\\ 0,&\textit{otherwise}\end{cases}$ (4)

The service constraints ensure that every service arc is visited once by at
least one agent. Thus every service arc is selected once in one of the
components of the resulting solution subgraph.

$\displaystyle\sum_{k,l}X_{a_{*}kl}=1,\ \forall a_{*}\in A_{*}$ (5)

The running-time constraints given by Equation (6) ensure that time spent by
an agent in any arc is at least as large as the time required to traverse the
underlying track section. Note that running-time for same track section varies
with choice of agents as well as the type of arc representation (service or
deadhead). If an arc is not occupied, the constraint is relaxed by restricting
the time difference between the two vertices of this arc by a large number.
For similar formulations in literature, this number $\tau$ is chosen very
large, generally equal to the service time limit, if available (Peng, 2013
[13]), or maximum of either the largest time entry in the unavailability list
or the sum of all running-time values over all arcs. Since the parameter
$\tau$ gets discarded on analytical simplification of the presented algorithm,
the value for $\tau$ is chosen as infinity.

$\displaystyle\tau X-\mathcal{B}^{\top}Y\leq\tau\mathbbm{1}-\mathcal{W}$ (6)

$\mathcal{B}^{\top}Y$ represents the difference in running-time, which must be
no smaller than $\mathcal{W}$, if corresponding arcs are occupied i.e.
$X_{akl}=1$. Note, ‘$\mathbbm{1}$’ and ‘$0$’ represent column vectors with
each entry being ones and zeros respectively.

The last set of constraints restrict the temporal values of $Y$ for each of
the arcs of directed multi-graph $a\in A$ based on its list of
unavailabilities $Z_{a}$. Construction of the either-or constraints, given by
Equation (7) with $v$ as the head vertex of arc $a$ (arc $a$ enters at vertex
$v$), using unavailability data is illustrated in Figure 4. Although the set
can be made convex using dummy $0$-$1$ variables, the presented algorithm
directly tackles this type of non-convex constraints, thus reducing the search
space. Notice that the size of unavailability list can be very large, however
with $\tau$ being infinity, this proposed algorithm can tackle infinite number
of unavailability constraints (say periodic unavailabilities that are
generated algorithmically).

$\displaystyle either,$ $\displaystyle\tau
X_{akl}+Y_{vkl^{\prime}}\leq\tau+\underline{\omega}-\mathcal{W}_{akl}$ (7)
$\displaystyle or,$ $\displaystyle\tau
X_{akl}-Y_{vkl^{\prime}}\leq\tau-\overline{\omega}$
$\forall\;(\underline{\omega},\overline{\omega})\in
Z_{a},(a,k,l)\in\mathcal{A},\mbox{ for any }0\leq l,l^{\prime}\leq\mathcal{L}$

$(a_{1},k,l)$ $(a_{2},k,l)$ $(a_{3},k,l)$ spacetime11 am12 pm1 pm2 pm3 pm4
pm$\underline{\omega}$$\overline{\omega}$

Figure 4: Three deadhead arcs of replicated graph having different
unavailability periods are shown. Based on the cruising time data of the agent
$k$ and the unavailability data of the arcs, infeasible regions (red) for each
of the tail vertices are computed. This constraint guarantees that for any
arc, if an agent leaves the tail vertex outside highlighted (red) region, then
the agent can reach the other end without entering the unavailable (blue) zone
and then wait as long as required.

For algorithmic advantages, the unavailability lists $Z_{a}$’s are sorted in
chronological order, which is also a practically reasonable format of
available dataset. Each element $z\in Z_{a}$ contributes a pair of
constraints, as shown in Equation (7). For simplicity in referencing, we name
the constraints of lower time limit for each unavailability period $z$ as
l-cuts while the counterpart is named as u-cuts. All such unavailability
constraints together are denoted by $Z$.

#### Master problem

The optimization problem given by Equations (3), (4), (5), (6) and (7) can be
written in a compact form as follows:

$\displaystyle\underset{X,Y}{\text{min}}$
$\displaystyle\mathcal{W}^{\top}X+\mathcal{C}^{\top}Y;$
$\displaystyle\mathcal{W}\in\mathbb{R}^{|\mathcal{A}|},\
\mathcal{C}\in\mathbb{R}^{|\mathcal{V}|}$ s.t. $\displaystyle
A_{31}X-{\mathcal{B}}^{\top}Y\leq b_{3};$
$\displaystyle\mathcal{A}_{31}\in\mathbb{R}^{|\mathcal{A}|\times|\mathcal{A}|},\
\mathcal{B}\in{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\\{-1,0,1\\}}^{|\mathcal{V}|\times|\mathcal{A}|},b_{3}\in\mathbb{R}^{|\mathcal{A}|}$
$\displaystyle A_{40}X+A_{41}Y\leq b_{41};$
$\displaystyle\mathcal{A}_{40}\in\mathbb{R}^{|\mathcal{Z}|\times|\mathcal{A}|},\
\mathcal{A}_{41}\in\mathbb{R}^{|\mathcal{Z}|\times|\mathcal{V}|},b_{41}\in\mathbb{R}^{|\mathcal{Z}|}$
$\displaystyle A_{40}X+A_{42}Y\leq b_{42};$
$\displaystyle\mathcal{A}_{42}\in\mathbb{R}^{|\mathcal{Z}|\times|\mathcal{V}|},b_{42}\in\mathbb{R}^{|\mathcal{Z}|}$
$\displaystyle X\in S_{X},Y\geq 0;$ $\displaystyle
Y\in\mathbb{R}^{|\mathcal{V}|}$ $\displaystyle where,$ $\displaystyle
S_{X}:=\\{X\in\\{0,1\\}^{|\mathcal{A}|}|\mathcal{B}X=b_{1},A_{2}X=\mathbbm{1}\\}$
$\displaystyle A_{31}=\tau\mathcal{I}$
$\displaystyle\mathcal{I}\rightarrow\textit{ identity matrix of suitable
size}$ $\displaystyle b_{3}=\tau\mathbbm{1}-\mathcal{W}$ $\displaystyle
A_{40},A_{41},\ \&\ A_{42}\textit{ are matrix
representation{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}s}}\textit{
for unavailabilities}$

Next section describes the proposed algorithm in detail, followed by
simulation results to validate the work.

##  IV Algorithm

The routing problem formulated in the previous section (Section III), has
binary spatial variables and linear temporal variables. Benders’ decomposition
is applied to this problem to separate the integer and linear decision
variables, such that the sub-problem becomes a linear programming (LP)
problem. The final algorithm is derived from Benders’ decomposition strategy
with modifications to tackle the LP sub-problem with disjunctive hyper-cube
constraints (defined in Section III) in polynomial time. This exact algorithm
has the NP-Hard part embedded in the spatial portion of the problem, which is
equivalent to a multi-agent rural postman problem (RPP) without unavailability
constraints. The final algorithm is iterative that solves the spatial and
temporal part separately. This work doesn’t indulge into completion,
existence, and convergence arguments which is well established in literature
for Benders’ decomposition. Instead this work shows equivalence relation with
the Benders steps such that disjunctive convex hyper-cube constraints are
tackled, and time of convergence is improved.

### Benders’ decomposition on MTRPP

Benders decomposition decomposes the problem at hand to a master and one or
more sub-problems, which improves computation time (see Section II-B).

On applying Benders’ decomposition directly on the master problem III, the
following sub-problem IV is achieved:

$\displaystyle\underset{Y}{\text{min}}$ $\displaystyle\mathcal{C}^{\top}Y$
s.t. $\displaystyle-{\mathcal{B}}^{\top}Y\leq b_{3}-A_{31}\bar{X}_{i}$
$\displaystyle A_{41}Y\leq b_{41}-A_{40}\bar{X}_{i}$ $\displaystyle
A_{42}Y\leq b_{42}-A_{40}\bar{X}_{i}$ $\displaystyle Y\geq 0$

where $\bar{X}_{i}$ ($\in S_{X}$) is a spatial solution, obtained at $i^{th}$
iteration.

For the sub-problem IV, the dual maximizer search problem IV and the dual
feasibility search problem IV are represented as follows:

$\displaystyle\underset{\lambda}{\text{max}}$
$\displaystyle\tilde{b}_{3}^{\top}\lambda+\tilde{b}_{41}^{\top}\gamma_{1}+\tilde{b}_{42}^{\top}\gamma_{2}$
s.t. $\displaystyle{\mathcal{B}}\lambda-
A_{41}^{\top}\gamma_{1}-A_{42}^{\top}\gamma_{2}\leq\mathcal{C}$
$\displaystyle\lambda,\gamma_{1},\gamma_{2}\geq 0$

$\displaystyle\underset{\lambda}{\text{max}}$
$\displaystyle\tilde{b}_{3}^{\top}\lambda+\tilde{b}_{41}^{\top}\gamma_{1}+\tilde{b}_{42}^{\top}\gamma_{2}$
s.t. $\displaystyle{\mathcal{B}}\lambda-
A_{41}^{\top}\gamma_{1}-A_{42}^{\top}\gamma_{2}\leq 0$
$\displaystyle\lambda,\gamma_{1},\gamma_{2}\geq 0$ $\displaystyle where,$
$\displaystyle\tilde{b}_{3}=A_{31}\bar{X}_{i}-b_{3},\
\tilde{b}_{41}=A_{41}\bar{X}_{i}-b_{41},\textit{ and}$
$\displaystyle\tilde{b}_{42}=A_{42}\bar{X}_{i}-b_{42}$

The reduced master problem IV of $i^{th}$ iteration contains only spatial
decision variables. $\bar{S}_{i}$ collects the dual optimizers obtained from
the dual sub-problem IV, obtained till the $i^{th}$ iteration. Similarly,
$\hat{S}_{i}$ collects all unbounded extreme rays obtained from dual sub-
problem (IV), till the $i^{th}$ iteration. Now, the reduced master problem IV
is given as follows:

$\displaystyle\underset{z,X}{\text{min}}$ $\displaystyle z$ s.t.
$\displaystyle z\geq\mathcal{W}^{\top}X$ $\displaystyle
z\geq\bar{\lambda}^{\top}(A_{31}X-b_{3})\ \
\forall\bar{\lambda}\in\bar{S}_{i}$
$\displaystyle\hat{\lambda}^{\top}(A_{31}X-b_{3})\geq 0\ \
\forall\hat{\lambda}\in\hat{S}_{i}$
$\displaystyle{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}X\in
S_{X}\textit{ and }z\geq 0}$

The constraints constructed using the sets $\bar{S}_{i}$ and $\hat{S}_{i}$ are
obtained from the dual sub-problem IV and IV in the usual way, using Benders’
decomposition steps, as follows:

  * •

If the spatial solution $\bar{X}_{j}$ ($1\leq j\leq i-1$) is feasible for
master problem, then the cost-cuts are derived from the dual maximizer
$\bar{\lambda}_{j}$, such that, at $i^{th}$ iteration,
$z\geq\mathcal{W}^{\top}\bar{X}_{j}+\mathcal{C}^{\top}\bar{Y}_{j}$ (the cost
produced by III), $\forall j\in\\{1,\dots,i-1\\}$.

  * •

If the spatial solution $\bar{X}_{j}$ ($1\leq j\leq i-1$) is infeasible for
master problem, then the feasibility-cuts are derived from the extreme rays
$\hat{\lambda}_{j}$, such that this solution is infeasible for reduced master
problem of further iterations.

  * •

The sub-problem is a minimization problem with positive cost coefficients
($\mathcal{C}\geq 0$) and non-negative variables. This bounds the cost by
zero, from below, and hence unbounded scenarios for IV are irrelevant.

Note that, since unavailability constraints form disjunctive hyper-cubes i.e.
the feasible region for IV is empty, there is _no feasible solution_.
Consequently, the above steps cannot be utilized, unless the either-or
constraints in Equation (7) are modeled in a different way, e.g. using dummy
binary variables.

### Proposed algorithm

The spatial part of the problem III is inherited in the reduced master problem
IV, which is evaluated once per iteration. Similar to Benders’ decomposition
algorithm, IV outputs a lower bound as minimum cost, alongside a spatial
solution $\bar{X}_{i}$ which is feasible with respect to the flow and service
constraints defined in set $S_{X}$. This reduced problem IV is given as

$\displaystyle\underset{z,X,d_{1},\dots,d_{\mathcal{K}}}{\text{min}}$
$\displaystyle z+\sum_{k\in\mathcal{K}}d_{k}$ s.t. $\displaystyle z\geq
2\mathcal{W}^{\top}X$ $\displaystyle d_{k}\geq\lambda^{\top}X\ \
\forall\lambda\in\tilde{S}_{k,i},k\in\mathcal{K}$
$\displaystyle\hat{\lambda}^{\top}(A_{31}X-b_{3})\geq 0\ \
\forall\hat{\lambda}\in\hat{S}_{i}$ $\displaystyle X\in S_{X}\textit{ and
}z,d_{k}\geq 0$

where $d_{k}$ approximates cost due to delay/waiting  for each
$k\in\mathcal{K}$, and $z$ approximates the cost due to movement. Further
details on variables $z$ and $\\{d_{k};\ \forall k\in\mathcal{K}\\}$ will be
described as we develop the constraints of IV, gradually in this section.

The first two constraints in IV lower bounds the variables $z$ and $d_{k}$,
and hence the cost. These constraints are introduced only if $\bar{X}_{j}$ and
$\bar{Y}_{j}$, for any previous iteration (i.e. $j\leq{1,\dots,i-1}$), is
feasible for the master problem. The sets $\tilde{S}_{k,i}$ are determined
from a polynomial time algorithm proposed in Theorem IV.5 to store the
constraint coefficients for each agent and any suitable previous iteration. On
the other hand, if the solution of previous iteration is infeasible for the
master problem then Lemma IV.2 motivates a simple search problem to determine
$\hat{\lambda}$, which are added to the set $\hat{\mathcal{S}}_{i}$ to
iteratively remove infeasible solutions from search space. In this work, the
constraints that lower bounds the cost of IV are termed as cost-cuts, while
the constraints that reduce the initial search space given by $S_{X}$ are
termed as feasibility-cuts.

In order to show the relation between our reduced master problem and sub-
problems, with Benders’ description, a theory is first established for the
temporal graph, followed by analysis and reduction of the Benders’ constraint
representations to the form given in IV.

A temporal graph-agent tuple $(G,Z,R)$ is said to be well-defined if and only
if:

  * •

every arc in $G$ is available for at least the corresponding running-time data
$R$ for any agent, in the period $[0,\infty]$,

  * •

availability period222time gap between two unavailability periods, derived
from $Z$, for every arc in $G$ is larger than the running-time data of that
arc for any of the agents (otherwise it is equivalent to being unavailable for
a subset of the agents), and

  * •

the deadhead graph $G^{\prime}:=(V,A\backslash A_{*},B^{\prime})$, must be
strongly connected i.e. there exist a vertex-disjoint path between any two
nodes in the deadhead graph.

Since replicated graph $\mathcal{G}$ is an assimilation of agent and network
data, the temporal graph $(\mathcal{G},Z)$ is also well-defined iff temporal
graph-agent tuple $(G,Z,R)$ is well-defined in the same interval.

In this work, a less conservative definition for connectedness is introduced:
a temporal graph $(\mathcal{G},Z)$ is called $T_{p}$-periodically connected if
every arc in the parent deadhead graph $G^{\prime}$ is available in
$T_{p}$-period interval. Note that union of all instances in any $T_{p}$
period will produce the replicated graph $\mathcal{G}$.

###### Theorem IV.1

If a temporal graph $(\mathcal{G},Z)$ is well-defined and $T_{p}$-periodically
connected
($T_{p}{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\geq
0}$), with cost of per unit waiting time being at most as large as cost of
equivalent deadhead traversal time (i.e. $\beta\geq 0$ in the cost function),
and every vertex having indefinite waiting permit; then the replicated graph
$\mathcal{G}$ with $|A_{*}|+1$ layers, contains an optimal solution.

Proof.  Every layer of replicated graph $\mathcal{G}$ is constructed from
deadhead graph $G^{\prime}$, hence we can map paths in $G^{\prime}$ to an
equivalent path in $\mathcal{G}$. For simplicity in notations, a sequence is
constructed in $G^{\prime}$ as
$\\{v_{a_{1}},a_{1},v_{a_{2}},a_{2},\dots,a_{r},v_{a_{r+1}}\\}$, with vertices
$v_{a_{1}},v_{a_{2}},\dots,v_{a_{r+1}}\in V$ and arcs
$a_{1},a_{2},\dots,a_{r}\in A\backslash A_{*}$. Let $a_{r+1}\in A_{*}$ be an
arbitrary choice of service arc, and this sequence is assumed to be a path
between depot vertex ($v_{a_{1}}\in D_{k}$) of $k^{th}$ agent, and tail vertex
of arc $a_{r+1}$ (denoted as $v_{a_{r}}$). tail vertices are represented with
$u$ Since every arc is available (implying traversable due to
$(\mathcal{G},Z)$ being well-defined) in any $T_{p}$ time interval, feasible
time-stamps can be assigned to all vertices of the path making it a feasible
trajectory with respect to unavailabilities. If this path has a cycle, then
all arcs of the cycle can be contracted to form a single vertex whose time-
stamp is the maximum among all time-stamps in the cycle, thus preserving
feasibility. Hence, the path
$P(:=\\{v_{a_{1}},a_{1},v_{a_{2}},a_{2},\dots,a_{r-1},v_{a_{r}}\\})$ can be
reduced to a simple vertex disjoint path preferring waiting instead of
deadhead cyclic traversal. Since cost of per unit waiting time is at most the
cost of per unit deadhead traversal time, there exists such a vertex-disjoint
optimal path between any two vertices in the deadhead graph $G^{\prime}$ i.e.
in one layer of replicated graph $\mathcal{G}$. Therefore given the order of
attending the service arcs in an optimal solution, one can construct an
optimal trajectory by joining the source vertex, the sink vertex, and the
service arcs by vertex-disjoint paths; such that each path can be contained in
a single layer of the replicated graph. In other words, an optimal solutions
is always present in this graph with finite number of layers. ❒

###### Corollary IV.1.1

An upper bound on the inspection time $T_{I}$ is $T_{p}(|A_{*}|+1)d$, where
$d$ is the diameter of the deadhead graph $G^{\prime}$.

Proof.  define vertex-disjoint path in section II; also what do you mean by
occupied arcs in this proof? By the definition of diameter, $d$ is the number
of arcs constituting the largest vertex-disjoint path between any two vertices
in $G^{\prime}$. Theorem IV.1 permits only vertex-disjoint path in each layer
of replicated graph $\mathcal{G}$, thus the number of arcs forming such a path
is bounded by $d$. Note that the time required to traverse each arc is at most
$T_{p}$, given the fact that the temporal graph $(\mathcal{G},\mathcal{Z})$ is
$T_{p}$-periodically connected. There are $|A_{*}|+1$ layers in the replicated
graph, thus the optimal time for completing the inspection task by any agent
is bounded by $T_{p}(|A_{*}|+1)d$. ❒

###### Remark 1

Removal of unavailability constraints doesn’t alter the spatial feasible
region given by $S_{X}$.

As a consequence of Remark 1, for any spatial solution of $i^{th}$ iteration
$\bar{X}_{i}$ one can determine a feasible temporal solution $Y_{i}$ which
satisfies the unavailability constraints, even if it doesn’t satisfy the
running-time constraints. In other words, all feasibility-cuts (constructed
from unbounded extreme rays; see Sec-II-B) are derived from the sub-problem
involving running-time constraints only.

###### Lemma IV.2

In case of empty unavailability list, each unbounded extreme ray in IV
corresponds to a unique cycle in the solution subgraph given by
$\bar{X}_{i}\in S_{X}$.

Proof.  In case of empty unavailability list, dual variables $\gamma_{1}$ and
$\gamma_{2}$ can be removed from the dual sub-problem IV. Moreover, from the
master problem III, described in Section III, $A_{31}$ and $b_{3}$ can be
replaced by $\tau\mathcal{I}$ and $\tau\mathbbm{1}-\mathcal{W}$ respectively.
This reduces IV to the following linear program:

$\displaystyle\underset{\lambda}{\text{max}}$
$\displaystyle(\tau\bar{X}_{i}-\tau\mathbbm{1}+\mathcal{W})^{\top}\lambda$
s.t. $\displaystyle{\mathcal{B}}\lambda\leq 0$ $\displaystyle\lambda\geq 0$

Note that, IV is only evaluated when its primal counterpart is infeasible i.e.
there exist an unbounded ray in the dual feasible region.

The proof begins by establishing a relation between the dual vectors $\lambda$
and flow in the graph given by incidence matrix $\mathcal{B}$ (which is
representation of the replicated graph $\mathcal{G}$). Each entry in $\lambda$
corresponds to one column of matrix $\mathcal{B}$ and hence one arc of graph
$\mathcal{R}$. The value of row in vector $\lambda$ can be considered a flow
value in corresponding arc; thus the equation $\mathcal{B}\lambda=0$ becomes a
flow conservation constraint (similar to Equation (4)). Since this product
$\mathcal{B}\lambda$ is evaluated against $0$ vector, the sum of incoming and
outgoing flow values at any vertex must be zero.

The if part of this lemma is proved using the following claims:

  * •

Claim 1: the sets $\\{\lambda\geq 0|\mathcal{B}\lambda\leq 0\\}$ and
$\\{\lambda\geq 0|\mathcal{B}\lambda=0\\}$ represent the same region,

  * •

Claim 2: feasible region represented by the set $\\{\lambda\geq
0|\mathcal{B}\lambda=0\\}$ is a cone,

  * •

Claim 3: any $\lambda\in\\{\lambda\geq 0|\mathcal{B}\lambda\leq 0\\}$, that
does not correspond to a cycle, is not an extreme ray.

Claim 1 is evident, if the vector $\mathcal{B}\lambda$ cannot have identically
non-positive entry for any $\lambda\geq 0$. This statement can be argued by
contradiction using two exhaustive cases: firstly, flow cannot originate at a
vertex as this implies positive right hand side corresponding to that vertex
(or row); and secondly, cyclic flow having imbalance implies there is a source
and sink for the surplus value, making atleast one entry positive and one
negative in the right hand side.

Claim 2 is true because any feasible vector $\lambda$ can be scaled by
$\alpha\in[0,\infty)$ and it still remains feasible. As a consequence, the
entire feasible region can be described using conic combination of extreme
rays (see Section II-A).

For claim 3, the decomposition property of flow is utilized. Consider
$\lambda\in\\{\lambda\geq 0|\mathcal{B}\lambda\leq 0\\}$, and constuct a sub-
graph having only those arcs that have non-zero entry in $\lambda$. If the in-
degree or out-degree for any vertex in this sub-graph is greater than $1$,
then neglecting all except one direction, selects a cycle of the graph (it
cannot be acyclic as our first claim holds). The minimum flow value of this
cycle (say $\lambda_{1}$) can be subtracted from all the flow values of this
cycle, thus reducing flow of some arcs to zero and eliminating them in the
process. Since both $\lambda$ and $\lambda_{1}$ are in null space of
$\mathcal{B}$, the resulting flow (say $\lambda_{2}:=\lambda-\lambda_{1}$)
will also belong to null space of $\mathcal{B}$, aside from being non-
negative. Hence, any feasible $\lambda$, that doesn’t correspond to a cycle,
can be represented using sum of two other vectors from the cone, implying that
its not an extreme ray. In other words sub-graph corresponding to an extreme
ray must be a cycle.

The other part of the proof can be established by showing that a cycle is not
only non-decomposable, but also related to a single valued positive flow that
satisfies flow constraints ($\mathcal{B}\lambda=0$). In other words,
$\lambda$, corresponding to any cycle, cannot be represented as a sum of any
two non-collinear elements of the cone, and hence it is an extreme ray.

Note that if an extreme ray, denoted by $\lambda$, relates to a cycle in the
solution subgraph given by $\bar{X}_{i}$, then the dual cost is unbounded in
that direction i.e. the dual cost
$(\tau\bar{X}_{i}-\tau\mathbbm{1}+\mathcal{W})^{\top}(\alpha\lambda)\rightarrow+\infty$,
as $\alpha\rightarrow+\infty$. Such unbounded extreme rays are binary $0$-$1$
vectors denoted by $\hat{\lambda}$. ❒

From Lemma IV.2, an equivalent representation of feasibility-cut is derived.
Equation (8) renders this equivalent representation where $\hat{\lambda}$’s
are unbounded extreme rays corresponding to cycles in solution subgraph of
$\bar{X}_{i}$, stored in iteratively increasing feasibility-cut set,
$\hat{S}_{i}=\hat{S}_{i-1}\cup\\{\textit{all }\hat{\lambda}\textit{ of
}i^{th}\textit{ iteration}\\}$; initialized by null set,
$\hat{S}_{0}=\\{\hbox{}\\}$. This is used as a constraint in the proposed
reduced master problem IV.

$\displaystyle\hat{\lambda}^{\top}X\leq(\mathbbm{1}^{\top}\hat{\lambda}-1)$
(8)

###### Remark 2

All cycles are restricted to only one layer of a component in replicated graph
$\mathcal{G}$. Also, each cycle has a counterpart in every layer of all
components in $\mathcal{G}$ due to its repetitive structure.

###### Corollary IV.2.1

In case of empty unavailability list, set of feasibility-cut
$\hat{\mathcal{S}}_{i}$ is always empty.

Proof.  From the argument presented in Theorem IV.1, any cycle, in any layer
of the replicated graph, can be collapsed to form a vertex. Contrary to the
method used in Theorem IV.1, the vertex is assigned the minimum time-stamp
among all time-stamps in the cycle. This yields a feasible solution of lower
cost. Again, since there are no inter-layer or inter-component cycles possible
(by Remark 2), solution of IV will never have cycles. ❒

All infeasibility scenarios in running time are eliminated by feasibility-
cuts, and the unavailability constraints don’t contribute to infeasibility
(see Remark 1). Hence, the cost-cuts are determined with an assumption that
there exists feasible temporal solution $Y$, given $\bar{X}_{i}$ (solution of
IV). Optimal temporal solution for $i^{th}$ iteration $\bar{Y}_{i}$ is
determined from the dual sub-problem IV in two steps: first considering no
unavailability case and solving for $\bar{Y}_{i}^{\prime}$ (Remark 3); then
increasing $\bar{Y}_{i}^{\prime}$ minimally to achieve feasibility with
respect to unavailability constraints (Theorem IV.5). This minimal increase in
temporal values results in optimal temporal solution of $i^{th}$ iteration,
$\bar{Y}_{i}$.

###### Lemma IV.3

In case of no unavailability, if the dual sub-problem IV, given
$\bar{X}_{i}\in S_{X}$, has a finite maximum then the dual maximizer
$\bar{\lambda}_{i}$ equals $\bar{X}_{i}$.

Proof.  The finite maximizer search problem is derived from IV, thus _no
unavailability_ scenario changes the feasible set in IV to $\\{\lambda\geq
0|\mathcal{B}\lambda\leq\mathcal{C}\\}$; while the cost vector remains same.
Since it is assumed that this optimization problem has a finite maximum, by
Lemma IV.2, the solution subgraph doesn’t contain cycles. Also, as discussed
earlier, positive cost is only possible if non-zero entires of $\lambda$ are
related to arcs of solution subgraph. Such a choice of $\lambda$ restricts the
$\lambda$-based subgraph to vertex-disjoint paths. Moreover, the constraint
set restricts the value of $\lambda$, at vertex $t\in E_{k}$, to be at most
$(\mathcal{C})_{tk0}$. By backtracking the sink vertices along the vertex-
disjoint paths, it is found out that all non-zero entires of $\lambda$ are at
most $(\mathcal{C})_{tk0}$.

Observe that $\bar{X}_{i}$ is an element of the set $S_{X}$, which implies
$(\mathcal{B}\lambda)_{tk0}\leq(\mathcal{C})_{tk0}:=1$ for all sink vertices
in graph $\mathcal{G}$, and $(\mathcal{B}\lambda)_{vkl}\leq 0$ for all other
vertices. Hence the choice of maximizer $\bar{\lambda}_{i}=\bar{X}_{i}$
produces maximum cost by selecting all arcs in the solution subgraph. ❒

The cost-cuts can be equivalently represented by Equation (9), where
$\bar{\lambda}_{i}$’s are picked from the iteratively increasing cost-cut set,
$\bar{S}_{i}=\bar{S}_{i-1}\cup\\{\bar{\lambda}_{i}\\}$;
($\bar{S}_{0}=\\{\phi\\}$).

$\displaystyle\\{z\geq\bar{\lambda}_{i}^{\top}(A_{31}X-b_{3}),z\geq 0\\}$ (9)
$\displaystyle\iff$ $\displaystyle\\{z\geq\bar{\lambda}_{i}^{\top}(\tau
X-\tau\mathbbm{1})+\bar{\lambda}_{i}^{\top}\mathcal{W},z\geq 0\\}$
$\displaystyle\iff$
$\displaystyle\\{z=\bar{\lambda}_{i}^{\top}\mathcal{W}\textit{ if
$X=\bar{X}_{i}$, }0\textit{ otherwise}\\}$

Note that the feasibility-cut in Equation (8) and the cost-cut in Equation (9)
are independent of $\tau$.

Based on above observations, following changes can be made in the standard
Benders decomposition based algorithm to transform the LP sub-problem (with
convex constraints) to a polynomial-time search:

  * •

Find all cycles in the $i^{th}$ solution subgraph using depth first search
(DFS) based greedy algorithm. The feasibility-cuts are expressed using the
cycle data, as given in Equation (8).

  * •

If there are no cycles, the cost-cuts as given in Equation (9) are introduced
to IV. This set of constraints (all cost-cuts along with
$z\geq\mathcal{W}^{\top}X$) is named as an approximator, more precisely an
under-approximator, i.e. for $\bar{X}_{j}$ ($\forall j|\ 1\leq j\leq i-1$) it
returns the corresponding master cost with minimum temporal values
$\bar{Y}_{i}^{\prime}$ (unavailability scenario is not considered yet), and
for any other $X$ it gives very small value (approaching negative infinity).
If the set of constraints termed as under-approximator (constructed
iteratively) is known apriori i.e. independent of iteration, we call it exact-
approximator. It maybe represented by a single equivalent constraint rather
than a set.

###### Remark 3

Temporal solution based on running-time constraint ($\bar{Y}_{i}^{\prime}$) is
computed in polynomial time, along the path obtained in the solution sub-graph
given by $\bar{X}_{i}$; if IV has a finite maximizer.

The exact approximator, for no unavailability scenario, is given by Equation
(10). For unavailability based scenario, this is a ‘good’ under-approximator
because it reduces the number of iterations, as will be shown in the following
Lemma. This constraint is employed in proposed reduced master problem IV to
approximate the cost due to running-time (no waiting) at the vertices.

$\displaystyle z\geq 2\mathcal{W}^{\top}X$ (10)

###### Lemma IV.4

Any under-approximator is suitable for feasibility-cut, while exact-
approximator being the best choice.

Proof.  Consider $z\geq e_{i}(X)$ be a general representation for the cost-
cuts of Equation (9), where $e_{i}(.)$ is a function of spatial variable
determined from $i^{th}$ iteration. Assume that the function $e_{i}(X)$ has
the following properties:

  * •

$e_{i}(\bar{X}_{j})$ outputs same value as
$\mathcal{W}^{\top}\bar{X}_{j}+\mathcal{C}^{\top}\bar{Y}_{j}$, if previous
iteration $j\in[0,i-1]$ had a feasible dual solution

  * •

$e_{i}(X)$ outputs at most $\mathcal{W}^{\top}X+\mathcal{C}^{\top}Y$, if
$X\neq\bar{X}_{i}$ (here $Y$ is feasible choice for given $X$)

  * •

$e_{i}(X)$ gives any arbitrary value, if iteration $i$ had an unbounded dual
solution

Let $(\bar{X},\bar{Y})$ be an optimal solution for the master problem. If such
an $\bar{X}$ exists, then it cannot be eliminated from the search space by a
feasibility-cut. Also, using the second property of $e_{i}(.)$, its evident
that $e_{i}(\bar{X})$ is at most the optimal master cost, for any $i$.

For any $X\in S_{X}$ the lower bound obtained on solving IV is always
$\mathcal{W}^{\top}X$. Hence, the cardinality of the set $S:=\\{X\in
S_{X}|\mathcal{W}^{\top}X\leq\mathcal{W}^{\top}\bar{X}+\mathcal{C}^{\top}\bar{Y}\\}$
determines the number of iterations in Benders’ algorithm. Once all feasible
choice of $X$ that has cost approximate333lower bound produced by the under-
approximator for given $X$, using set of all cost-cuts of previous iterations
less than the optimal have been evaluated, the minimum value of $z$ attains
the optimal cost. Thus the duality gap becomes zero and the algorithm outputs
an optimal solution of the master problem. If such a $\bar{X}$ doesn’t exist,
then IV returns an empty solution in the first iteration. Consequently, for
exact approximator, this set $S$ which determines the number of iterations,
includes only the infeasible (master) solution with cost less than optimal.
This ensures duality gap to become zero at the exact iteration when the
optimal spatial solution is first evaluated i.e. as soon as a feasible $Y$ is
found for the spatial solution. ❒

The cost contribution due to running-time (and no waiting) is analytically
evaluated into a single constraint, while the unavailability based update can
be determined in polynomial time. Since search of feasibility-cut is also a
polynomial time algorithm, the following can be stated:

###### Theorem IV.5

The sub-problem IV is solvable in polynomial-time.

Proof.  Presenting a polynomial time algorithm for updating temporal solution
$\bar{Y}_{i}^{\prime}$ (obtained in Remark 3) suffices this claim. Here,
$\bar{Y}_{i}^{\prime}$ represents temporal values evaluated using running-time
constraints, without unavailability list based correction. This claim is shown
using the disjunctive hyper-cubes, described by constraints in Equation (7),
where the l-cuts consists of upper bounds while the u-cuts consists of lower
bounds on temporal variables. Since the sub-problem IV is a minimization
problem, only the lower bound given in u-cuts can be tight.

Taking the base case to be the portion of trajectory leading from source
vertex to the depot vertex; this is always feasible for all agents
$k\in\mathcal{K}$. Let a portion of trajectory be directed from $(s,k,0)$ to
$(v^{\prime},k,l^{\prime})$ via $(v,k,l)$ for any arbitrary choice of $k$.
Assume that the portion of trajectory from $(s,k,0)$ to $(v,k,l)$ (excluding
$(v,k,l)$) is feasible for all $k\in\mathcal{K}$ with respect to the
unavailability constraints, and the cost of this portion of trajectory is
minimum. If $(v,k,l)$ is minimally increased by $\Delta Y$ such that $Y_{vkl}$
is feasible, as shown in Figure 5, then the portion of trajectory from
$(s,k,0)$ to $(v,k,l)$ is feasible with respect to unavailability constraints,
and the cost is also minimum by inductive hypothesis. In short, correction of
temporal variables (given the spatial solution) doesn’t affect the temporal
data in variables occurring earlier in the trajectory.

Figure 5: The top figure shows a sample vertex disjoint solution trajectory
for agent $k$. The bottom figure shows that the time value (marked with cicle)
at vertex $Y_{vkl}$ is infeasible because it lies in the unavailable region
given by $[\underline{\omega}_{2}-\mathcal{W}_{akl},\overline{\omega}_{2}]$;
where $\mathcal{W}_{akl}$ is the running time of the outgoing occupied arc
(see Figure 4 to observe the effect of running-time on lower limit of
unavailable region). Hence minimum increment of $\Delta Y$ is required to
correct it. Since $Y_{v^{\prime}kl^{\prime}}$ appears later in this solution
sequence, correction is also reflected in this variable.

Thus this simple greedy algorithm (Algorithm 1: polYcut) corrects the temporal
variable $Y_{i}$ in ascending order in $\mathbf{O}(|Y|)$ steps, as there are
at most $|Y|$ temporal variables. There can only be a finite number of
unavailabilities (say $n$) within the inspection time bound (given in
Corollary IV.1.1) because the availability gaps have finite non-zero lower
bound. In other words, since temporal-graph agent data
$(\mathcal{G},\mathcal{Z})$ is well defined, every availability gap of an arc
is at least as large as the running-time value for that arc. If the smallest
running-time value is given as $\mathcal{W}_{min}$
$(:=\min_{vkl}\\{W_{vkl}\\})$, and the number of availability gaps is at most
$n$, then the product (given as $nW_{min}$) must be upper bounded by the
inspection time ($T_{I}$). Thus the number of comparisons with entries of
unavailability list will be at most $\frac{T_{I}}{\mathcal{W}_{min}}$ per
step. ❒

input : Temporal solution $\bar{Y}_{i}^{\prime}$ of $i^{th}$ iteration
evaluated without unavailability constraints,

: Replicated graph $\mathcal{G}$,

: Vector $O_{Y}$, showing order in which temporal variables appear in solution
path \- all variables of one agent are clubbed together, in increasing order
of the agent numbers,

: Unavailability list for each arc $Z_{a}$, and

: Set of previous cost-cut coefficients due to unavailability
$\tilde{S}_{k,i-1}$

output : Temporal solution $\bar{Y}_{i}$ of $i^{th}$ iteration evaluated with
unavailability constraints, and

: Set of cost-cut coefficients due to unavailability $\tilde{S}_{k,i}$

while  _$O_{Y}$ is not empty_ do

$o_{y}\leftarrow$ Get-smallest-index ($O_{Y}$) ;

// Outputs the index of smallest value in $O_{Y}$ i.e. index of the first
occurring temporal variable (e.g. index of the source vertex of agent 1 is
output of the first step)

$O_{Y}\leftarrow$ _remove $o_{y}$ from $O_{Y}$_;

$y\leftarrow\bar{Y}_{i}^{\prime}[o_{y}]$ ;

$[a,k]\leftarrow$ _determine incoming arc and agent number related to the
vertex of $y$_;

$chk\leftarrow$ Not-in-availability-gap ($y$, $\mathcal{G}$, $Z_{a}$) ;

// Outputs $true$ if $y$ lies in an unavailability period defined in $Z_{a}$
for arc $a$, otherwise $false$

if  _$chk$_ then

$\Delta Y\leftarrow$ _determine minimum increase to achieve feasibility_ ;

$\bar{Y}_{i}^{\prime}\leftarrow$ Correct-values ($\bar{Y}_{i}^{\prime}$,
$O_{Y}$, $k$, $\Delta Y$) ;

// All temporal values in $O_{Y}$ related to agent $k$ is incremented by
$\Delta Y$, and retured as output

$\tilde{S}_{k,i}\leftarrow$ Add-cost-cut ($\tilde{S}_{k,i-1}$, $\Delta Y$,
$o_{y}$, $\mathcal{G}$) ;

// Outputs the constraints given by Equation (11)

end if

end while

$\bar{Y}_{i}\leftarrow\bar{Y}_{i}^{\prime}$

Algorithm 1 polYcut; Polynomial-time algorithm for correcting the temporal
values at $i^{th}$ iteration, and obtaining cost-cuts due to unavailability
constraints simultaneously

The concept of under-approximator is utilized in constructing cost-cuts for
unavailabilities. Following the properties required from an under-
approximator, we propose a cost-cut given in Equation (11), whose coefficient
vectors are stored in $\tilde{S}_{k,i}$ for every feasible Benders iteration.

$\displaystyle d_{k}$ $\displaystyle\geq\Delta
Y\big{(}X_{a^{\prime}kl}-\sum_{a\in\mathcal{N}(P_{vkl})}X_{akl}\big{)}\
\forall k\in\mathcal{K}$ (11) $\displaystyle\textit{and }d_{k}$
$\displaystyle\geq 0\ \forall k\in\mathcal{K}$

$P_{vkl}$ is a path to vertex $(v,k,l)$, having unavailability violated,
$\mathcal{N}(P_{vkl})$ are set of arcs incoming at the path but not at vertex
$(v,k,l)$, $a^{\prime}$ is the occupied arc incoming at the vertex $v$, and
$\Delta Y$ is computed based on the delay caused at vertex $(v,k,l)$ due to
all unavailabilities in $P_{vkl}$.

These cost-cuts due to unavailabilities are constructed such that even if a
part of the path causing unavailability is repeated in any solution,
corresponding delay gets reflected in the cost. Hence it produces a lower
bound for the delay variables $d_{k}$ in the reduced master problem IV.

Finally, the cost-cuts are determined for unavailability scenarios, and
repeated across all components if the agent set is homogeneous i.e. running-
time for all agents is same. Moreover, Remark 2 shows that the feasibility-
cuts are invariant across layers and components, hence they can be repeated
accordingly to save many redundant iterations. The final algorithm for solving
MTRPP is shown in Algorithm 2: mtrpp. For homogeneous agents, i.e. agents with
identical running times and unavailablities; feasibility-cuts and cost-cuts
can be copied accross layers to reduce the number of iterations.

input : Graph-agent data $(G,Z,R)$, where $G=(V,A,B)$;

: Set of service arcs $A_{*}$; and

: Costs $\mathcal{W},\mathcal{C}$;

output : Optimal spatial and temporal solution, $\bar{X}_{i},\bar{Y}_{i}$;

$\mathcal{G}\leftarrow$ _Create replicated graph using $G$ and $R$_;

$Z_{a}\leftarrow$ _Create unavailablity list for each arc $a\in A$_ ;

// $|A|$ lists

$i\leftarrow 0$;

$\mathcal{K}\leftarrow$ _Array $1$ through $\frac{|R|}{|A|}$_ ;

$\hat{S}_{i},S_{\hat{\lambda}}\leftarrow$ _empty_ ;

$\tilde{S}_{k,i}\leftarrow$ _empty for all $k\in\mathcal{K}$_;

while  _$true$_ do

$z,d_{k},\bar{X}_{i}\leftarrow$ _Solve MILP problem IV_ ;

// $|\mathcal{K}|$ variables are generated for storing waiting time in $d_{k}$

$i\leftarrow i+1$;

$S_{\hat{\lambda}}\leftarrow$ _Construct unbounded extreme rays (
$\hat{\lambda}$s), and over-write $S_{\hat{\lambda}}$_ ;

// Depth First Search (DFS) to identify infeasible cycles (Lemma IV.2)

$chk\leftarrow$ Is-Not-Empty ($S_{\hat{\lambda}}$) ;

// Returns a true if the set is not empty

if  _$chk$_ then

$\hat{S}_{i}\leftarrow\hat{S}_{i-1}\cup S_{\hat{\lambda}}$ ;

// Store all feasibility-cuts

else

$\bar{Y}_{i}^{\prime},O_{Y}\leftarrow$ _Compute temporal values and their
order of occurrence in vertex-disjoint path(s) in solution subgraph given by
$\bar{X}_{i}$_ ;

// (Remark 3)

$\bar{Y}_{i}\leftarrow$ polYcut
($\bar{Y}_{i}^{\prime},\mathcal{G},O_{Y},Z_{a},\tilde{S}_{k,i-1}$) // Correct
the temporal solution, and find cost-cuts occuring due to unavailability
constraints, in polynomial time (Theorem IV.1)

end if

if  _$z$ equals $\mathcal{W}^{\top}\bar{X}_{i}+\mathcal{C}^{\top}\bar{Y}_{i}$_
then

_break loop_ // Stop if LB equals UB (a small tolerance maybe used to
compensate for numerical errors)

end if

end while

Algorithm 2 mtrpp; Benders’ decomposition based algorithm to solve Multi-agent
Temporal Rural Postman Problem

Lastly, to ensure that optimal solution for the actual problem statement is
achieved by the proposed methodology, the following claim is proved:

###### Theorem IV.6

Optimal solution of MTRPP can be achieved using the proposed methodology of
replicated graph based formulation and algorithm (see Algorithm 2: mtrpp).

Proof.  All spatial solution from the set $\mathcal{S}_{X}$ are either a path
from source to sink with cycles that are connected or disconnected to the
path. Since Lemma IV.2 eliminates formation of cycles, the feasible solution
set in the proposed master problem formulation III contains only vertex-
disjoint solutions. Thus this methodology finds a vertex-disjoint optimal
solution, represented by arcs and vertices of replicated graph.

Also, by Theorem IV.1, set of vertex disjoint paths, in each layer of the
replicated graph, can be concatenated to form an optimal solution for MTRPP;
without repeating any arc of the directed multi-graph $G$ more than
$|A_{*}|+1$ times. Hence an optimal solution can be found using the proposed
methodology. ❒

###### Remark 4

Vertex disjoint path in replicated graph $\mathcal{G}$ doesn’t mean vertex
disjoint path in $G$.

##  V Simulation examples

###  V-A Case studies

In this work, two simulation case studies are presented to validate the
effectiveness of the proposed algorithm.

The first case is demonstrated using an example picked from article by
(Lannez, 2015 [5]) - a graph with three vertices and six arcs. This example is
slightly different from literature as some side constraints are removed and
details are modified to represent it in the problem setting of this work. This
also supports comparison between Time-Space Network TSN and the proposed
approach. TSN approach is implemented by discretizing time, and expanding
temporal graph to time dimension, thus making it a static graph in time-space
dimension. The approach is detailed in (Peng, 2011 [17]) and (Lannez, 2015
[5]). This approach is also called static expansion of temporal graph
(Michail, 2016 [18]).

The second case is a real application of the proposed work, based on a portion
of Mumbai sub-urban rail network.

#### Ex-I

The first example is a directed multi-graph with three vertices and six arcs,
shown earlier in Figure 2. A TSN based decision graph, constructed from the
same temporal graph-agent tuple $(G,Z,R)$, is illustrated in Figure 6. It’s
replicated graph based representation is shown in Figure 3.

Figure 6: TSN representation of a graph with three vertices and six arcs (see
Figure 2). The x-axis represents time stamps discretized with interval of $1$
minute. The three vertices are shown vertically for each time-stamp (x-axis),
and they are connected by arcs from left to right based on the running-time of
that arc. An arc with running-time $3$ minutes connects vertices at $t$ and
$t+3$. The waiting arcs connect same vertex of consecutive time-stamps, thus
displayed by horizontally drawn arcs. Based on unavailability list (see Table
I, Ex-I, Case I), some arcs are not present at certain times. Arc to the sink
vertex (rightmost vertex) is introduced for every time-stamp such that an
agent can end the trajectory anytime, and not be penalized for waiting.
Finally, the optimal trajectory for Case I is highlighted. (The graph is
generated using Matlab)

Unavailability list for different cases is shown in Table I, along with
comparison results. Notice that the concept of inspection time bound/limit
(see Corollary IV.1.1) is not assumed anywhere in the proposed method, unlike
TSN where algorithm convergence is highly sensitive to the number of discrete
intervals. This is observed by the change in algorithmic runtime of the two
methods as time limit is increased (see Table I, Ex-IIa and Ex-IIb). Time
limit is designer’s choice; for simplicity we have chosen a value using prior
knowledge of the solution to obtain finite number of time-space decision
variables for TSN.

TABLE I: Results comparing algorithm runtime, number of variables and
constraints of TSN and the proposed method

Unavailability list Time limit No. of agents No. of variables No. of
constraints Algorithm runtime Arcs Lower limit ($\underline{\omega}$) Upper
limit ($\overline{\omega}$) for TSN (secs) (minutes) (minutes) TSN Proposeda
TSN Proposedb TSN Proposed Ex-Ia $a_{1}$ 4 5 8 1 65 15L+25I 31 68 0.8 0.7
$a_{2}$ 4 5 $a_{5}$ 7 8 $a_{6}$ 7 8 Ex-Ib $a_{4}$ 4 5 20 1 168 15L+25I 67 68
0.8 0.8 $a_{5}$ 5 14 Ex-Ic $a_{3}$ 4.1 4.9 12 2 190 30L+50I 84 134 0.8 0.9
$a_{3}$ 8.1 10.9 $a_{4}$ 4.1 4.9 $a_{5}$ 5.5 14 Ex-IIa $719$ unavailabilities
222 2 29,532 724L 16,069 2901 255 251 Ex-IIb $719$ unavailabilities 296 2
42,434 \+ 1084I 21,397 2609 251 Ex-IIc $4727$ unavailabilities 1554 2 - - -
260 Ex-IId No unavailability - 2 - - - 1.0 a The subscript $L$ indicates
linear and $I$ integer. b The constraints due to unavailabilities are not
shown these constraints are not used directly (see Algorithm 1).

In Table I, unavailability list for Ex-Ia is also directly picked from the
corresponding article by (Lannez, 2015 [5]), while unavailabilty list for Ex-
Ib and Ex-Ic is constructed to show different trajectory possibilities. For
Ex-Ia, the agent visits service arc $a_{2}$ first, followed by arc $a_{5}$
(see Figure 6); while the sequence is reversed for Ex-Ib i.e. arc $a_{5}$ is
visited first, followed by arc $a_{2}$. In Ex-Ic, optimal solution is improved
due to presence of an extra agent. The performance of TSN approach in Ex-Ic is
better, showing that the improvement is not justified in small graphs with
small unavailability list.

#### Ex-II

A directed multi-graph is constructed to resemble the Kurla-Thane-Vashi-Kurla
sub-urban railway network located in Mumbai, India. Since entire data is not
available with us, following assumptions are made: all stations considered in
this network can park any number of agent for indefinite period of time,
unavailability schedules are periodic, and its possible for an agent to move
between any two stations (while servicing or cruising) in any $74$ minute time
interval i.e. $T_{p}$-periodic connectivity property is guaranteed; where
$T_{p}$ is $74$ minutes. Figure 7 shows the directed multi-graph for the
Kurla-Thane-Vashi-Kurla sub-urban railway. Such periodic schedule that repeats
every $74$ minutes is not actually employed by Mumbai sub-urban railways,
however its not completely fictitious as Dutch has periodic schedule (also
called cyclic timetable) that repeats every hour (Peeters, 2003 [2])

Figure 7: This graph portrays 36 stations all modeled as vertices - indicated
using circles of two sizes. There are 45 deadhead arcs in this model, of which
9 are to be serviced. As shown in the enlarged figure in right, the service
arcs marked in bold have a duplicate deadhead arc in parallel.

The results are presented in Table I, and compared with TSN using four
settings:

  * •

In Ex-IIa, prior knowledge of the total inspection time taken by optimal
solution ($207$ minutes) is assumed. Hence the least integral multiple of
$T_{p}$ that can contain the optimal solution i.e. $3$ times $T_{p}$ ( = $222$
minutes) is considered as time limit for formulating TSN.

  * •

In Ex-IIb, effect on increasing the total inspection time is observed i.e. $4$
times $T_{p}$ (= $296$ minutes)444for single agent scenario, the best deadhead
trajectory based solution (i.e. solution with empty unavailability list) gets
delayed to $301$ minutes on considering unavailabilities is considered for
formulating TSN.

  * •

In Ex-IIc, using inspection time bound derived from Remark III \- diameter of
deadhead graph is $21$ arcs long, hence inspection time bound is $21*74=1554$
minutes.

  * •

In Ex-IId, simulation without unavailability scenario is presented.

Results of TSN for Ex-IIc is not shown in Table I because Matlab responded
with an out of memory error. For the special case (Ex-IId), there are existing
RPP algorithms to provide a comparable, if not better, solution. Since the two
agents are identical, the trajectories can be interchanged if there is an
overtaking occurrence in the schedule.

#### Observation

From the results shown in Table I, it can be observed that the proposed
algorithm has less sensitivity to size of unavailability list compared to TSN.
The most important factor contributing to this improvement is the polynomial
time algorithm (Algorithm 1). The number of variables and constraints don’t
change for the proposed algorithm, however the computation time varies mostly
due to the structure of unavailabilities.

###  V-B Comparison with TDRPP

The comparison study presented in this work is not an exact one because the
underlying data and scenario cannot be made identical with the TDRPP
(Calogiuri, 2019 [4]). However, it gives an insight on performance of the
proposed algorithm relative to algorithms handling similar problems from
literature. Note that the proposed problem setting is more general compared to
TDRPP due to relaxation of FIFO assumption; and traversablity in railway based
temporal graph is not implied by connectedness of the underlying graph.
Moreover, since the presented work deals with servicing (two arc
representation; deadhead-service arc pair) instead of congestion (which can be
represented with one arc), the total number of arcs in each comparison
instance is $|A_{*}|$ more compared to TDRPP. Even though there are these
differences, both problems are temporal variant of the well known rural
postman problem (RPP). Hence the two different model-algorithm pairs (one
proposed in TDRPP article and another proposed in this work) might be adjusted
to suit the other problem setting and improve the final results.

For comparison, all data are picked from the results presented for TDRPP such
that the problem is identical. Three cases of vertex set sizes are considered:
$20$, $40$, and $60$. Further, three arc sizes are put in comparison table
having number of arcs equal $1.2$, $1.6$, and $2$ times the number of vertices
(chosen such that the number of arcs are integers), as shown in Table II. The
arcs are selected in two stages: firstly, a random Hamiltonian tour is chosen
with random running time weights for $|V|$ number of arcs, followed by random
arc generation for the remaining number or arcs; to ensure strongly connected
underlying graph. Three cases of service arc set sizes are considered, given
by $round(\beta|A|)$, where $\beta$ takes values from the set
$\\{0.3,0.5,0.7\\}$. These arcs are added to the existing set, such that the
service time and deadhead time are represented by different arcs. For each
node sizes, arc sizes, and service arc sizes (given by $\beta$); $30$ random
graph instances are generated, and the result is averaged out before
displaying in Table II. If all $30$ instances are not solved to optimality,
another row is added to show the duality gap achieved within a pre-specified
time interval. For the temporal part, a periodic schedule of train is assumed
to traverse in the arcs forming Hamiltonian tour, with randomly generated
traversal times, resulting in unavailabilities occurring in $|V|$ arcs only.

TABLE II: Comparison with TDRPP

TDRPP Proposed $\beta$ $|V|$ $|A|$ OPT NODES TIME GAPi GAPf OPT ITER TIME GAPi
GAPf $0.3$ $20$ $24$ $30$ $11$ $0.40$ $4.37$ $0.00$ $30$ $8$ $0.61$ $10.04$
$0.00$ $32$ $30$ $1,173$ $30.57$ $4.83$ $0.00$ $30$ $18$ $3.22$ $8.26$ $0.00$
$40$ $29$ $835$ $22.81$ $4.26$ $0.00$ $30$ $25$ $10.08$ $11.06$ $0.00$ $-$ $-$
$5.19$ $0.30$ $40$ $48$ $30$ $127$ $3.80$ $4.39$ $0.00$ $29$ $41$ $28.15$
$15.44$ $0.00$ $250$ $>700$ $32.15$ $4.29$ $64$ $27$ $5,212$ $148.00$ $4.49$
$0.00$ $23$ $42$ $51.53$ $17.53$ $0.00$ $-$ $-$ $5.47$ $0.57$ $175$ $>700$
$29.90$ $7.09$ $80$ $26$ $8,415$ $216.75$ $4.25$ $0.00$ $3$ $50$ $238.06$
$9.95$ $0.00$ $-$ $-$ $5.42$ $0.54$ $89$ $>700$ $21.10$ $7.16$ $60$ $72$ $30$
$1,029$ $29.05$ $4.73$ $0.00$ $30$ $68$ $34.23$ $22.99$ $0.00$ $96$ $24$
$1,927$ $57.87$ $4.61$ $0.00$ $2$ $46$ $367.37$ $35.174$ $0.00$ $-$ $-$ $5.78$
$0.51$ $131$ $>1300$ $28.24$ $11.38$ $120$ $28$ $2,513$ $82.47$ $1.94$ $0.00$
$0$ $-$ $-$ $5.57$ $0.43$ $73$ $>1300$ $28.35$ $12.52$ $0.5$ $20$ $24$ $30$
$34$ $0.98$ $3.43$ $0.00$ $30$ $7$ $0.64$ $9.52$ $0.00$ $32$ $30$ $2,721$
$67.72$ $4.17$ $0.00$ $28$ $79$ $54.28$ $14.54$ $0.00$ $-$ $-$ $-$ $-$ $428$
$>700$ $13.49$ $0.89$ $40$ $22$ $22,743$ $614.20$ $3.81$ $0.00$ $25$ $72$
$99.31$ $16.67$ $0.00$ $-$ $-$ $4.21$ $0.78$ $245$ $>700$ $20.93$ $3.21$ $40$
$48$ $30$ $402$ $10.34$ $3.54$ $0.00$ $28$ $73$ $136.4$ $20.00$ $0.00$ $-$ $-$
$-$ $-$ $156$ $>700$ $18.54$ $6.15$ $64$ $12$ $16,729$ $453.01$ $2.97$ $0.00$
$1$ $34$ $89.53$ $20.11$ $0.00$ $-$ $-$ $3.85$ $0.75$ $127$ $>700$ $30.84$
$8.30$ $80$ $9$ $10,158$ $260.85$ $2.89$ $0.00$ $4$ $28$ $390.41$ $26.32$
$0.00$ $-$ $-$ $3.68$ $0.62$ $43$ $>700$ $32.81$ $14.33$ $60$ $72$ $30$
$4,050$ $116.01$ $3.69$ $0.00$ $29$ $43$ $145.24$ $32.615$ $0.00$ $130$
$>1300$ $36.54$ $4.13$ $96$ $9$ $30,815$ $810.22$ $3.27$ $0.00$ $1$ $29$
$610.37$ $30.54$ $0.00$ $-$ $-$ $3.78$ $0.54$ $64$ $>1300$ $35.11$ $13.11$
$120$ $6$ $8,926$ $293.97$ $3.23$ $0.00$ $0$ . . . . $-$ $-$ $3.77$ $0.53$ . .
. . $0.7$ $20$ $24$ $30$ $46$ $1.41$ $2.80$ $0.00$ $30$ $8$ $1.48$ $13.39$
$0.00$ $32$ $30$ $11,420$ $287.06$ $3.10$ $0.00$ $29$ $105$ $116.62$ $16.75$
$0.00$ $359$ $>700$ $12.33$ $1.17$ $40$ $5$ $15,285$ $426.92$ $2.26$ $0.00$
$15$ $106$ $292.02$ $23.92$ $0.00$ $-$ $-$ $2.94$ $1.00$ $193$ $>700$ $20.04$
$3.40$ $40$ $48$ $30$ $889$ $22.72$ $2.63$ $0.00$ $23$ $52$ $204.13$ $2.18$
$0.01$ $41$ $>700$ $33.33$ $13.65$ $64$ $2$ $17,516$ $481.34$ $2.52$ $0.00$
$0$ $-$ $-$ $3.22$ $0.88$ $23$ $>700$ $34.95$ $17.76$ $80$ $2$ $10,641$
$258.96$ $2.57$ $0.00$ $0$ . . . . $-$ $-$ $2.99$ $0.77$ . . . . $60$ $72$
$27$ $21,554$ $554.05$ $2.88$ $0.00$ $0$ . . . . $-$ $-$ $3.17$ $0.13$ . . . .
$96$ $1$ $43,817$ $1029.50$ $2.63$ $0.00$ $0$ . . . . $-$ $-$ $2.99$ $0.80$ .
. . . $120$ $1$ $51,691$ $1283.04$ $1.94$ $0.00$ $0$ . . . . $-$ $-$ $3.05$
$0.81$ . . . . ‘$-$’ data unavailable from the source (Calogiuri, 2019 [4])
how do you classify it as garbage? ‘.’ algorithm failed to generate this data

Data for TDRPP is directly copied from the article by (Calogiuri, 2019 [4]),
and compared with the one generated using the proposed algorithm. The terms in
the table indicate:

OPT

Number of instances (out of 30) that reached the optimal solution within a
pre-specified time ($700$ seconds for vertex size $20$ and $40$, while $1300$
seconds for vertex size $60$),

NODES

Number of nodes explored in the branch-and-bound while solving TDRPP
(Calogiuri, 2019 [4]),

ITER

Average number of iterations (rounded) of the proposed algorithm,

TIME

Average computation time (in seconds) using the indicated method,

GAPi

Average duality gap during the initial phase, and

GAPf

Average duality gap during the final phase.

From the comparison result in Table II, generated using $810$ random graphs,
it can be observed that the proposed algorithm performs well for smaller and
sparse graphs. The time taken per iteration is observed to increase, which is
justified because an Integer Linear Program (ILP) is solved in each iteration
of the proposed algorithm. Matlab uses a branch and bound technique internally
to solve the ILP at each iteration. Hence incorporating branch and bound
technique in the proposed algorithm, such that each iteration is a relaxed
linear programming problem, can further improve the performance results.

For graph with $20$ vertices, the overall performance by the proposed
algorithm is better compared to TDRPP. Its difficult to judge for graphs with
$40$ vertices and small number of service arcs ($\beta\in\\{0.3,0.5\\}$).
However, as $\beta$ is increased to $0.7$, the proposed algorithm struggles
more than TDRPP approach. A slightly decreasing trend is observed in the
performance of proposed algorithm for graphs with $60$ vertices and
$\beta\in\\{0.3,0.5\\}$. For $\beta=0.7$, graph with $60$ vertices demands lot
of storage resources, forcing Matlab to throw out of memory error. Although
our results are generated in a HP Probook (Core i7, 8th generation) and 8 GB
of memory - which is a slightly better system compared to MacBook Pro with a
2.33-GHz Intel Core 2 Duo processor and 4 GB of memory used for TDRPP - yet an
improved performance is expected for smaller and sparser graphs with identical
dataset. Furthermore, TDRPP is solved using CPLEX solver, as compared to
Matlab for proposed method.

##  VI Conclusion

Routing and scheduling problem with time-dependent unavailability scenarios
can be planned optimally using the proposed methodology. With inspection
schedules being sandwiched between train schedules, more frequent monitoring
of railway track is possible, thus improving reliability of the railway
services. Moreover, in a general setting, this work contributes an approach to
obtain exact solution of multi-agent temporal RPP. To consider practical
scenarios, application oriented MILP formulation of this inspection planning
problem with side constraints can be adapted from work presented by (Peng,
2011 [17]).

The backbone of this work is established by Theorem IV.1 which shows existence
of optimal result in replicated graph based formulation. In particular, the
assumptions of well-defined and $T_{p}$-periodic connectivity in the graph-
agent data, dictate the usability of presented results. The clause, insisting
the graph-agent data to be well-defined, governs compatibility between network
and agents in the optimization problem i.e. an agent that moves too slow (e.g.
ultrasonic inspection vehicles) will not be able to traverse most tracks of a
sub-urban railway network, in any availability gap, during peak hours. Having
a well-defined graph-agent data, or producing one before attempting the
optimization problem, is a wise choice to reduce computational burden.
Similarly, most schedules are periodic, hence $T_{p}$-periodic connectivity is
an intrinsic property in many applications. The theorem has another assumption
on indefinite waiting permit at vertices which restricts our result to cases
with small number of agents in many applications, like railways. Since large
agent swarm is not practical in railways, and considering the fact that such
problems are solved regionally by railway undertakings, portions of railway
networks can be modeled with small and sparser graphs, making our results
suitable for regional implementation. Swarm cases on large and dense graphs
can be considered a future work to inspect complex railway networks on large
scale.

The proposed work is more suited to on-demand routing and scheduling tasks,
however it can also be applied for periodic inspection in a sub-optimal
manner, while ensuring no disruption of the train schedule. This novel
methodology, supported by simulation and comparison study, serves as a bridge
between the theory developed on temporal graphs and its practical application
in railway track inspection routing and scheduling problems. This work also
addresses MTRPP which adds exact solution, developed from Benders’
decomposition, for yet another variant of rural postman problems.

Apart from research on optimizing agent swarms, an application specific
optimization, to address challenges arising from practical side-constraints,
would also be another research direction for future. Moreover, addressing the
limitation, arising due to headway and overtaking scenarios in heterogeneous
agents, will further strengthen the practicality of this work. Research
directed towards robust optimization may be attempted to handle yet another
challenge of schedule uncertainty.

## Acknowledgment

One of the authors would like to thank IITB-Monash Research Academy for
providing with scholarship for this research.

## References

  * [1] T. Schlechte, “Railway track allocation: Models and algorithms,” Ph.D. dissertation, der Technischen Universit at Berlinzur Erlangung des akademischen Grades, 01 2012.
  * [2] L. Peeters, “Cyclic railway timetable optimization,” Ph.D. dissertation, Erasmus University Rotterdam, Jun. 2003. [Online]. Available: http://hdl.handle.net/1765/429
  * [3] J.-F. Cordeau, G. Ghiani, and E. Guerriero, “Analysis and branch-and-cut algorithm for the time-dependent travelling salesman problem,” _Transportation Science_ , vol. 48, no. 1, pp. 46–58, 2014. [Online]. Available: https://doi.org/10.1287/trsc.1120.0449
  * [4] T. Calogiuri, G. Ghiani, E. Guerriero, and R. Mansini, “A branch-and-bound algorithm for the time-dependent rural postman problem,” _Computers $\&$ Operations Research_, vol. 102, pp. 150 – 157, 2019. [Online]. Available: http://www.sciencedirect.com/science/article/pii/S0305054818302004
  * [5] S. Lannez, C. Artigues, J. Damay, and M. Gendreau, “A railroad maintenance problem solved with a cut and column generation matheuristic,” _Networks_ , vol. 66, no. 1, pp. 40–56, 2015. [Online]. Available: https://onlinelibrary.wiley.com/doi/abs/10.1002/net.21605
  * [6] A. J. Conejo, E. Castillo, R. Mínguez, and R. Garcia-Bertrand, _Decomposition Techniques in Mathematical Programming_. Springer-Verlag Berlin Heidelberg, 2006.
  * [7] H. A. Eiselt, M. Gendreau, and G. Laporte, “Arc routing problems, part i: The chinese postman problem,” _Operations Research_ , vol. 43, no. 2, pp. 231–242, 1995. [Online]. Available: https://doi.org/10.1287/opre.43.2.231
  * [8] H. A. Eiselt, M. Gendreau, and G. Laporte, “Arc routing problems, part ii: The rural postman problem,” _Operations Research_ , vol. 43, no. 3, pp. 399–414, 1995. [Online]. Available: https://doi.org/10.1287/opre.43.3.399
  * [9] S. M. Pour, J. H. Drake, L. S. Ejlertsen, K. M. Rasmussen, and E. K. Burke, “A hybrid constraint programming/mixed integer programming framework for the preventive signaling maintenance crew scheduling problem,” _European Journal of Operational Research_ , vol. 269, no. 1, pp. 341 – 352, 2018. [Online]. Available: http://www.sciencedirect.com/science/article/pii/S0377221717307646
  * [10] G. Budai, D. Huisman, and R. Dekker, “Scheduling preventive railway maintenance activities,” _Journal of the Operational Research Society_ , vol. 57, no. 9, pp. 1035–1044, Sep 2006. [Online]. Available: https://doi.org/10.1057/palgrave.jors.2602085
  * [11] C. Dao, R. Basten, and A. Hartmann, “Maintenance scheduling for railway tracks under limited possession time,” _Journal of Transportation Engineering, Part A: Systems_ , vol. 144, no. 8, p. 04018039, 2018. [Online]. Available: https://ascelibrary.org/doi/abs/10.1061/JTEPBS.0000163
  * [12] A. Higgins, L. Ferreira, and M. Lake, “Scheduling rail track maintenance to minimise overall delays,” in _The 14th International Symposium on Transportation and Traffic Theory_ , Jerusalem, Israel, 1999. [Online]. Available: https://eprints.qut.edu.au/2531/
  * [13] F. Peng, Y. Ouyang, and K. Somani, “Optimal routing and scheduling of periodic inspections in large-scale railroad networks,” _Journal of Rail Transport Planning $\&$ Management_, vol. 3, no. 4, pp. 163 – 171, 2013. [Online]. Available: http://www.sciencedirect.com/science/article/pii/S221097061400016X
  * [14] E. Resendiz, J. M. Hart, and N. Ahuja, “Automated visual inspection of railroad tracks,” _IEEE Transactions on Intelligent Transportation Systems_ , vol. 14, no. 2, pp. 751–760, June 2013.
  * [15] M. Papaelias, M. Lugg, C. Roberts, and C. Davis, “High-speed inspection of rails using acfm techniques,” _NDT $\&$ E International_, vol. 42, no. 4, pp. 328 – 335, 2009. [Online]. Available: http://www.sciencedirect.com/science/article/pii/S0963869508001527
  * [16] G. Nicholson and C. Davis, “Modelling of the response of an ACFM sensor to rail and rail wheel RCF cracks,” _NDT $\&$ E International_, vol. 46, pp. 107 – 114, 2012. [Online]. Available: http://www.sciencedirect.com/science/article/pii/S0963869511001757
  * [17] F. Peng, “Scheduling of track inspection and maintenance activities in railroad networks,” Ph.D. dissertation, Civil Engineering in the Graduate College of the University of Illinois at Urbana-Champaign, 2011. [Online]. Available: http://hdl.handle.net/2142/24517
  * [18] O. Michail, “An introduction to temporal graphs: An algorithmic perspective,” _Internet Mathematics_ , vol. 12, no. 4, pp. 239–280, 2016. [Online]. Available: https://doi.org/10.1080/15427951.2016.1177801

left black intentionally

## Appendix

### Notations

All notations are listed in the order of their occurrence.

Graph notations used in this work:

$G$

($:=(V,A,B)$) A directed multi-graph, where $V$ is the vertex set, $A$ is the
arc set, and $B$ is incidence matrix describing vertex-arc interconnections.
In this work, for notational simplicity, arcs are not represented as tuples,
instead they are represented as single variables (say $a_{1}\in A$), and thus
the relationship between arcs and vertices are expressed using the incidence
matrix.

$G^{\prime}$

($:=(V,A\backslash A_{*},B^{\prime})$) A sub-graph of $G$, called deadhead
sub-graph, that doesn’t have service arcs ($A_{*}$); hence $B^{\prime}$ is an
incidence matrix without the columns for any service arc.

$\mathcal{G}$

($:=(\mathcal{V},\mathcal{A},\mathcal{B},\mathcal{W})$) Replicated graph
constructed from directed multi-graph $G$ and agent running time
characteristics; $\mathcal{V}$ is the vertex set, $\mathcal{A}$ is the arc
set, $\mathcal{B}$ is the incidence matrix, and $\mathcal{W}$ denotes column
vector of arc weights.

List of variables is given as follows:

$X_{akl}$

binary spatial variable for each arc in $\mathcal{G}$

$X$

Column vector concatenated using spatial variables $X_{akl}$

$X_{i}$

A spatial solution vector of $i^{th}$ iteration

$\bar{X}_{i}$

Spatial optimal solution vector with respect to the reduced master problem of
the $i^{th}$ iteration

$\bar{X}$

Spatial solution vector which is optimal with respect to the main master
problem

$Y_{vkl}$

($\in\mathbb{R}$) temporal variable for each vertex in $\mathcal{G}$

$Y$

Column vector concatenated using temporal variables $Y_{vkl}$

$Y_{i}$

A temporal solution vector at $i^{th}$ iteration

$\bar{Y}_{i}$

Temporal optimal solution vector with respect to the master problem given
$\bar{X}_{i}$; this is also the solution vector of sub-problem at $i^{th}$
iteration

$\bar{Y}_{i}^{\prime}$

Temporal optimal solution vector with respect to the master problem given
$\bar{X}_{i}$, but without unavailability constraints; this is also the
solution of sub-problem of $i^{th}$ iteration with no unavailability

$\bar{Y}$

Temporal solution vector which is optimal with respect to the main master
problem

$\lambda$

Dual variables related to running-time constraints

$\hat{\lambda}$

Dual extreme rays related to cycles in subgraph constructed from $\bar{X}_{i}$

$\bar{\lambda}_{i}$

Dual optimizer related to the sub-problem at $i^{th}$ iteration, without
unavailability constraints

$\gamma_{1}\&\gamma_{2}$

Dual variables related to unavailability constraints

$z$

Variable used to approximate master cost due to running (and not due to
waiting)

$d_{k}$

Variable used to approximate cost due to delay only

Few important notations are listed below:

$A_{*}$

($\subset A$) Set of service arcs - each service arc has an identical arc
(connecting the same set of vertices in same direction) not present in $A_{*}$

$\mathcal{K}$

Set of agents, each agent $k$ ($\in\mathcal{K}$)

$\mathcal{L}$

($:=|A_{*}|+1$) Number of layers

$D_{k}$

($\subset V$) Set of start vertices corresponding to each agent
$k\in\mathcal{K}$

$E_{k}$

($\subset V$) Set of terminating vertices corresponding to each agent
$k\in\mathcal{K}$

$R$

($\in\mathbb{R}^{|\mathcal{K}||A|}$) Vector comprising of running times of
each agent-arc pair

$Z_{a}$

List of unavailabilities per arc $a\in A$

$(G,Z,R)$

Tuple describing graph, unavailabilities and agents running-time; termed as
temporal graph-agent tuple

$i\&j$

$i$ indicates current iteration, while $j$ indicates a previous iteration

$\mathcal{C}$

Column vector of dimension $|\mathcal{V}|\times 1$, denoting weighting factor
on the temporal variables

$(.)_{vkl}$

An entry of a vector related to $v^{th}$ vertex ($\in V$), $k^{th}$ agent and
$l^{th}$ layer; e.g. $(\mathcal{C})_{vkl}$ denotes entries in vector
$\mathcal{C}$ while $(\mathcal{B}X)_{vkl}$ denotes entries in vector formed by
the product of incidence matrix $\mathcal{B}$ and vector $X$

$\beta$

Penalty factor on the cost of deadhead traversals

$\mathbbm{1}$

Column vector with all entries as $1$

$0$

Column vector with all entries as zero

$\mathcal{I}$

Denotes identity matrix of dimension $|\mathcal{A}|\times|\mathcal{A}|$

$e_{i}(.)$

A function of spatial variables, termed as under-approximator

$S_{X}$

Set of spatial solutions that satisfy the flow and service constraints

$\bar{S}_{i}$

Set of feasibility-cuts, formed by dual vectors at the $i^{th}$ iteration

$\hat{S}_{i}$

Set of cost-cuts, formed by dual vectors at the $i^{th}$ iteration, due to
running-time constraints

$\tilde{S}_{k,i}$

Set of cost-cuts, formed by coefficient vectors at the $i^{th}$ iteration, due
to unavailability constraints