# Exact resolution of a simultaneous vehicle routing and crew scheduling
problem in long-haul transport

Mauro Lucci
<EMAIL_ADDRESS>
CONICET - Universidad Nacional de Rosario, Argentina.
Daniel Severín
<EMAIL_ADDRESS>
Universidad Nacional de Rosario, Argentina  Paula Zabala
<EMAIL_ADDRESS>
CONICET - Universidad de Buenos Aires, Argentina

###### Abstract

This work focuses on exact methods for a Simultaneous Vehicle Routing and Crew
Scheduling Problem in long-haul transport. Pickup-and-delivery requests with
time windows must be fullfiled over a multi-day planning horizon. Unlike some
classic approaches, the correspondence between trucks and drivers is not fixed
and they can be exchanged in some locations and at any time. Drivers can also
travel for free as truck passengers or take external taxis for an additional
cost. The objective is to minimise the truck and taxi travel costs and the
penalties for late deliveries. Routes for trucks and drivers are represented
separately as directed paths in certain digraphs and then synchronised in time
and space. Three compact Integer Linear Programming formulations are proposed
and many families of valid inequalities are described. Extensive computational
experiments are conducted on randomly generated instances. The formulations
are experimentally compared and the effectiveness of the proposed valid
inequalities as cutting planes in a branch-and-cut algorithm is evaluated.

Keywords: Integer Programming; Valid inequalities; Long-haul transport;
Vehicle routing; Crew scheduling.

##  1 Introduction

Vehicle Routing Problems (VRPs) are among the most popular in transportation
and logistics. Applications go from supplying fuel to gas stations, collecting
fresh produce from farms, long-haul and last-mile transport, home health care,
and dial-a-ride services, to many others (see Toth and Vigo,, 2014). The
opportunity to minimise the high costs they entail in the supply chain and
their notorious hardness from the point of view of combinatorial optimization
explains the broad interest they have aroused. New requirements that
constantly emerge in real-world applications often limit or prohibit the
practical use of commercial software for general-purpose VRPs, motivating the
development of specific approaches.

In recent decades, there has been growing interest in considering hours-of-
service regulations for drivers during vehicle routing. They are generally
regulated by labour laws, collective labour agreements, and internal company
policies. In this work, we focus on long-haul transport, whose requirements
differ greatly from other transport sectors such as airways, railways, and
urban mass transit, since trips are not timetabled and can be interrupted for
drivers to take breaks. Considering that the truck travel times can extend to
several days, taking driver rest periods into account is essential to building
feasible routes. There is extensive literature in this regard, e.g. for United
States law see Rancourt et al., (2013); Goel and Vidal, (2014); Koç et al.,
(2018), for European Union law see Goel, (2009); Kok et al., (2010); Prescott-
Gagnon et al., (2010); Goel and Vidal, (2014), and for exact methods see Goel
and Irnich, (2017); Tilk and Goel, (2020). For crew scheduling in other
transport sectors see Deveci and Çetin Demirel, (2018); Heil et al., (2020);
Ibarra-Rojas et al., (2015). It is also worth mentioning the work of Goel et
al., (2021), who compare the advantages of team vs. single driving in European
road freight transport. For that social legislation, different laws apply to
team driving, e.g. the total daily driving time is doubled and a driver can
take a break while the other is driving. Generally, multiple objective
functions are optimised simultaneously though the most common are the number
of vehicles and drivers used, the total travel distance of the vehicles, and
the total working time of the drivers.

All the above works for long-haul transport assume a fixed correspondence
between trucks and drivers. That is, the same driver travels the entire route
of a truck, thus remaining unnecessarily idle during her/his rest periods.
There is room for improvement when drivers are allowed to change trucks over
time; e.g. when a driver needs to take a break, another can relieve her/him
and continue driving without delaying transport. This motivates some recent
works focused on Simultaneous Vehicle Routing and Crew Scheduling Problems
(SVRCSPs), where such a fixed correspondence is relaxed. For these variants,
any route plan where drivers do not change trucks remains feasible, but more
efficient use of resources can lead to new solutions to improve overall costs.
As a counterpart, the combinatorics grows since routes must also be built for
drivers, which must be synchronised in time and space with the truck routes.
Besides, sinchronisation constraints introduce interdependence, meaning that
changes to one route can affect the feasibility of the others (see Drexl,,
2012). In the brief literature on SVRCSPs, the most common resolution method
consists of a two-stage sequential decomposition, e.g. see Drexl et al.,
(2013); Mendes and Iori, (2020); Lucci et al., (2022). In the first stage,
restrictions on drivers are completely or partially ignored to build truck
routes that serve all customers. After that, drivers are assigned to segments
of the previous routes in compliance with labour laws.

Although decomposition is effective in solving large instances in reasonable
times, there is no guarantee of global optimality; e.g. crew scheduling could
be too costly or even impossible for a given input (truck routes).
Furthermore, decomposition indirectly introduces a hierarchical order between
the objectives, which makes it difficult to prioritise the objectives of the
second stage over those considered in the first. Exact methods that address
the problem in a single stage can overcome these drawbacks, but they are
usually more complex and more time-consuming. To our knowledge, very few works
in the literature deal with exact methods for SVRCSPs. Lam et al., (2015)
consider an application in humanitarian and military logistics, where crews
can change vehicles in different locations and also travel as passengers, and
the objective is to minimise a weighted sum that considers the number of
vehicles and crews used and the total vehicle and crew travel distances. The
authors propose a constraint programming formulation that is solved using a
LNS, a mixed integer program that is directly optimised by a general-purpose
ILP solver, and a two-stage method. From computational results on instances
with up to 96 requests, they conclude that a single-stage method produces
considerable beneﬁts over a sequential one, the constraint program scales
better than the mixed integer program, and decoupling crews from vehicles is
crucial to obtaining high-quality solutions. Domínguez-Martín et al., (2018)
consider an application inspired by local air traffic operations in the Canary
Islands, which involves two depots. Vehicle routes start and end at different
bases, but since crews always return to their starting base, drivers can
switch vehicles in some exchange locations and travel as passengers. The
authors present an ILP formulation to minimise the cost of driver routes,
introduce some valid inequalities, and develop separation procedures. They are
then incorporated into a branch-and-cut algorithm to solve instances with up
to 30 customer locations.

Both mathematical models proposed by Lam et al., (2015) and Domínguez-Martín
et al., (2018) limit vehicle exchanges to occur only when customers are being
served and also restrict drivers to single shifts. Some recent works seem to
reveal that greater savings can be achieved when drivers have more
opportunities to switch vehicles. Despite not involving a VRP, it is worth
mentioning the work of Ammann et al., (2023), who address a driver routing and
scheduling problem in long-distance bus networks, where drivers are allowed to
exchange buses at arbitrary intermediate stops en route. They define a MIP
formulation based on a time-expanded multi-digraph and propose a destructive-
bound-enhanced matheuristic. A computational study with real-world data
concludes that considering driver exchanges at intermediate stops in addition
to regular stops allows additional cost savings of 43%.

This work seeks to contribute to the study of exact single-stage methods for
solving SVRCSPs, where the planning horizon extends over multiple days, trucks
can have up to two drivers, and drivers can change trucks at any time in a set
of locations. A case study is followed that involves performing a set of
pickup-and-delivery transport requests with multiple time windows over long
distances, and the objective is to minimise a weighted sum that considers the
travel costs and penalties for late deliveries. We propose three digraphs-
based ILP formulations that model the movement of trucks and drivers over
time, evaluate alternatives for some of the constraints, study some families
of valid inequalities and their incorporation into branch-and-cut (B&C)
algorithms, and perform extensive computational experiments to compare the
algorithms and evaluate their limits. Some preliminaries were already
presented in the conference paper Lucci et al., (2022). A coffee distribution
company inspires the case study and was also addressed in Lucci et al., (2023)
through a two-stage sequential decomposition. In that previous work, hybrid
metaheuristics were proposed to solve both stages, including a LNS with SA
stopping criteria for the first stage and a GRASP$\times$ILS for the second.
Several generated instances were used to test the algorithms with up to 3000
requests distributed in 15 Argentine cities and a planning horizon of 1 to 4
weeks.

The rest of the paper is structured as follows. Section 2 presents the study
case. Section 3 describes the representation chosen for truck and driver
routes and how they can be synchronised. The ILP formulations are defined in
Section 4 and some families of valid inequalities in Section 5. Section 6
presents the computational experiments. Finally, the conclusion and future
work appear in Section 7.

##  2 Problem description

This work considers a case study that involves the planning of a homogeneous
fleet of trucks and drivers to transport pickup-and-delivery requests with
multiple time windows over long distances. The goal is to minimise the total
travel cost of the trucks and taxis and the penalties for customer
dissatisfaction caused by delivery delays. More details are provided below and
Fig. 1 summarises the notation/definition list.

$H\in\mathbb{N}_{0}$:

Number of days in the planning horizon.

$I\in\mathbb{N}$:

Number of time instants per day, under certain discretization.

[$n$]$\doteq\\{0,\ldots,n\\}$:

Set of integers from 0 to $n$.

$\mathcal{I}\doteq$[$IH$]:

Set of time instants in the planning horizon.

$L,V,D,R$:

Set of locations, trucks, drivers, and requests, respectively.

$l_{v},l_{d}\in L$:

Start locations of $v\in V$ and $d\in D$, respectively.

$L_{V}\doteq\\{l_{v}:v\in V\\}$:

Subset of truck start locations.

$L_{D}\doteq\\{l_{d}:d\in D\\}$:

Subset of driver start locations.

$l_{r}^{\textsc{p}},l_{r}^{\textsc{d}}\in L$:

Pickup and delivery locations of $r\in R$, respectively.

$j_{r}^{\textsc{p}},j_{r}^{\textsc{d}}\in$[$H-1$]:

Pickup and delivery start days of $r\in R$, respectively.

$a_{r}^{\textsc{p}},b_{r}^{\textsc{p}}\in$[$I-1$]:

Start and end times of the pickup time window of $r\in R$, respectively.

$a_{r}^{\textsc{d}},b_{r}^{\textsc{d}}\in$[$I-1$]:

Start and end times of the delivery time window of $r\in R$, respectively.

$s_{r}^{\textsc{p}},s_{r}^{\textsc{d}}\in\mathbb{N}_{0}$:

Pickup and delivery service times of $r\in R$, respectively.

$w_{r}\in\mathbb{N}_{0}$:

Penalty for each day of delay in the delivery of $r\in R$.

$\operatorname{day}:\mathcal{I}\to$[$H-1$]:

Given $i\in\mathcal{I}$, $\operatorname{day}(i)$ is the quotient of dividing
$i$ by $I$.

$\operatorname{time}:\mathcal{I}\to$[$I-1$]:

Given $i\in\mathcal{I}$, $\operatorname{time}(i)$ is the remainder of dividing
$i$ by $I$.

$\operatorname{len^{\textsc{v}}},\operatorname{len^{\textsc{t}}}:L\times
L\to\mathbb{N}$:

Truck and taxi travel time functions, respectively.

$\operatorname{cost^{\textsc{v}}},\operatorname{cost^{\textsc{t}}}:L\times
L\to\mathbb{N}$:

Truck and taxi travel cost functions, respectively.

$\mathcal{I}_{r}^{\textsc{p}},\mathcal{I}_{r}^{\textsc{d}}\subset\mathcal{I}$:

Subsets of instants where the pickup and delivery of $r\in R$ can begin, resp.

Figure 1: Notation/Definition list.

Each request involves picking up (or loading) cargo in a location,
transporting it to a different location, and delivering (or unloading) it
there. All requests must be delivered before the end of the planning horizon,
typically 1 week. The same truck that picks up a request must be in charge of
delivering it, i.e. cargo transshipment is not allowed, and each truck can
only transport a single request at a time. Loading/unloading cargo into/from a
truck takes 1 hour, known as service time, and the driver must remain in the
truck in the meantime. Each request has a pickup start day and a (hard) pickup
time window indicating the day and time the loading is allowed to begin, but
may end after hours. Analogously, it also has a delivery start day and a
(hard) delivery time window. These time constraints represent the day from
which the customer expects service and its hours of operation. Each time
window remains open from a start time to an end time (which can be on the next
day as long as the 24-hour limit is not exceeded) and repeats this behaviour
every day from its start day to the end of the planning horizon. Fig. 2a shows
an example of a time window that opens from day 0 during the time interval
[8:00, 20:00], where each box represents 1 hour. Instead, the one in Fig. 2b
does so during [20:00, 8:00], but when this happens (closing on the next day),
it must also open partially during [0:00, 8:00] on its start day (day 0) and
during [20:00, 0:00] on the last day of the planning horizon (day 1).
Deliveries have no deadlines, but each request has a penalty for each day of
delay from its delivery start day. For example, for delivery time windows as
Fig. 2, there is no penalty on day 0, a 1-day penalty on day 1, and so on.

Figure 2: Examples of time windows in a two-day planning horizon.

Each truck and driver has a start location, which represents their position at
the beginning of planning. They do not have to return there at the end of the
planning horizon, but their last location becomes the start of the next
planning. There is a set of locations where trucks and drivers can stop, and
they are not allowed to do so anywhere else, for example en route to some
location. This set contains all pickup and delivery locations of the requests
and the start locations of the trucks and drivers. The travel distance between
locations is known and the travel time is constant throughout the day.

Drivers are the only ones authorised to drive the trucks, but they can also
request an external taxi service to travel for an extra cost. Therefore, crew
scheduling becomes more flexible as drivers can take a taxi to a location
other than their current position without using a truck. There is an unlimited
number of taxis available in every location and at every time, but each can
carry one driver per ride. A driver is considered to be working while
travelling by vehicle (truck or taxi) and during service times and is supposed
to be resting the rest of the time. A rest begins when the driver gets out of
a vehicle (truck or taxi) and ends as soon as he/she gets in the same or a
different one. Drivers must rest at least 12 hours in each 24-hour interval
and must have at least 1 day off in each 7-day interval. As mentioned, drivers
can change trucks arbitrarily at any authorised location and at any time
(before, during, or after service). For example, a driver who needs to rest
can divert the truck to the nearest location, switch with another driver, and
immediately continue transport without further interruptions. Drivers are not
required to take a break after getting out of a vehicle, e.g. they can
immediately switch to another.

Optionally, each truck can carry an extra driver as a passenger, who is also
considered to be working. Thus, crews can be composed of 1 or 2 drivers, and
each one can be relieved independently. Although this increases combinatorics
even more, it may avoid taxis in certain situations. For example, a driver can
divert the truck to pick up a second driver and leave her/him in another
location before continuing the route, instead of taking a taxi. The choice
between travelling as a passenger or taking a taxi will depend on the routes
and the associated costs; e.g. it may be too expensive or impossible for a
driver to travel as a passenger if all the trucks are far away at the time and
need to travel long distances to her/his location.

Before proceeding to the mathematical models and more specific details, a toy
instance of the problem is exemplified and two possible route plans are
presented.

######  Example 1.

Consider an instance with a planning horizon of one day discretised in 3-hour
instants, i.e. $H=1$ and $I=8$. There two locations $L=\\{l_{1},l_{2}\\}$, two
vehicles $V=\\{v_{1},v_{2}\\}$, and two drivers $D=\\{d_{1},d_{2}\\}$. The
start location of $v_{1}$ and $d_{1}$ is $l_{1}$ and that of $v_{2}$ and
$d_{2}$ is $l_{2}$, i.e. $l_{v_{1}}=l_{d_{1}}=l_{1}$ and
$l_{v_{2}}=l_{d_{2}}=l_{2}$. The travel time between both locations is one
instant, i.e.
$\operatorname{len^{\textsc{v}}}(l_{1},l_{2})=\operatorname{len^{\textsc{v}}}(l_{2},l_{1})=\operatorname{len^{\textsc{t}}}(l_{1},l_{2})=\operatorname{len^{\textsc{t}}}(l_{2},l_{1})=1$.
There are two requests $R=\\{r_{1},r_{2}\\}$, whose attributes appear in Table
1.

Table 1: Request attributes in Example 1 Request | Pickup | Delivery
---|---|---
| Location | Start day | Time w. | Service t. | Location | Start day | Time w. | Service t.
$r\in R$ | $l^{\textsc{p}}_{r}$ | $j^{\textsc{p}}_{r}$ | $[a^{\textsc{p}}_{r},b^{\textsc{p}}_{r}]$ | $s^{\textsc{p}}_{r}$ | $l^{\textsc{d}}_{r}$ | $j^{\textsc{d}}_{r}$ | $[a^{\textsc{d}}_{r},b^{\textsc{d}}_{r}]$ | $s^{\textsc{d}}_{r}$
$r_{1}$ | $l_{1}$ | 0 | $[0,2]$ | 1 | $l_{2}$ | 0 | $[6,2]$ | 1
$r_{2}$ | $l_{2}$ | 0 | $[3,5]$ | 1 | $l_{1}$ | 0 | $[3,5]$ | 1

Fig. 3a depicts possible routes for the trucks and drivers. Recall that the
delivery of $r_{1}$ is allowed to begin at 6:00 since its time window opens
during [0,2] = [0:00,6:00] and [6,8] = [18:00,0:00] on day 0, as explained in
Fig. 2b. Fig. 3b depicts a better route for $d_{2}$, where he/she travels as a
passenger.

Figure 3: Routes of the Example 1. The tasks performed by truck $v_{1}$ are
orange and those by $v_{2}$ are blue.

##  3 Route representation

This section defines digraphs that represent truck and driver routes as
directed paths. The arcs specify the actions performed over time, which can be
resting, travelling, loading, or unloading. However, there are some categories
of directed paths that do not correspond to valid routes and should be
prohibited. In the end, the conditions that routes must meet to be
synchronised in space and time are identified.

###  3.1 Truck routes

#### Location–Time digraph

Let $G_{\textsc{lt}}\doteq(N,E,c)$ be a weighted multidigraph constructed as
follows. For each $l\in L$ and $i\in\mathcal{I}$, there is a node $n_{l,i}$ in
$N$, representing that a truck is in location $l$ at instant $i$. There are
also two distinguished nodes in $N$, the source $n_{s}$ and the sink
$\tilde{n}_{s}$. The multiset $E$ is the union of the following sets.

  * •

$E^{\textsc{rest}}$ has an arc $(n_{l,i},n_{l,i+1})$ for each $l\in L$ and
$i\in[|\mathcal{I}|-1]$, representing a rest of one time instant.

  * •

$E^{\textsc{trip}}$ has an arc $(n_{l_{1},i},n_{l_{2},i+\Delta})$ for each
pair $(l_{1},l_{2})\in L\times L$ of adjacent locations and
$i\in[|\mathcal{I}|-\Delta]$ with
$\Delta=\operatorname{len^{\textsc{v}}}(l_{1},l_{2})$, representing that a
truck travels.

  * •

For each $r\in R$, $E^{\textsc{p},r}$ has an arc
$(n_{l^{\textsc{p}}_{r},i},n_{l^{\textsc{p}}_{r},i+s^{\textsc{p}}_{r}})$ for
each $i\in\mathcal{I}^{\textsc{p}}_{r}$, representing that a truck loads $r$.

  * •

For each $r\in R$, $E^{\textsc{d},r}$ has an arc
$(n_{l^{\textsc{d}}_{r},i},n_{l^{\textsc{d}}_{r},i+s^{\textsc{d}}_{r}})$ for
each $i\in\mathcal{I}^{\textsc{d}}_{r}$, representing that a truck unloads
$r$.

  * •

$E^{\textsc{sour}}$ has an arc $(n_{s},n_{l,0})$ for each $l\in L_{V}$ and
$E^{\textsc{sink}}$ has an arc $(n_{l,|\mathcal{I}|},\tilde{n}_{s})$ for each
$l\in L$.

The weight of an arc $e=(n_{l_{1},i},n_{l_{2},i+\Delta})\in E^{\textsc{trip}}$
is the truck travel cost between $l_{1}$ and $l_{2}$, i.e.
$c(e)=\operatorname{cost^{\textsc{v}}}(l_{1},l_{2})$, the weight of
$e^{\prime}=(n_{l^{\textsc{d}}_{r},i},n_{l^{\textsc{d}}_{r},i+s^{\textsc{d}}_{r}})\in
E^{\textsc{d},r}$ is the product between the penalty of $r$ and the number of
days of delay, i.e.
$c(e^{\prime})=w_{r}(\operatorname{day}(i)-j^{\textsc{d}}_{r})$, and the
remaining are zero-weighted arcs. Fig. 4 shows the resulting construction for
the instance of Example 1. The position of the nodes on the vertical and
horizontal axis determines the location (first subscript) and the time instant
(second subscript), respectively. The arcs in $E^{\textsc{rest}}$ are parallel
to the horizontal axis and those in $E^{\textsc{trip}}$ are oblique and not
incident on the source/sink. For each $r\in R$, the ones in $E^{\textsc{p},r}$
and $E^{\textsc{d},r}$ are curved with the label “${\uparrow}r$” and
“${\downarrow}r$”, respectively, to distinguish them from other parallel arcs.
The orange $(n_{s},\tilde{n}_{s})$-directed path corresponds to the route of
the truck $v_{1}$ and the blue one to $v_{2}$.

$n_{s}$$\tilde{n}_{s}$$n_{l_{2},0}$$n_{l_{1},0}$$n_{l_{2},|\mathcal{I}|}$$n_{l_{1},|\mathcal{I}|}$${\uparrow}r_{1}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{2}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\downarrow}r_{2}$
Figure 4: Digraph $G_{\textsc{lt}}$ for Example 1.

When there are ambiguities, the node and arc sets of $G_{\textsc{lt}}$ will
also have the subscript “lt”, e.g. $N_{\textsc{lt}}$, $E_{\textsc{lt}}$,
$E_{\textsc{lt}}^{\textsc{rest}}$, etc., as well as the weight vector
$c_{\textsc{lt}}$. To avoid referring to particular requests, the multisets
$E^{\textsc{p}}$ and $E^{\textsc{d}}$ will denote ${\textstyle\bigcup_{r\in
R}E^{\textsc{p},r}}$ and ${\textstyle\bigcup_{r\in R}E^{\textsc{d},r}}$,
respectively. By construction, $G_{\textsc{lt}}$ is an acyclic digraph where
every truck route corresponds to a $(n_{s},\tilde{n}_{s})$-directed path.
However, the converse is not necessarily true. Let $P$ be a
$(n_{s},\tilde{n}_{s})$-directed path, four categories of forbidden paths can
be identified based on the requirements of the case study.

  * •

Repeated services. Since each request must be picked up and delivered exactly
once, then for all $r\in R$, $P$ cannot have more than one arc of
$E^{\textsc{p},r}$ and $E^{\textsc{d},r}$.

  * •

Unpaired services. Since each request must be fully served by the same truck,
then for all $r\in R$, the total number of arcs of $E^{\textsc{p},r}$ and
$E^{\textsc{d},r}$ in $P$ cannot differ. Two of these forbidden paths for the
instance in Example 1 are depicted in Fig. 5(a).

  * •

Disordered services. Since each request must be picked up before being
delivered, then for all $r\in R$ and $P^{\prime}$ subdirected path of $P$ from
$n_{s}$ such that its last arc belongs to $E^{\textsc{d},r}$, the total number
of arcs of $E^{\textsc{p},r}$ in $P^{\prime}$ cannot be zero. Fig. 5(b) shows
an example of such a forbidden path that also has no unpaired services.

  * •

Excess capacity. Since each truck can service only one request at a time, then
for all $P^{\prime}$ subdirected path of $P$ from $n_{s}$ such that its last
arc belongs to $E^{\textsc{p}}$, the difference between the total number of
arcs of $E^{\textsc{p}}$ and $E^{\textsc{d}}$ in $P^{\prime}$ cannot be
greater than one. Fig. 5(c) shows an example of such a forbidden path that
also has no unpaired or disordered services.

(a)$n_{s}$$\tilde{n}_{s}$$n_{l_{2},0}$$n_{l_{1},0}$$n_{l_{2},|\mathcal{I}|}$$n_{l_{1},|\mathcal{I}|}$${\uparrow}r_{1}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{2}$${\uparrow}r_{1}$${\downarrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{1}$

(b)$n_{s}$$\tilde{n}_{s}$$n_{l_{2},0}$$n_{l_{1},0}$$n_{l_{2},|\mathcal{I}|}$$n_{l_{1},|\mathcal{I}|}$${\uparrow}r_{1}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{2}$${\uparrow}r_{1}$${\downarrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{1}$

(c)$n_{s}$$\tilde{n}_{s}$$n_{l_{2},0}$$n_{l_{1},0}$$n_{l_{2},|\mathcal{I}|}$$n_{l_{1},|\mathcal{I}|}$${\uparrow}r_{1}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{2}$${\uparrow}r_{1}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{1}$

Figure 5: Examples of forbidden paths in $G_{\textsc{lt}}$.

In the remainder of the subsection, digraphs with greater structure are
proposed that prevent the existence of certain types of forbidden paths.

#### Location–Time–Cargo digraph

Let $G_{\textsc{ltc}}\doteq(N,E,c)$ be a weighted multidigraph constructed as
follows. For each $l\in L$, $i\in\mathcal{I}$, and $q\in\\{\circ,\bullet\\}$,
there is a node $n_{l,i}^{q}$ in $N$, representing that the truck is empty in
location $l$ at instant $i$ if $q=\circ$ and is carrying cargo if $q=\bullet$.
There are also two distinguished nodes in $N$, the source $n_{s}$ and the sink
$\tilde{n}_{s}$. The multiset $E$ follows the construction of
$E_{\textsc{lt}}$, taking into account that the pickup and delivery arcs are
the only ones that modify cargo and that trucks must begin and end their route
empty. Formally, $E$ is the union of the following sets.

  * •

$E^{\textsc{rest}}$ has an arc $(n_{l,i}^{q},n_{l,i+1}^{q})$ for each $l\in
L$, $i\in[|\mathcal{I}|-1]$, and $q\in\\{\circ,\bullet\\}$.

  * •

$E^{\textsc{trip}}$ has an arc $(n_{l_{1},i}^{q},n_{l_{2},i+\Delta}^{q})$ for
each pair $(l_{1},l_{2})\in L\times L$ of adjacent locations,
$i\in[|\mathcal{I}|-\Delta]$ with
$\Delta=\operatorname{len^{\textsc{v}}}(l_{1},l_{2})$, and
$q\in\\{\circ,\bullet\\}$.

  * •

For each $r\in R$, $E^{\textsc{p},r}$ has an arc
$(n_{l^{\textsc{p}}_{r},i}^{\circ},n_{l^{\textsc{p}}_{r},i+s^{\textsc{p}}_{r}}^{\bullet})$
for each $i\in\mathcal{I}^{\textsc{p}}_{r}$.

  * •

For each $r\in R$, $E^{\textsc{d},r}$ has an arc
$(n_{l^{\textsc{d}}_{r},i}^{\bullet},n_{l^{\textsc{d}}_{r},i+s^{\textsc{d}}_{r}}^{\circ})$
for each $i\in\mathcal{I}^{\textsc{d}}_{r}$.

  * •

$E^{\textsc{sour}}$ has an arc $(n_{s},n_{l,0}^{\circ})$ for each $l\in L_{V}$
and $E^{\textsc{sink}}$ has an arc
$(n_{l,|\mathcal{I}|}^{\circ},\tilde{n}_{s})$ for each $l\in L$.

The definition of the weight vector $c$ is omitted since its adaptation is
straightforward. Fig. 6 shows the resulting construction for the instance of
Example 1. The two lower rows of nodes have the superscript “$\circ$” and the
two upper ones have “$\bullet$”. The pickup and delivery arcs are now oblique
but, unlike travel arcs, they are incident on nodes with different
superscripts. For the sake of simplicity, the given definition of
$G_{\textsc{ltc}}$ omitted some technicalities. There may be nodes with the
superscript “$\bullet$” through which it is impossible for a
$(n_{s},\tilde{n}_{s})$-directed path to pass; e.g. nodes before the end time
of the earliest pickup arc or after the start time of the latest delivery arc.
Such nodes and the arcs incidents on them, which are drawn with a light dashed
line, can be erased from the digraph.

$n_{s}$$\tilde{n}_{s}$$n^{\circ}_{l_{2},0}$$n^{\circ}_{l_{1},0}$$n^{\bullet}_{l_{2},0}$$n^{\bullet}_{l_{1},0}$$n^{\circ}_{l_{2},|\mathcal{I}|}$$n^{\circ}_{l_{1},|\mathcal{I}|}$$n^{\bullet}_{l_{2},|\mathcal{I}|}$$n^{\bullet}_{l_{1},|\mathcal{I}|}$${\uparrow}r_{1}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{2}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\downarrow}r_{2}$
Figure 6: Digraph $G_{\textsc{ltc}}$ for Example 1.

The greater structure of $G_{\textsc{ltc}}$ guarantees that every
$(n_{s},\tilde{n}_{s})$-directed path already alternates arcs of
$E^{\textsc{p}}$ and $E^{\textsc{d}}$. Consequently, there cannot be excess
capacity but services may still be unpaired or disordered. In the worst case,
$G_{\textsc{ltc}}$ has $|\mathcal{I}|.|L|$ more nodes and
$|E^{\textsc{rest}}_{\textsc{lt}}|+|E^{\textsc{trip}}_{\textsc{lt}}|$ more
edges than $G_{\textsc{lt}}$. The following digraph presents an even more
structured approach to modelling truck routes.

#### Location–Time–Request digraph

The weighted (simple) digraph $G_{\textsc{ltr}}\doteq(N,E,c)$ is constructed
as follows. For each $l\in L$, $i\in\mathcal{I}$, and $r\in R^{\,\circ}\doteq
R\cup\\{\circ\\}$, there is a node $n_{l,i}^{r}$ in $N$, representing that the
truck is empty in location $l$ at instant $i$ if $r=\circ$ and is carrying
request $r$ otherwise. There are also two distinguished nodes in $N$, the
source $n_{s}$ and the sink $\tilde{n}_{s}$. The arc set $E$ follows the
construction of $E_{\textsc{ltc}}$ and is the union of the following sets.

  * •

$E^{\textsc{rest}}$ has an arc $(n_{l,i}^{r},n_{l,i+1}^{r})$ for each $l\in
L$, $i\in[|\mathcal{I}|-1]$, and $r\in R^{\,\circ}$.

  * •

$E^{\textsc{trip}}$ has an arc $(n_{l_{1},i}^{r},n_{l_{2},i+\Delta}^{r})$ for
each pair $(l_{1},l_{2})\in L\times L$ of adjacent locations,
$i\in[|\mathcal{I}|-\Delta]$ with
$\Delta=\operatorname{len^{\textsc{v}}}(l_{1},l_{2})$, and $r\in R^{\,\circ}$.

  * •

For each $r\in R$, $E^{\textsc{p},r}$ has an arc
$(n_{l^{\textsc{p}}_{r},i}^{\circ},n_{l^{\textsc{p}}_{r},i+s^{\textsc{p}}_{r}}^{r})$
for each $i\in\mathcal{I}^{\textsc{p}}_{r}$.

  * •

For each $r\in R$, $E^{\textsc{d},r}$ has an arc
$(n_{l^{\textsc{d}}_{r},i}^{r},n_{l^{\textsc{d}}_{r},i+s^{\textsc{d}}_{r}}^{\circ})$
for each $i\in\mathcal{I}^{\textsc{d}}_{r}$.

  * •

$E^{\textsc{sour}}$ has an arc $(n_{s},n_{l,0}^{\circ})$ for each $l\in L_{V}$
and $E^{\textsc{sink}}$ has an arc
$(n_{l,|\mathcal{I}|}^{\circ},\tilde{n}_{s})$ for each $l\in L$.

For $r\in R^{\,\circ}$, we will refer to $E^{\textsc{rest},r}$ (resp.
$E^{\textsc{trip},r}$) as the subset of rest arcs (resp. trip arcs) connecting
nodes with the superscript $r$. Once again, the adaptation of the weight
vector $c$ is omitted. Fig. 7 shows the resulting construction for the
instance of Example 1. Now, the two rows of nodes in the middle have the
superscript “$\circ$”, the upper ones have “$r_{1}$”, and the lower ones have
“$r_{2}$”. Again, the definition given for $G_{\textsc{ltr}}$ simplifies some
technicalities in the construction. Particularly, for all $r\in R$, the nodes
with the superscript “$r$” before the end time of the earliest arc in
$E^{\textsc{p},r}$ and the ones after the start time of the latest arc in
$E^{\textsc{d},r}$ can be erased along with the arcs incidence on them.

$n_{s}$$\tilde{n}_{s}$$n^{r_{2}}_{l_{2},0}$$n^{r_{2}}_{l_{1},0}$$n^{\circ}_{l_{2},0}$$n^{\circ}_{l_{1},0}$$n^{r_{1}}_{l_{2},0}$$n^{r_{1}}_{l_{1},0}$$n^{r_{2}}_{l_{2},|\mathcal{I}|}$$n^{r_{2}}_{l_{1},|\mathcal{I}|}$$n^{\circ}_{l_{2},|\mathcal{I}|}$$n^{\circ}_{l_{1},|\mathcal{I}|}$$n^{r_{1}}_{l_{2},|\mathcal{I}|}$$n^{r_{1}}_{l_{1},|\mathcal{I}|}$${\uparrow}r_{1}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{2}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\downarrow}r_{2}$
Figure 7: Digraph $G_{\textsc{ltr}}$ for Example 1.

The major benefit of $G_{\textsc{ltr}}$ is that its construction already
guarantees that every $(n_{s},\tilde{n}_{s})$-directed path alternates arcs of
$E^{\textsc{p}}$ and $E^{\textsc{d}}$ associated with the same request. Thus,
the services are already paired and well-ordered and there can be no excess
capacity. The main disadvantage is associated with its size, which has
$|R|.|\mathcal{I}|.|L|$ more nodes and
$|R|.(|E^{\textsc{rest}}_{\textsc{lt}}|+|E^{\textsc{trip}}_{\textsc{lt}}|)$
more arcs that $G_{\textsc{lt}}$ in the worst case. The following subsection
presents a digraph that can be utilised to represent the driver routes.

###  3.2 Driver routes

#### Location–Time digraph with taxi arcs

Let $G_{\textsc{ltx}}\doteq(N,E,c)$ be the weighted multidigraph that has the
same node set as $G_{\textsc{lt}}$, i.e. $N=N_{\textsc{lt}}$, and whose arc
set is the union of the following sets.

  * •

$E^{\textsc{rest}}$, $E^{\textsc{trip}}$, $E^{\textsc{p}}$, $E^{\textsc{r}}$,
and $E^{\textsc{sink}}$, as defined in $G_{\textsc{lt}}$.

  * •

$E^{\textsc{taxi}}$ has an arc $(n_{l_{1},i},n_{l_{2},i+\Delta})$ for each
pair $(l_{1},l_{2})\in L\times L$ of adjacent locations and
$i\in[|\mathcal{I}|-\Delta]$ with
$\Delta=\operatorname{len^{\textsc{t}}}(l_{1},l_{2})$, representing that a
driver takes a taxi.

  * •

$E^{\textsc{sour}}$ has an arc $(n_{s},n_{l,0}^{\circ})$ for each $l\in
L_{D}$.

The weight of an arc $e=(n_{l_{1},i},n_{l_{2},i+\Delta})\in E^{\textsc{taxi}}$
is the taxi travel cost between $l_{1}$ and $l_{2}$, i.e.
$c(e)=\operatorname{cost^{\textsc{t}}}(l_{1},l_{2})$, and the remaining are
zero-weighted arcs. Fig. 8 shows the resulting construction for the instance
of Example 1 and the driver routes of Fig. 3(a). Since this example considers
travel times to be identical for trucks and taxis, the arcs in
$E^{\textsc{taxi}}$ are parallel to those in $E^{\textsc{trip}}$. To
distinguish them, taxi arcs are drawn curved with the label “t”. The green
$(n_{s},\tilde{n}_{s})$-directed path corresponds to the route of $d_{1}$ and
the pink one to $d_{2}$.

$n_{s}$$\tilde{n}_{s}$$n_{l_{2},0}$$n_{l_{1},0}$ttttttttttttttt${\uparrow}r_{1}$${\uparrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\downarrow}r_{1}$${\uparrow}r_{2}$${\uparrow}r_{2}$${\downarrow}r_{2}$${\downarrow}r_{2}$${\uparrow}r_{1}$${\uparrow}r_{2}$${\downarrow}r_{1}$t${\downarrow}r_{2}$
Figure 8: Digraph $G_{\textsc{ltx}}$ for Example 1.

The construction of $G_{\textsc{ltx}}$ guarantees that every driver route
corresponds to a $(n_{s},\tilde{n}_{s})$-directed path, but it is also
necessary to forbid some paths for the converse to hold. Let $P$ be a
$(n_{s},\tilde{n}_{s})$-directed path, there are two categories of forbidden
paths according to the labour laws of the case study.

  * •

Daily rest. Since drivers must rest at least 12 hours in each 24-hour
interval, then for all $i\in[I(H-1)]$ and $P^{\prime}$ subdirected path of $P$
from nodes $n_{l_{1},i}$ to $n_{l_{2},i+I}$ with $l_{1},l_{2}\in L$, the total
number of arcs of $E^{\textsc{rest}}$ in $P^{\prime}$ cannot be lower than
$\frac{I}{2}$.

  * •

Weekly rest. Since drivers must have at least 1 day off in each 7-day
interval, then for all $j\in[H-7]$ and $P^{\prime}$ subdirected path of $P$
from nodes $n_{l_{1},Ij}$ to $n_{l_{2},I(j+7)}$ with $l_{1},l_{2}\in L$, there
must be a subdirected path of $P^{\prime}$ with all its arcs in
$E^{\textsc{rest}}$ from nodes $n_{l,Ij^{\prime}}$ to $n_{l,I(j^{\prime}+1)}$
for some $l\in L$ and $j^{\prime}\in\\{j,\ldots,j+6\\}$.

Unlike the approach followed for trucks, no further representations are
proposed for driver routes, since enriching the nodes with cargo information
has no clear advantage for them. In what follows, we will discuss how to
synchronise the routes of the trucks and drivers.

###  3.3 Synchronisation of routes

Synchronization is of utmost importance to ensure that trucks and drivers
match in time and space. The selected representation for the routes has a
significant benefit as it enables synchronization to be expressed in a
simplified manner. The most direct case is when truck routes are represented
with $G_{\textsc{lt}}$ and driver routes with $G_{\textsc{ltx}}$, since the
arc set of the former is contained in the latter. For each truck route that
passes through a trip, pickup, or delivery arc (which are the ones that
require synchronization), there must be some driver route that passes through
that same arc, and vice-versa. However, since drivers can travel as
passengers, at most two of them can share the same arc of a truck route.
Formally, for each $e\in E^{\textsc{trip}}_{\textsc{ltx}}\cup
E^{\textsc{p}}_{\textsc{ltx}}\cup E^{\textsc{d}}_{\textsc{ltx}}$, if $V_{e}$
and $D_{e}$ are the numbers of trucks and drivers whose routes pass through
$e$, respectively, then it must hold $V_{e}\leq D_{e}\leq 2V_{e}.$

For the other representations, the arc correspondence between the digraphs is
less direct and no longer one-to-one. For example, a driver can travel in
either an empty or a full truck, but $G_{\textsc{ltc}}$ and $G_{\textsc{ltr}}$
identify such trip with different arcs (varying the superscript of the nodes
they connect). To overcome this, functions $\operatorname{arc_{\textsc{ltc}}}$
and $\operatorname{arc_{\textsc{ltr}}}$ are defined, that return the subset of
arcs of $E_{\textsc{ltc}}$ and $E_{\textsc{ltr}}$ that match a given $e\in
E_{\textsc{ltx}}$, respectively. For the sake of brevity, we present only the
relevant cases of $\operatorname{arc_{\textsc{ltc}}}$, even though the
definition can be extended to cover its entire domain, as well as extended to
$\operatorname{arc_{\textsc{ltr}}}$.

$\displaystyle\operatorname{arc_{\textsc{ltc}}}(e)\doteq$$\displaystyle\begin{cases}\\{e^{\prime}\in
E^{\textsc{trip}}_{\textsc{ltc}}:e^{\prime}=(n^{q}_{l_{1},i_{1}},n^{q}_{l_{2},i_{2}})\wedge
q\in\\{\circ,\bullet\\}\\}&\text{if }e=(n_{l_{1},i_{1}},n_{l_{2},i_{2}})\in
E^{\textsc{trip}}_{\textsc{ltx}}.\\\ \\{e^{\prime}\in
E^{\textsc{p},r}_{\textsc{ltc}}:e^{\prime}=(n^{\circ}_{l,i_{1}},n^{\bullet}_{l,i_{2}})\\}&\text{if
}e=(n_{l,i_{1}},n_{l,i_{2}})\in E^{\textsc{p},r}_{\textsc{ltx}}.\\\
\\{e^{\prime}\in
E^{\textsc{d},r}_{\textsc{ltc}}:e^{\prime}=(n^{\bullet}_{l,i_{1}},n^{\circ}_{l,i_{2}})\\}&\text{if
}e=(n_{l,i_{1}},n_{l,i_{2}})\in E^{\textsc{d},r}_{\textsc{ltx}}.\\\
\end{cases}$

Then, the previous inequality can be stated as
$\sum_{e^{\prime}\in\operatorname{arc_{\textsc{ltc}}}(e)}V_{e^{\prime}}\leq
D_{e}\leq
2\sum_{e^{\prime}\in\operatorname{arc_{\textsc{ltc}}}(e)}V_{e^{\prime}}.$

##  4 Integer linear programs

Before presenting the ILP formulations, some additional notations must be
introduced. Given a node $n\in N$, $\Gamma^{-}(n)\subset E$ (resp.
$\Gamma^{+}(n)$) is the subset of incoming (resp. outgoing) arcs of $n$. Given
$n_{l,i}\in N$ and $e\in\Gamma^{+}(n_{l,i})$ (resp.
$e\in\Gamma^{-}(n_{l,i})$), then $\operatorname{L}(e)\doteq l$ and
$\operatorname{\mathcal{I}}(e)\doteq i$ (resp.
$\operatorname{L^{\prime}}(e)\doteq l$ and
$\operatorname{\mathcal{I}^{\prime}}(e)\doteq i$). Given $v\in V$ and $d\in
D$, $E^{v}\subset E$ and $E^{d}\subset E$ are the subsets of possible arcs for
$v$ and $d$ according to their start location, resp., i.e. $E^{v}\doteq
E\setminus\\{e\in E^{\textsc{sour}}:\operatorname{L^{\prime}}(e)\neq l_{v}\\}$
and $E^{d}\doteq E\setminus\\{e\in
E^{\textsc{sour}}:\operatorname{L^{\prime}}(e)\neq l_{d}\\}$.

#### LT formulation

The first ILP formulation is based on $G_{\textsc{lt}}$ and $G_{\textsc{ltx}}$
to represent truck and driver routes, respectively. For each $v\in V$ and
$e\in E^{v}_{\textsc{lt}}$, there is a binary variable $X_{ve}$ such that
$X_{ve}=1$ if and only if $v$ passes through $e$. For each $d\in D$ and $e\in
E^{d}_{\textsc{ltx}}$, there is a binary variable $Y_{de}$ such that
$Y_{de}=1$ if and only if $d$ passes through $e$. For each $d\in D$ and
$j\in[H-1]$, there is a binary variable $W_{dj}$ such that $W_{dj}=1$ if $d$
has day $j$ off.

(LT) $\displaystyle\min\sum_{v\in V}\sum_{e\in
E^{v}_{\textsc{lt}}}c_{\textsc{lt}}(e)X_{ve}+\sum_{d\in D}\sum_{e\in
E^{d}_{\textsc{ltx}}}c_{\textsc{ltx}}(e)Y_{de}$ (1) s.t.
$\displaystyle\sum_{\mathclap{e\in\Gamma^{-}(n)\cap
E^{v}_{\textsc{lt}}}}\hskip 7.0ptX_{ve}=\hskip
7.0pt\sum_{\mathclap{e\in\Gamma^{+}(n)\cap E^{v}_{\textsc{lt}}}}\hskip
5.0ptX_{ve}$ $\displaystyle\forall v\in V,\ n\in
N_{\textsc{lt}}\setminus\\{n_{s},\tilde{n}_{s}\\},$ (2)
$\displaystyle\sum_{\mathclap{e\in\Gamma^{-}(n)\cap
E^{d}_{\textsc{ltx}}}}\hskip 7.0ptY_{de}=\hskip
7.0pt\sum_{\mathclap{e\in\Gamma^{+}(n)\cap E^{d}_{\textsc{ltx}}}}\hskip
5.0ptY_{de}$ $\displaystyle\forall d\in D,\ n\in
N_{\textsc{ltx}}\setminus\\{n_{s},\tilde{n}_{s}\\},$ (3)
$\displaystyle\sum_{v\in V}\sum_{e\in E^{\textsc{p},r}_{\textsc{lt}}}X_{ve}=1$
$\displaystyle\forall r\in R,$ (4) $\displaystyle\sum_{\mathclap{e\in
E^{\textsc{p},r}_{\textsc{lt}}}}X_{ve}=\sum_{e\in
E^{\textsc{d},r}_{\textsc{lt}}}X_{ve}$ $\displaystyle\forall v\in V,\ r\in R,$
(5) $\displaystyle\sum_{v\in V}\hskip
8.0pt\sum_{\mathclap{\begin{subarray}{c}e^{\prime}\in
E^{\textsc{p},r}_{\textsc{lt}}:\\\
\operatorname{\mathcal{I}^{\prime}}(e^{\prime})\leq\operatorname{\mathcal{I}}(e)\end{subarray}}}\hskip
3.0ptX_{ve^{\prime}}\geq\sum_{v\in V}X_{ve}$ $\displaystyle\forall r\in R,\
e\in E^{\textsc{d},r}_{\textsc{lt}},$ (6) $\displaystyle\sum_{r\in R}\hskip
4.0pt\Bigl{(}\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{p},r}_{\textsc{lt}}:\\\ \operatorname{\mathcal{I}^{\prime}}(e)\leq
i\end{subarray}}}X_{ve}-\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{d},r}_{\textsc{lt}}:\\\ \operatorname{\mathcal{I}^{\prime}}(e)\leq
i\end{subarray}}}X_{ve}\Bigr{)}\leq 1$ $\displaystyle\forall v\in V,\
i\in{\textstyle\bigcup_{r\in R}\mathcal{I}^{\textsc{p}}_{r}},$ (7)
$\displaystyle\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{rest}}_{\textsc{ltx}}:\\\
i\leq\operatorname{\mathcal{I}}(e)<i+I\end{subarray}}}\hskip 5.0ptY_{de}\geq
I/2$ $\displaystyle\forall d\in D,\ i\in[I(H-1)],$ (8)
$\displaystyle\sum_{\mathclap{j\leq j^{\prime}\leq j+6}}\hskip
5.0ptW_{dj^{\prime}}\geq 1$ $\displaystyle\forall d\in D,\ j\in[H-7],$ (9)
$\displaystyle\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{desc}}_{\textsc{ltx}}:\\\
Ij\leq\operatorname{\mathcal{I}}(e)<I(j+1)\end{subarray}}}\hskip
5.0ptY_{de}\geq IW_{dj}$ $\displaystyle\forall d\in D,\ j\in[H-1],$ (10)
$\displaystyle\sum_{v\in V}X_{ve}\leq\sum_{d\in D}Y_{de}\leq 2\sum_{v\in
V}X_{ve}$ $\displaystyle\forall e\in E^{\textsc{trip}}_{\textsc{ltx}}\cup
E^{\textsc{p}}_{\textsc{ltx}}\cup E^{\textsc{d}}_{\textsc{ltx}},$ (11)
$\displaystyle X_{ve}\in\\{0,1\\}$ $\displaystyle\forall v\in V,\ e\in
E^{v}_{\textsc{lt}},$ (12) $\displaystyle Y_{de}\in\\{0,1\\}$
$\displaystyle\forall d\in D,\ e\in E^{c}_{\textsc{ltx}},$ (13) $\displaystyle
W_{dj}\in\\{0,1\\}$ $\displaystyle\forall d\in D,\ j\in[H-1].$ (14)

The objective function (1) minimises total travel costs and penalties for late
deliveries. Constraints (2)–(3) are the flow conservation equations for trucks
and drivers, respectively, and they guarantee that each route is empty or a
$(n_{s},\tilde{n}_{s})$-directed path. Then, (4)–(7) and (8)–(10) prohibit
directed paths that do not correspond to truck and driver routes according to
the categories of Section 3.1 and 3.2, respectively. Observe that (4) also
prevent a request from being executed by more than one truck and adding
similar constraints for the deliveries is not necessary since they are
implied. When considering individual truck routes, (4) and (6) disaggregated
by truck are also implied. While (9) force each driver to have at least one
day off, (10) guarantee that their route contains only rest arcs on such day.
Constraints (11) synchronise the routes. Finally, (12)–(14) are the
integrality constraints.

#### LTC formulation

In this ILP formulation, $G_{\textsc{ltc}}$ is used to represent the truck
routes. In addition to the above variables $Y$ and $W$, for each $v\in V$ and
$e\in E^{v}_{\textsc{ltc}}$, there is a binary variable $X_{ve}$ such that
$X_{ve}=1$ if and only if $v$ passes through $e$.

(LTC) $\displaystyle\min\sum_{v\in V}\sum_{e\in
E^{v}_{\textsc{ltc}}}c_{\textsc{ltc}}(e)X_{ve}+\sum_{d\in D}\sum_{e\in
E^{d}_{\textsc{ltx}}}c_{\textsc{ltx}}(e)Y_{de}$ (15) s.t.
$\displaystyle(\ref{form.LT.LTX.Driver.Flow}),(\ref{form.LT.LTX.Driver.Daily}),(\ref{form.LT.LTX.Driver.Weekly.I}),(\ref{form.LT.LTX.Driver.Weekly.II}),(\ref{form.LT.LTX.Range.II}),(\ref{form.LT.LTX.Range.III}),$
$\displaystyle\sum_{\mathclap{e\in\Gamma^{-}(n)\cap
E^{v}_{\textsc{ltc}}}}\hskip 7.0ptX_{ve}=\hskip
7.0pt\sum_{\mathclap{e\in\Gamma^{+}(n)\cap E^{v}_{\textsc{ltc}}}}\hskip
5.0ptX_{ve}$ $\displaystyle\forall v\in V,\ n\in
N_{\textsc{ltc}}\setminus\\{n_{s},\tilde{n}_{s}\\},$ (16)
$\displaystyle\sum_{v\in V}\sum_{e\in
E^{\textsc{p},r}_{\textsc{ltc}}}X_{ve}=1$ $\displaystyle\forall r\in R,$ (17)
$\displaystyle\sum_{\mathclap{e\in
E^{\textsc{p},r}_{\textsc{ltc}}}}X_{ve}=\sum_{e\in
E^{\textsc{d},r}_{\textsc{ltc}}}X_{ve}$ $\displaystyle\forall v\in V,\ r\in
R,$ (18) $\displaystyle\sum_{v\in V}\hskip
8.0pt\sum_{\mathclap{\begin{subarray}{c}e^{\prime}\in
E^{\textsc{p},r}_{\textsc{ltc}}:\\\
\operatorname{\mathcal{I}^{\prime}}(e^{\prime})\leq\operatorname{\mathcal{I}}(e)\end{subarray}}}\hskip
3.0ptX_{ve^{\prime}}\geq\sum_{v\in V}X_{ve}$ $\displaystyle\forall r\in R,\
e\in E^{\textsc{d},r}_{\textsc{ltc}},$ (19) $\displaystyle\sum_{v\in V}\hskip
3.0pt\sum_{\mathclap{\begin{subarray}{c}e^{\prime}\in\\\
\operatorname{arc_{\textsc{ltc}}}(e)\end{subarray}}}X_{ve^{\prime}}\leq\sum_{d\in
D}Y_{de}\leq 2\sum_{v\in V}\hskip
3.0pt\sum_{\mathclap{\begin{subarray}{c}e^{\prime}\in\\\
\operatorname{arc_{\textsc{ltc}}}(e)\end{subarray}}}X_{ve^{\prime}}$
$\displaystyle\forall e\in E^{\textsc{trip}}_{\textsc{ltx}}\cup
E^{\textsc{p}}_{\textsc{ltx}}\cup E^{\textsc{d}}_{\textsc{ltx}},$ (20)
$\displaystyle X_{ve}\in\\{0,1\\}$ $\displaystyle\forall v\in V,\ e\in
E^{v}_{\textsc{ltc}}.$ (21)

#### LTR formulation

Finally, the ILP formulation is presented that considers $G_{\textsc{ltr}}$ to
represent the truck routes. In addition to the above variables $Y$ and $W$,
for each $e\in E_{\textsc{ltr}}$, there is an integer variable $X_{e}$ whose
value is the number of trucks that passes through $e$. The range of $X_{e}$ is
upper bounded by the capacity of $e$. Particularly, the number of trucks in
each start location defines the capacity for each outgoing arc of the source
node, the capacity for the arcs where the truck is loaded is one to avoid
repeated services, otherwise the capacity is limited by the total number of
trucks.

(LTR) $\displaystyle\min\sum_{e\in
E_{\textsc{ltr}}}c_{\textsc{ltr}}(e)X_{e}+\sum_{d\in D}\sum_{e\in
E^{d}_{\textsc{ltx}}}c_{\textsc{ltx}}(e)Y_{de}$ (22) s.t.
$\displaystyle(\ref{form.LT.LTX.Driver.Flow}),(\ref{form.LT.LTX.Driver.Daily}),(\ref{form.LT.LTX.Driver.Weekly.I}),(\ref{form.LT.LTX.Driver.Weekly.II}),(\ref{form.LT.LTX.Range.II}),(\ref{form.LT.LTX.Range.III}),$
$\displaystyle\sum_{\mathclap{e\in\Gamma^{-}(n)}}X_{e}=\sum_{\mathclap{e\in\Gamma^{+}(n)}}X_{e}$
$\displaystyle\forall n\in
N_{\textsc{ltr}}\setminus\\{n_{s},\tilde{n}_{s}\\},$ (23)
$\displaystyle\sum_{\mathclap{e\in E^{\textsc{p},r}_{\textsc{ltr}}}}X_{e}=1$
$\displaystyle\forall r\in R,$ (24)
$\displaystyle\sum_{\mathclap{\begin{subarray}{c}e^{\prime}\in\operatorname{arc_{\textsc{ltr}}}(e)\end{subarray}}}X_{e^{\prime}}\leq\sum_{d\in
D}Y_{de}\leq
2\sum_{\mathclap{\begin{subarray}{c}e^{\prime}\in\operatorname{arc_{\textsc{ltr}}}(e)\end{subarray}}}X_{e^{\prime}}$
$\displaystyle\forall e\in E^{\textsc{trip}}_{\textsc{ltx}}\cup
E^{\textsc{p}}_{\textsc{ltx}}\cup E^{\textsc{d}}_{\textsc{ltx}},$ (25)
$\displaystyle X_{e}\in\\{0,\ldots,\operatorname{cap}(e)\\}$
$\displaystyle\forall e\in E_{\textsc{ltr}}.$ (26)

Where the function $\operatorname{cap}$ is defined as
$\operatorname{cap}(e)\doteq|\\{v\in V:l_{v}=l\\}|$ if $e\in
E^{\textsc{sour}}_{\textsc{ltr}}\text{ and }\operatorname{L^{\prime}}(e)=l$,
$\operatorname{cap}(e)\doteq 1$ if, for some $r\in R,e\in
E^{\textsc{p},r}_{\textsc{ltr}}\cup E^{\textsc{d},r}_{\textsc{ltr}}\cup
E^{\textsc{rest},r}_{\textsc{ltr}}\cup E^{\textsc{trip},r}_{\textsc{ltr}}$,
and $\operatorname{cap}(e)\doteq|V|$ otherwise. In this case, translating a
feasible solution to truck routes is not straightforward. It involves
partitioning the flow given by $X$ into $(n_{s},\tilde{n}_{s})$-directed
paths, which can be achieved in polynomial time.

##  5 Alternative constraints and valid inequalities

This section describes some alternative constraints and valid inequalities
that cut fractional solutions, which will later be used as cutting planes. For
simplicity, they are expressed for the LT formulation, but they can be adapted
to the others, except in specific cases that will be commented on in due
course.

#### Precedence inequalities

They are inspired by an alternative definition of the forbidden paths
mentioned in Section 3.1. Let $P$ be a $(n_{s},\tilde{n}_{s})$-directed path,

  * •

Disordered services (alternative). For all $r\in R$ and $P^{\prime}$
subdirected path of $P$ from $n_{s}$, the difference between the total number
of arcs of $E^{\textsc{d},r}$ and $E^{\textsc{p},r}$ in $P^{\prime}$ cannot be
greater than zero.

This alternative definition now excludes, for example, a path that has an arc
of $E^{\textsc{p},r}$ followed by two arcs of $E^{\textsc{d},r}$. Despite
being more restrictive than the original, they coincide when $P$ has no
repeated services. At the same time, it is unnecessary to consider all arcs of
$E^{\textsc{p},r}$ in $P^{\prime}$, but only those that give the truck enough
time to trip to the delivery location. Then, the following valid inequalities
arise:

$\displaystyle\hskip 10.0pt\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{p},r}_{\textsc{lt}}:\\\\[2.0pt]
\operatorname{\mathcal{I}^{\prime}}(e)\leq
i-\operatorname{len^{\textsc{t}}}(l^{\textsc{p}}_{r},l^{\textsc{d}}_{r})\end{subarray}}}\hskip
10.0pt\sum_{v\in V}X_{ve}\geq\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{d},r}_{\textsc{lt}}:\\\\[2.0pt] \operatorname{\mathcal{I}}(e)\leq
i\end{subarray}}}\hskip 10.0pt\sum_{v\in V}X_{ve}$ $\displaystyle\forall r\in
R,\ i\in\mathcal{I}^{\textsc{d}}_{r}.$ (PREC)

An additional advantage is that they can replace constraints (6), yielding a
tighter ILP formulation.

#### Synchronisation constraints

The synchronisation of routes requires one or two drivers to be present during
truck activities, as proposed in Section 3.3. While it may seem logical for
drivers to avoid taking taxis by travelling as passengers, there is no
advantage to sharing a truck in other circumstances. Observe that during
loading or unloading, one of the drivers may always be forced to remain in the
truck and the other to rest outside, which preserves the cost and feasibility
of the solution. The following constraints can achieve this restriction on the
solution space:

$\displaystyle\sum_{v\in V}X_{ve}$ $\displaystyle=\sum_{d\in D}Y_{de}$
$\displaystyle\forall e\in E^{\textsc{p}}_{\textsc{ltx}}\cup
E^{\textsc{d}}_{\textsc{ltx}}.$ (SYNC1)

The opposite approach involves allowing an unlimited number of drivers in the
truck during loading and unloading. In this case, this relaxation is achieved
by the following constraints:

$\displaystyle\sum_{v\in V}X_{ve}$ $\displaystyle\leq\sum_{d\in D}Y_{de}$
$\displaystyle\forall e\in E^{\textsc{p}}_{\textsc{ltx}}\cup
E^{\textsc{d}}_{\textsc{ltx}}.$ (SYNC2)

The constraints (11) on pickup and delivery arcs can be replaced by any of the
above. In particular, SYNC1 eliminate some feasible solutions but keep some
optimal ones, and SYNC2 introduce additional feasible solutions but an optimal
one can always be recovered for the original problem.

#### Pickup-and-delivery trip inequalities

Their motivation comes from observing that trucks always travel from the
pickup location to the delivery location (not necessarily directly) for each
request they serve. The first family of inequalities forces a truck route to
have as many trip arcs with origin (resp. destination) in a given location as
requests it loads (resp. unloads) there.

$\displaystyle\sum_{\begin{subarray}{c}e\in
E^{\textsc{trip}}_{\textsc{lt}}:\\\
\operatorname{L}(e)=l\end{subarray}}X_{ve}$
$\displaystyle\geq\sum_{\begin{subarray}{c}r\in R:\\\
l^{\textsc{p}}_{r}=l\end{subarray}}\hskip 2.0pt\sum_{e\in
E^{\textsc{p},r}_{\textsc{lt}}}X_{ve}$ $\displaystyle\forall l\in L,\ v\in V.$
(PD1)

To avoid making the description too long, similar inequalities for the
deliveries are also included. However, it is important to note that such trips
always occur after loading (resp. before unloading). By restricting the arcs
according to their start time, the following alternative inequalities are
obtained.

$\displaystyle\hskip 5.0pt\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{trip}}_{\textsc{lt}}:\\\ \operatorname{L}(e)=l,\
\operatorname{\mathcal{I}}(e)\geq i+s^{\textsc{p}}_{r}\end{subarray}}}\hskip
5.0ptX_{ve}$ $\displaystyle\geq\sum_{\begin{subarray}{c}r\in R:\\\
l^{\textsc{p}}_{r}=l\end{subarray}}\hskip 2.0pt\sum_{\begin{subarray}{c}e\in
E^{\textsc{p},r}_{\textsc{lt}}:\\\ \operatorname{\mathcal{I}}(e)\geq
i\end{subarray}}X_{ve}$ $\displaystyle\forall l\in L,\
i\in\bigcup_{\begin{subarray}{c}r\in R:\\\
l^{\textsc{p}}_{r}=l\end{subarray}}\mathcal{I}^{\textsc{p}}_{r},\ v\in V.$
(PD2)

Observe that PD2 include PD1 when $i$ is the earliest time instant. The last
family of these inequalities requires the truck to travel from the pickup
location (resp. to the delivery location) both after loading and before
unloading a given request. Instead of considering a single truck, the
following inequalities are generalised over a subset $V^{\prime}$ of trucks.

$\displaystyle\sum_{v\in V^{\prime}}\hskip
5.0pt\Bigl{(}\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{p},r}_{\textsc{lt}}:\\\
\operatorname{\mathcal{I}}(e_{1})\leq\operatorname{\mathcal{I}}(e)\end{subarray}}}X_{ve}+\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{d},r}_{\textsc{lt}}:\\\
\operatorname{\mathcal{I}^{\prime}}(e)\leq\operatorname{\mathcal{I}^{\prime}}(e_{2})\end{subarray}}}X_{ve}\Bigr{)}-1$
$\displaystyle\leq\sum_{v\in V^{\prime}}\hskip
7.0pt\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{trip}}_{\textsc{lt}}:\\\ \operatorname{L}(e)=l^{\textsc{p}}_{r},\\\
\operatorname{\mathcal{I}^{\prime}}(e_{1})\leq\operatorname{\mathcal{I}}(e),\
\operatorname{\mathcal{I}^{\prime}}(e)\leq\operatorname{\mathcal{I}}(e_{2})\end{subarray}}}X_{ve}$
$\displaystyle\forall~{}r\in R,\ e_{1}\in E^{\textsc{p},r},\ e_{2}\in
E^{\textsc{d},r},\ V^{\prime}\subset V.$ (PD3)

They should be taken into consideration only when
$\operatorname{\mathcal{I}^{\prime}}(e_{1})+\operatorname{len^{\textsc{t}}}(l^{\textsc{r}}_{p},l^{\textsc{e}}_{p})\leq\operatorname{\mathcal{I}}(e_{2})$.
The validity follows from the fact that $r$ can be served by at most one truck
in $V^{\prime}$ since repeated services are prohibited, thus one is an upper
bound for the left-hand side.

In their adaptations to the LTC formulation, a significant optimization can be
applied. Since the trip from the pickup location (resp. to the delivery
location) must occur while the truck is loaded, the set
$E^{\textsc{trip}}_{\textsc{ltc}}$ in the summation can be tightened to
$E^{\textsc{trip},{\bullet}}_{\textsc{ltc}}$, i.e. those trip arcs connecting
nodes with the superscript $\bullet$. Regarding the LTR formulation, PD1 and
PD2 can be adapted only for instances with a single truck and PD3 when
$V^{\prime}=V$, since the variables $X$ are no longer indexed by trucks.

#### Sequencing inequalities

There exists a minimum time for a truck to deliver a given set of requests.
Although it is a naive estimation since the restrictions on drivers are
ignored, it can be used as a lower bound on the route duration. For example,
the truck needs to (i) travel from its start location to the pickup location
of the first request, (ii) wait until the pickup time window opens, (iii) load
the cargo, (iv) travel to the delivery location, (v) wait until the delivery
time window opens, (vi) unload the cargo, (vi) travel to the pickup location
of the second request, and so on until delivering the last request. Observe
that all possible sequences of requests must be considered since their order
is not known in advance. Despite the computation of this parameter involves
solving a Pickup-and-Delivery Problem with Time Windows (see Dumas et al.,,
1991), few requests are tractable even by brute-force algorithms.

Given $R^{\prime}\subset R$ with $|R^{\prime}|\geq 2$ and $v\in V$, let
$\text{dur}^{\textsc{d}}_{0}(R^{\prime},v)$ be the minimum time needed by $v$
to deliver all requests in $R^{\prime}$. By noting that $v$ can only deliver
$|R^{\prime}|-1$ requests of $R$ in a time less than
$\text{dur}^{\textsc{p}}_{0}(R^{\prime},v)$, the following inequalities are
obtained.

$\displaystyle\sum_{r\in R^{\prime}}\hskip
10.0pt\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{d},r}_{\textsc{lt}}:\\\
\operatorname{\mathcal{I}^{\prime}}(e)<\text{dur}^{\textsc{d}}_{0}(R^{\prime},v)\end{subarray}}}X_{ve}$
$\displaystyle\leq|R^{\prime}|-1$ $\displaystyle\forall\ R^{\prime}\subseteq
R,\ v\in V.$ (SEC1)

Again, to avoid making the description too long, similar inequalities are
included for the pickups (i.e. considering the minimum time needed by $v$
until picking up the last request of $R^{\prime}$). Notice that this approach
spans from the beginning of the planning horizon. Alternative inequalities can
be proposed when considering an arbitrary $i\in\mathcal{I}$ as the beginning.
In this case, let $\text{dur}^{\textsc{d}}(R^{\prime},i)$ be the minimum time
a truck needs to deliver all requests in $R^{\prime}$ beginning at $i$. It is
important to remark that, as the location of a truck at a particular instant
is unknown in advance, it does not count the initial trip to the pickup
location of the first request, and as a consequence, is independent of the
truck.

$\displaystyle\sum_{r\in R^{\prime}}\hskip
10.0pt\sum_{\mathclap{\begin{subarray}{c}e\in
E^{\textsc{d},r}_{\textsc{lt}}:\\\ \operatorname{\mathcal{I}}(e)\geq i,\\\
\operatorname{\mathcal{I}^{\prime}}(e)<i+\text{dur}^{\textsc{d}}(R^{\prime},i)\end{subarray}}}X_{ve}$
$\displaystyle\leq|R^{\prime}|-1$ $\displaystyle\forall\ R^{\prime}\subseteq
R,\ v\in V,\ i\in\mathcal{I}.$ (SEC2)

Again, they can be adapted to the LTR formulation only for instances with a
single truck.

##  6 Computational experiments

In this section, the performance of the previous ILP formulations is compared
and the incorporation of the proposed valid inequalities as cutting planes is
evaluated, using a B&C algorithm from a commercial ILP solver. Experiments are
performed by a computer equipped with Ubuntu 18.04 64bits, 6 GB of memory, and
a core i7-9700 at 3.00 GHz. It is common to use the benchmark program _fmax_
and the benchmark instance _r500.5_ (both available in Trick, (2002)) to allow
future comparisons of our results with different machines. In our case, 2.86 s
(user time) were spent. The algorithms are implemented in C++11 and call CPLEX
Optimization Studio 20.1.0 through the Concert Technology API (IBM,, 2021).
The CPLEX parameters have the default values, except a time limit of 2 h per
instance, a single thread, and a deterministic search mode.

For this purpose, the following sets of random instances are generated, based
on the road network depicted in Fig. 9. The generation considers a 1-hour
discretisation of the planning horizon, a speed of 90 km/h for trucks and
taxis, a travel cost for trucks equal to the travel time, and a travel cost
for taxis twice that of trucks (since taxis must return to their origin).

  * •

S1: 15 instances with 3 locations ($L=\\{1,2,3\\}$), $|P|\in\\{4,5,6\\}$,
$H=7$, $|V|=1$ and $|D|=2$.

  * •

S2: 15 instances with $|P|\in\\{7,8,9\\}$ and the other parameters as S1.

  * •

S3: 15 instances with 3 locations, $|P|\in\\{8,9,10\\}$, $|V|=2$, and $|D|=4$.

  * •

S4: 15 instances with 6 locations and the other parameters as S1.

  * •

S5: 4 instances with 3 locations, $|P|\in\\{40,42\\}$, $H\in\\{40,42\\}$,
$|V|=1$, and $|D|=3$.

Figure 9: Fragment of the Argentine road network; extracted from Google Maps.

#### First experiment

For the first experiment, the ILP formulations, namely LT, LTC, and LTR,
presented in Section 4, are directly tackled with CPLEX. The results are shown
in Table 2, where the first columns report the name of the formulation, the
average number of variables and constraints, and the total number of solved
instances (to optimality). After that, the average execution time in seconds
over all instances (7200 s is considered for unsolved instances) and over the
solved instances are displayed. Then, the average relative gap over all
instances (0 % is considered for solved instances) and over the unsolved
instances are reported. The last column shows the average number of explored
nodes over the solved instances. The performance profile of the algorithms is
depicted in Fig. 10(a), using the execution time as the performance measure
(see Dolan and Moré,, 2002). Fig. 10(b) shows, for each instance, the relative
gap between the optimal value of the ILP formulation and its linear relaxation
(LR). It is worth mentioning that the optimal value is known for all the
instances of S1, as will be shown in the second experiment. The average
execution time of the LRs is 2.17 s (LT), 2.88 s (LTC), and 5.17 s (LTR).

Table 2: Comparison of ILP formulations when solving instances of S1 (part I). Form. | Var. | Cons. | Solved | Time (s) | Gap (%) | Nodes
---|---|---|---|---|---|---
|  |  |  | All | Sol. | All | Unsol. | Sol.
LT | 8682 | 5813 | 1 | 7109 | 5828 | 49 | 53 | 139372
LTC | 10166 | 6151 | 8 | 4597 | 2319 | 14 | 33 | 9154
LTR | 16102 | 7828 | 14 | 894 | 444 | 0.8 | 11 | 1915

(a) Performance profiles of the algorithms.

(b) Gap between the optimal value of the formulations and their LRs.

Figure 10: Comparison of ILP formulations when solving instances of S1 (part
II).

A big difference is observed in the behaviour of the algorithms, mainly in the
number of instances solved, where LTR far outperforms the others and LTC far
outperforms LT. The results show that the formulations are more competitive as
the digraphs have greater structure, despite the increase in the number of
variables and constraints. At the same time, the LRs are tighter as the
digraphs have greater structure, which might explain the difference observed.
Consequently, the next experiments focus on improving the best two algorithms
(based on LTC and LTR) to solve larger instances.

#### Second experiment

The next experiment evaluates the alternative formulations for LTC and LTR
mentioned at the beginning of Section 5. Their names are augmented with the
family of constraints being replaced, e.g. in “LTC + PREC + SYNC1”,
constraints (6) and those in (11) for pickup and delivery arcs are replaced
with PREC and SYNC1, respectively. It is worth remembering that the former
family is already implied in LTR, so LTR + PREC does not make sense in the
current experiment. The results are shown in Table 3, where the average number
of variables can be seen from the previous table. Fig. 11(a) presents the
performance profiles, where some of the less performant formulations are not
shown to facilitate reading. The LRs are analysed in Fig. 11(b), where
constraints SYNC1 are SYNC2 are omitted because they preserve the gaps. The
average execution time of the LRs is 2.57 s for LTC + PREC.

Table 3: Comparison of alternative ILP formulations when solving instances of S1 (part I). Form. | Cons. | Solved | Time (s) | Gap (%) | Nodes
---|---|---|---|---|---
|  |  | All | Sol. | All | Unsol. | Sol.
LTC | 6151 | 8 | 4597 | 2319 | 14 | 33 | 9154
LTC + PREC | 6151 | 14 | 1274 | 851 | 1.5 | 22 | 7151
LTC + SYNC1 | 5404 | 8 | 4674 | 2463 | 12 | 31 | 105495
LTC + SYNC2 | 5404 | 8 | 4597 | 2319 | 7.3 | 16 | 168612
LTC + PREC + SYNC1 | 5404 | 14 | 1019 | 577 | 2.6 | 39 | 4611
LTC + PREC + SYNC2 | 5404 | 12 | 2201 | 951 | 0.6 | 3.1 | 54125
LTR | 7828 | 14 | 894 | 444 | 0.8 | 11 | 1915
LTR + SYNC1 | 7082 | 15 | 917 | 917 | 0 | – | 2318
LTR + SYNC2 | 7082 | 14 | 2057 | 1690 | 1.3 | 19 | 56299

(a) Performance profiles of the algorithms.

(b) Gap between the optimal value of the formulations and their LRs.

Figure 11: Comparison of alternative ILP formulations when solving instances
of S1 (part II).

A remarkable improvement is noticed, mainly in the number of solved instances,
when LTC uses PREC, which greatly increases its competitiveness against LTR.
Again, the tighter LRs might explain this behaviour. The best algorithms are
based on LTC + PREC + SYNC1, LTR, and LTR + SYNC1. Focusing on solving the
greatest number of instances, we decided to continue working with LTC + PREC +
SYNC1 and LTR + SYNC1, which from now on will directly be called LTC and LTR,
respectively.

#### Third experiment

The following experiment analyses the impact of performing a warm start, i.e.
providing the ILP solver with an external initial solution. In our
implementation, we use the metaheuristics developed in Lucci et al., (2023) to
rapidly generate a heuristic solution, which is translated into a feasible
solution for the integer linear programs. For each instance of S1, the average
execution time of this initialisation is 7 s, and for the other sets it
requires 18 s (S2 and S3), 8 s (S4), and 464 s (S5).

The results obtained for S1 are discussed below. The warm start does not
modify the number of solved instances but impacts the execution time. In
particular, LTR is notably faster when performing a warm start, whose average
execution time over the solved instances decreases from 917 s to 647 s. This
is less noticeable for LTC, from 577 s to 484 s, but the gap for the only
unsolved instance decreases from 39% to 12%. To accompany the analysis, Fig.
12(a) depicts the performance profiles, and Fig. 12(b) shows, for each
instance, the gap between the optimal value and the value of the initial
solution provided. Since the warm start has proved to be useful, it will
continue to be enabled in the remainder of this work.

(a) Performance profiles of the algorithms.

(b) Gap between the optimal value and the initial value.

Figure 12: Impact of a warm start when solving instances of S1.

#### Fourth experiment

It remains to evaluate the impact of the valid inequalities presented in
Section 5 for the tightening of the LRs. In this experiment, we focus on the
root node of the B&C tree, while the remaining nodes are left for later. The
proposed methodology considers each family separately and incorporates all the
valid inequalities into the models before solving their LRs.

It is worth noting that PD3, SEC1, and SEC2 grow exponentially with the size
of the input. For this reason, a maximum size is introduced for the subsets of
trucks ($V^{\prime}$) and requests ($R^{\prime}$). The names of the
inequalities are augmented accordingly, e.g., for $k\in\mathbb{N}$,
“PD3-$V_{k}$” considers all the inequalities of PD3 such that
$|V^{\prime}|\leq k$ and similarly “SEC1-$R_{k}$” demands $|R^{\prime}|\leq
k$. However, this restriction might not be sufficient, and additional
subfamilies are considered for PD3 by limiting the initial time for loading
and unloading. In “PD3-$V_{k}$-A”, the start time of $e_{1}$ must coincide
with the opening of a pickup time window and that of $e_{2}$ with the closing
of a delivery time window. For example, for time windows as Fig. 2a, $e_{1}$
needs to start at 8:00 (on day 0 or 1) and $e_{2}$ at 20:00 (on day 0 or 1).
On the other hand, “PD3-$V_{k}$-B” is similar but with the conjunction of both
conditions, e.g. $e_{1}$ at 8:00 or $e_{2}$ at 20:00.

Tables 4 and 5 show the results for the instances of S2 and S3, which are
larger than S1. The columns report the average number of valid inequalities,
the average execution time, and the average relative gap with respect to the
value of the best incumbent solution for the integer programs (which will be
discovered later). Families that have reduced the gap little or nothing (less
than 0.2 percentage points) are excluded, e.g. PREC and PD for LTR in Table 4,
SEC for LTC in Table 5. Families that cannot be adapted are also omitted, e.g.
most of the families for LTR in Table 5. It is evident that the PD
inequalities are highly effective for LTC. In particular, the gap achieved by
LTC + PD2 is only 0.5 percentage points worse than LTR for the instances of
S2. Fig. 13(b) exhibits a more detailed comparison of the gaps obtained for
each instance of S1. When the instances have more than one truck, the best
value for $k$ in PD3-$V_{k}$ seems to be $k=|V|$, but they are still
outperformed by PD2. Regarding the SEC inequalities, larger values of $k$
obtain better gaps, but the decrease is less abrupt from $k>4$, while the
number of inequalities and the execution times increase in counterpart.

Table 4: Impact of valid inequalities on the root node for instances of S2 (part I). Form. | #Ineq. | Time (s) | Gap (%) | Algo. | #Ineq. | Time (s) | Gap (%)
---|---|---|---|---|---|---|---
LTC | 0 | 4.0 | 30.4 | LTR | 0 | 8.1 | 6.7
LTC + PD1 | 6 | 7.4 | 8.0 | LTR + SEC1-$R_{2}$ | 30 | 8.7 | 6.7
LTC + PD2 | 686 | 9.5 | 7.2 | LTR + SEC1-$R_{3}$ | 107 | 9.4 | 6.1
LTC + PD3-$V_{1}$-A | 349 | 7.1 | 9.8 | LTR + SEC1-$R_{4}$ | 223 | 8.2 | 5.5
LTC + PD3-$V_{1}$-B | 8531 | 54.3 | 8.2 | LTR + SEC1-$R_{5}$ | 335 | 8.3 | 5.4
LTC + SEC1-$R_{2}$ | 30 | 3.9 | 30.2 | LTR + SEC1-$R_{6}$ | 406 | 8.5 | 5.4
LTC + SEC1-$R_{3}$ | 107 | 4.1 | 29.4 | LTR + SEC2-$R_{2}$ | 6934 | 13.7 | 6.2
LTC + SEC1-$R_{4}$ | 223 | 4.5 | 28.9 | LTR + SEC2-$R_{3}$ | 23233 | 30.7 | 6.0
LTC + SEC1-$R_{5}$ | 335 | 4.6 | 28.7 | LTR + SEC2-$R_{4}$ | 45670 | 33.4 | 5.2
LTC + SEC1-$R_{6}$ | 406 | 4.9 | 28.7 | LTR + SEC2-$R_{5}$ | 65948 | 36.9 | 5.1
LTC + SEC2-$R_{2}$ | 6934 | 6.5 | 29.8 | LTR + SEC2-$R_{6}$ | 65948 | 36.9 | 5.1
LTC + SEC2-$R_{3}$ | 23233 | 11.8 | 29.3 |
LTC + SEC2-$R_{4}$ | 45670 | 14.1 | 28.4 |  |  |  |
LTC + SEC2-$R_{5}$ | 65948 | 19.0 | 28.2 |  |  |  |
LTC + SEC2-$R_{6}$ | 78140 | 22.6 | 28.1 |  |  |  |
Table 5: Impact of valid inequalities on the root node for instances of S3 (part II). Form. | Ineq. | Time (s) | Gap (%)
---|---|---|---
LTC | 0 | 7.4 | 39
LTC + PD1 | 12 | 14.6 | 24
LTC + PD2 | 1403 | 24.9 | 5.6
LTC + PD3-$V_{1}$-A | 789 | 10.3 | 38
LTC + PD3-$V_{1}$-B | 18273 | 142.3 | 38
LTC + PD3-$V_{2}$-A | 1184 | 17.4 | 15
LTC + PD3-$V_{2}$-B | 27409 | 147.7 | 12

Next, some of the best-performing families at the root node are selected and
evaluated on the entire tree. Small families are managed as an initial part of
the constraint set of the models, thus all the valid inequalities are present
in each node. Otherwise, larger families are grouped into a pool of cuts and
added to each node by the ILP solver on demand.

The results for PD inequalities are reported in Table 6 and Fig. 13(a). In
addition, the column “Mode” tells how the valid inequalities are handled, i.e.
as part of the initial constraint set (“Init.”) or through a pool of cuts
(“Pool”), and “Cuts” is the average number of cuts added from the pool per
instance. It is clear that PD2, PD3-$V_{2}$-$A$, and PD3-$V_{2}$-$B$ conduct
to a higher number of solved instances, lower execution times, and tighter
gaps; in particular, the former slightly outperforms the last two.

Table 6: Impact of PD inequalities on the entire tree for instances of S2, S3, and S4 (part I). Algo. | Mode | Ineq. | Cons. | Solved | Time (s) | Gap (%) | Nodes | Cuts
---|---|---|---|---|---|---|---|---
|  |  |  |  | All | Sol. | All | Unsol. | Sol. |
LTC | – | 0 | 10100 | 17 | 5474 | 2632 | 18 | 29 | 8415 | –
LTC + PD1 | Init. | 10 | 10110 | 18 | 5330 | 2525 | 13 | 21 | 42325 | –
LTC + PD2 | Init. | 912 | 11013 | 26 | 4459 | 2456 | 6.3 | 15 | 25916 | –
LTC + PD3-$V_{2}$-A | Init. | 587 | 10687 | 22 | 4632 | 1946 | 9.5 | 19 | 6206 | –
LTC + PD3-$V_{2}$-B | Pool | 13837 | 10100 | 23 | 4614 | 2140 | 8.5 | 17 | 5065 | 379

(a) Performance profiles of the algorithms.

(b) Gap between the optimal value of the formulations and their LRs.

Figure 13: Impact of PD inequalities on the entire tree for instances of S2,
S3, and S4 (part II).

While keeping PD2 as part of the constraint set of LTC, SEC inequalities are
now evaluated on the entire tree. The results are reported in Table 7 for the
instances of S2 and S4 (S3 is not considered since these inequalities were not
effective on the root node when there is more than one truck). The performance
profiles of the algorithms based on LTC and LTR can be seen in Fig. 14(a) and
14(b), respectively. The results are inconclusive for LTC as no algorithm
seems to outperform the others, e.g. the number of solved instances decreases.
Instead for LTR, SEC1-$R_{4}$ (Pool) seems superior as the number of solved
instances increases and the performance profile has the most wins within a
ratio of 1.5 or higher. It is interesting to mention that, although
SEC1-$R_{4}$ (Init) has the most wins within a ratio of 1, the performance
decays notoriously as the ratio increases. In addition, very few cuts are
added from the pool.

Table 7: Impact of SEC inequalities on the entire tree for instances of S2 and S4 (part I). Algo. | Mode | Ineq. | Var. | Cons. | Solved | Time (s) | Gap (%) | Nodes | Cuts
---|---|---|---|---|---|---|---|---|---
|  |  |  |  |  | All | Sol. | All | Unsol. | Sol. |
LTC | – | 0 | 18906 | 10788 | 20 | 3989 | 2384 | 7.4 | 22 | 28398 | –
LTC + SEC1-$R_{4}$ | Init. | 129 | 18906 | 10917 | 19 | 3882 | 1961 | 7.7 | 21 | 13447 | –
LTC + SEC1-$R_{4}$ | Pool | 129 | 18906 | 10788 | 19 | 4201 | 2464 | 7.9 | 22 | 19197 | 0.5
LTC + SEC2-$R_{4}$ | Pool | 26427 | 18906 | 10788 | 19 | 4049 | 2225 | 8.1 | 22 | 13601 | 3.1
LTR | – | 0 | 32664 | 13477 | 23 | 2652 | 1268 | 3.6 | 15 | 4834 | –
LTR + SEC1-$R_{4}$ | Init. | 129 | 32664 | 13606 | 22 | 2963 | 1422 | 2.9 | 11 | 42406 | –
LTR + SEC1-$R_{4}$ | Pool | 129 | 32664 | 13477 | 24 | 2452 | 1265 | 2.4 | 12 | 4008 | 0.5
LTR + SEC2-$R_{4}$ | Pool | 26427 | 32664 | 13477 | 22 | 2548 | 857 | 3.2 | 12 | 3537 | 2.3

(a) Performance profiles of the algorithms based on LTC.

(b) Performance profiles of the algorithms based on LTR.

Figure 14: Impact of SEC inequalities on the entire tree for instances of S2
and S4 (part II).

#### Fifth experiment

As a final experiment, the best algorithms for LTC and LTR are compared
against each other, considering the instances of S2, S3, and S4. To summarise
the results, LTR shows a clear superiority in the number of solved instances,
37 vs. 26. Every instance solved by LTC is also solved by LTR and more quickly
(except two of them), with an average execution time of 734s vs. 2250s.
Furthermore, when none of them prove optimality, LTR reports better gaps
(except for one instance), with an average gap of 12% vs. 19%. In 9 of the 19
unsolved instances by LTC, the final lower bound is at least 5% worse than
LTR, and in 3 of them (all from S4), at least 24%.

Then, the limits of the algorithms are analised using the instances of S5,
which have the greatest number of requests and the longest planning horizon.
In this opportunity, we migrate to a more powerful computer, using the 16
cores of an i7-10700 at 2.9GHz (parallel optimization), 6 GB of memory, and a
time limit of 10 hours. On the one hand, LTC has an average of 158k variables
and 95k constraints. The LRs are successfully solved on the root node but no
other nodes are explored within the time limit. The LRs are tightened with 3
to 5 rounds of general-purpose cutting planes and the final gaps vary from
6.7% to 18.5%. On the other hand, LTR runs out of memory while writing the
models for 3 instances and during the presolve for the other, with 485k
variables and 186k constraints.

Finally, we reconsider the first algorithm based on LTC, without alternative
constraints or valid inequalities, but keeping the warm start (with the same
initial solutions) for the instances of S5. Again, only the LRs are solved on
the root node and there is not enough time to branch. Although more rounds of
cutting planes are performed, 7 to 11, the final gaps deteriorate greatly,
varying from 69.8% y 79.4%.

##  7 Conclusions

In this work, we addressed a SVRCSP in long-haul transport. The goal was to
minimise the travel costs and the delay penalties involved in transporting
pickup-and-delivery requests with time windows over a multi-day planning
horizon. Many real-world requirements were considered, such as the hour-of-
service regulation of the drivers. Unlike most approaches in the literature,
the correspondence between trucks and drivers was not fixed and could be
exchanged in some locations. Different digraphs were defined to represent
truck and driver routes separately as directed paths, where nodes were indexed
by location and time. Some of them incorporated greater structure (additional
nodes and arcs) to ensure that certain forbidden paths could not exist by
construction. Three compact ILP formulations were proposed to model the
problem based on these digraphs, namely LT, LTC, and LTR. These models allowed
truck and driver routes to be easily synchronised in time and space. Among the
most notable differences, the number of variables of LT and LTC strongly
depended on the number of trucks, whereas LTR on the number of requests. Many
families of valid inequalities were presented and some of them were
exponential in size.

Extensive computational experiments were performed on random instances to
compare the ILP formulations and evaluate specific valid inequalities as
cutting planes. One of the most interesting conclusions was that, as the
digraphs gained structure and the formulations contained fewer families of
constraints, the linear relaxations were tighter and the number of solved
instances increased. Furthermore, tightening LTC with precedence inequalities
and pickup-and-delivery trip inequalities produced notable performance
improvements, primarily in the number of solved instances. In contrast, the
algorithm based on LTR, which initially performed best, supported notable
reductions in execution time with the addition of a warm start and a pool of
cuts with the sequencing inequalities. A final comparison revealed that the
latter significantly outperformed the former in instances with 3-6 locations,
7-10 requests, 1-2 trucks, 2-4 drivers, and a 1-week planning horizon.
However, the algorithm based on LTC stood out for providing relative gaps for
larger instances, with 40-42 requests and a planning horizon of 40-42 days.

In conclusion, this study has contributed to the development of exact methods
for solving SVRCSPs. To our knowledge, this is the first time that solutions
with proven optimality have been obtained for the case study. We believe that
these methods can be extended to other similar studies, but potential
adaptations might be necessary to handle different hour-of-service
regulations. This research also suggested some topics that would be
interesting to study in the future. Relaxing some of the requirements might
achieve further cost improvements, e.g. allowing drivers to share taxis.
Another line of research is the definition of new valid inequalities, which
might help to reduce the differences in the performance of the algorithms,
considering the limitation of LTR when distinguishing trucks. Regarding the
implementation part, it would be interesting to develop specific separation
routines for some families of valid inequalities, as well as primal
heuristics. Finally, since the models generally have many more variables than
constraints, a branch-and-price scheme could be competitive to address larger
instances.

## Acknowledgments

This work was partially supported by grants PIP-1900 (CONICET) and
PICT-2020-03032 (ANPCyT).

## References

  * Ammann et al., (2023) Ammann, P., Kolisch, R., and Schiffer, M. (2023).  Driver routing and scheduling with synchronization constraints.  Transportation Research Part B: Methodological, 174:102772.
  * Deveci and Çetin Demirel, (2018) Deveci, M. and Çetin Demirel, N. (2018).  A survey of the literature on airline crew scheduling.  Engineering Applications of Artificial Intelligence, 74:54–69.
  * Dolan and Moré, (2002) Dolan, E. and Moré, J. (2002).  Benchmarking optimization software with performance profiles.  Mathematical Programming, 91(2):201–213.
  * Domínguez-Martín et al., (2018) Domínguez-Martín, B., Rodriguez-Martin, I., and Salazar-Gonzalez, J.-J. (2018).  The driver and vehicle routing problem.  Computers & Operations Research, 92:56–64.
  * Drexl, (2012) Drexl, M. (2012).  Synchronization in Vehicle Routing - A Survey of VRPs with Multiple Synchronization Constraints.  Transportation Science, 46(3):297–316.
  * Drexl et al., (2013) Drexl, M., Rieck, J., Sigl, T., and Press, B. (2013).  Simultaneous vehicle and crew routing and scheduling for partial- and full-load long-distance road transport.  Business Research, 6(2):242–264.
  * Dumas et al., (1991) Dumas, Y., Desrosiers, J., and Soumis, F. (1991).  The pickup and delivery problem with time windows.  European Journal of Operational Research, 54(1):7–22.
  * Goel, (2009) Goel, A. (2009).  Vehicle scheduling and routing with drivers’ working hours.  Transportation Science, 43(1):17–26.
  * Goel and Irnich, (2017) Goel, A. and Irnich, S. (2017).  An exact method for vehicle routing and truck driver scheduling problems.  Transportation Science, 51(2):737–754.
  * Goel and Vidal, (2014) Goel, A. and Vidal, T. (2014).  Hours of service regulations in road freight transport: An optimization-based international assessment.  Transportation Science, 48(3):391–412.
  * Goel et al., (2021) Goel, A., Vidal, T., and Kok, A. L. (2021).  To team up or not: single versus team driving in European road freight transport.  Flexible Services and Manufacturing Journal, 33(4):879–913.
  * Heil et al., (2020) Heil, J., Hoffmann, K., and Buscher, U. (2020).  Railway crew scheduling: Models, methods and applications.  European Journal of Operational Research, 283(2):405–425.
  * Ibarra-Rojas et al., (2015) Ibarra-Rojas, O., Delgado, F., Giesen, R., and Muñoz, J. (2015).  Planning, operation, and control of bus transport systems: A literature review.  Transportation Research Part B: Methodological, 77:38–75.
  * IBM, (2021) IBM (2021).  ILOG CPLEX Optimization Studio 20.1.0.  https://www.ibm.com/products/ilog-cplex-optimization-studio.  Accessed 29 January 2024.
  * Koç et al., (2018) Koç, Ç., Jabali, O., and Laporte, G. (2018).  Long-haul vehicle routing and scheduling with idling options.  Journal of the Operational Research Society, 69(2):235–246.
  * Kok et al., (2010) Kok, A. L., Meyer, C. M., Kopfer, H., and Schutten, J. M. J. (2010).  A dynamic programming heuristic for the vehicle routing problem with time windows and european community social legislation.  Transportation Science, 44(4):442–454.
  * Lam et al., (2015) Lam, E., Van Hentenryck, P., and Kilby, P. (2015).  Joint vehicle and crew routing and scheduling.  In Pesant, G., editor, Principles and Practice of Constraint Programming, pages 654–670, Cham. Springer International Publishing.
  * Lucci et al., (2022) Lucci, M., Severin, D., and Zabala, P. (2022).  Integer programs for a simultaneous vehicle routing and crew scheduling problem.  In Joint ALIO/EURO, pages 28–32.
  * Lucci et al., (2023) Lucci, M., Severín, D., and Zabala, P. (2023).  A metaheuristic for crew scheduling in a pickup-and-delivery problem with time windows.  International Transactions in Operational Research, 30(2):970–1001.
  * Mendes and Iori, (2020) Mendes, N. F. and Iori, M. (2020).  A decision support system for a multi-trip vehicle routing problem with trucks and drivers scheduling.  In ICEIS (1), pages 339–349.
  * Prescott-Gagnon et al., (2010) Prescott-Gagnon, E., Desaulniers, G., Drexl, M., and Rousseau, L.-M. (2010).  European driver rules in vehicle routing with time windows.  Transportation Science, 44(4):455–473.
  * Rancourt et al., (2013) Rancourt, M.-E., Cordeau, J.-F., and Laporte, G. (2013).  Long-haul vehicle routing and scheduling with working hour rules.  Transportation Science, 47(1):81–107.
  * Tilk and Goel, (2020) Tilk, C. and Goel, A. (2020).  Bidirectional labeling for solving vehicle routing and truck driver scheduling problems.  European Journal of Operational Research, 283(1):108–124.
  * Toth and Vigo, (2014) Toth, P. and Vigo, D. (2014).  Vehicle routing: problems, methods, and applications.  SIAM.
  * Trick, (2002) Trick, M. (2002).  Color02/03/04: Graph coloring and its generalizations.  https://mat.gsia.cmu.edu/COLOR04/.  Accessed 29 January 2024.