- Email: [email protected]
A branch-and-price based heuristic for the stochastic vehicle routing problem with hard time windows
Available online at www.sciencedirect.com
Electronic Notes in Discrete Mathematics 52 (2016) 325–332 www.elsevier.com/locate/endm
A branch-and-price based heuristic for the stochastic vehicle routing problem with hard time windows G. Andreatta 1, M. Casula 2, C. De Francesco 3, L. De Giovanni 4 Dipartimento di Matematica Universit` a degli Studi di Padova Padova, Italy
Abstract In the Vehicle Routing Problem with Hard Time Windows and Stochastic travel times, a disruption occurs if, due to stochastic events, a vehicle arrives too late at a customer. In this case a recourse action is required such that the service starts within the time window, and a relevant penalty cost is incurred. Despite the problem has been inspired by a real-life application in airport ground handling optimization, it has never been addressed before in literature, to the best of our knowledge. We discuss how the expected penalty cost can be evaluated and how this computation can be integrated in a branch-and-price procedure to obtain heuristic solutions. Preliminary tests on literature instances show the effectiveness of the approach. Keywords: Vehicle routing, Hard time windows, Stochastic travel times, Branch-and-price
1 2 3 4
Email: Email: Email: Email:
[email protected] [email protected] [email protected] [email protected]
http://dx.doi.org/10.1016/j.endm.2016.03.043 1571-0653/© 2016 Elsevier B.V. All rights reserved.
326
1
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
Introduction
The Vehicle Routing Problem (VRP) is a well known network optimization problem where a set of customers is located at the nodes of a graph and a set of node-disjoint routes have to be determined in order to visit all the customers. Routes have to start and end at a same node called depot and travel distances between nodes, as well as node service times, are given. The objective is to minimize the total cost, normally related to the total distance traveled (see [13], among others). One of the possible extensions of this problem is the VRP with time windows (VRPTW), where the service of each customer has to start within a specific time interval. Time windows may be soft (they can be violated at some penalty cost), or hard (early vehicles have to wait at no additional cost, while late vehicles are prohibited to serve). Real-life applications often require that uncertainty is taken into account so that stochastic variants of the VRP have been introduced (see e.g. [9] for a review). We consider the VRP with Hard Time Windows and Stochastic travel and service times (S-VRPHTW). Due to the uncertainty, a customer originally included in a route may be reached after the end of its time window: we call disruption this event, and we assume that a recourse action is taken such that the customer is served within its time window and a large fixed penalty cost is incurred. As a consequence, the objective is to minimize a weighted sum of the total distance traveled and the expected number of disruptions. The study of the S-VRPHTW has been inspired by a real application related to the optimization of ground handling operations in airport aprons, where the customers are the aircraft that need specific operations before/after departure/landing, and the vehicles are the specific Ground Service Equipment (GSE) required by the same operations (see [3]). If a GSE cannot serve an aircraft on time (for example because the previous aircraft is late, or the delay cumulated along the service route is sufficiently large, or because of a breakdown or other unforeseen events), a recourse GSE takes its place in the originally planned route, causing an additional relevant penalty cost. Among the stochastic variants of the VRP, the ones involving stochastic travel times are less studied and literature is concerned with soft time windows: we cite some heuristic approaches based on genetic algorithms or tabu search (see [2,10] among others), and an exact algorithm based on branch-and-price [12]. Also the method proposed in [6] for the VRP with stochastic demands may be adapted to soft time windows, under suitable restrictions on delay distributions. The S-VRPHTW is concerned with hard time windows and, to the best of our knowledge, it has not been yet considered in literature.
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
327
In this paper, we propose a first heuristic approach for the S-VRPHTW. Theoretical and practical issues related to the computation of the expected number of disruptions are addressed in Section 2, proposing an approximation and discussing how it can be integrated in a branch-and-price procedure to obtain heuristic solutions. Preliminary results on literature benchmarks show that the proposed approach is able to provide good solutions, trading off between computational efficiency and accuracy of the stochastic approximations, as outlined in Section 3, together with some lines for future research.
2
The evaluation of the expected number of disruptions and the branch-and-price procedure
Let G = (N, A) be a directed graph where N is the set of customers. For each arc (i, j) ∈ A, let tij be the travel time from customer i to customer j and cij the related cost. For each node v ∈ N, let rv be the service time and [av , bv ], with av ≤ bv the time window. The depot is represented by node v0 . Let p = v0 − v1 − . . . − vn be an elementary path on G, that is, a path without node repetitions. With relation to p, the starting time si of a node vi is recursively defined as follows: s0 = 0; si = max{ai , si−1 + ri−1 + ti−1 i } for every i = 1...n. Path p is feasible if vn = v0 and ai ≤ si ≤ bi , for every i = 0...n. We assume that travel and service times are subject to independent random delays δi j , for each (i, j) ∈ A. If P is the penalty cost associated to a disruption, the cost of p is defined as c(p) =
n−1
cvi vi+1 + P · E(μp )
(1)
i=0
where μp is a random variable representing the number of disruptions on p. We consider the following set-partitioning formulation for S-VRPHTW [7]: c(p)xp (2) min p∈Ω
s.t.
xp = 1
∀v∈N
(3)
∀p∈Ω
(4)
p∈Ω(v)
xp ∈ {0, 1}
where Ω is the set of all feasible paths, Ω(v) is the set of feasible paths containing v, and xp takes value one if p is included in the solution, 0 otherwise. For each path p we have to evaluate c(p), and hence E(μp ). If d(i) is the event
328
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
“a disruption occurs at the node vi of p”, then E(μp ) = ni=1 P [d(i)]. Clearly, P [d(i)] depends on the delay cumulated up to node vi . We notice that, because of the hard time windows, the recourse action is intended to absorb the effect of delays in case of disruption (for example, a recourse vehicle replaces the delayed one on the originally planned route). Let σi be the random variable denoting the random starting time for node vi . We have σ0 = 0 and, for i = 1...n mi = max{ai , σi−1 + ri−1 + ti−1 i + δi−1 i } if mi ≤ bi σi = (5) ui otherwise where ui is any value within the time window. We assume that the recourse action allows the service of customer i to start at the originally planned starting time, that is ui = si . By denoting with e(k, i) (for 0 ≤ k < i) the event “the first disruption after node vk occurs in vi (independently from the occurrence of other disruptions on the route from v0 to vk )”, we have P [d(i)] = P [e(0, i)] +
i−1
P [e(k, i)|d(k)] · P [d(k)]
(6)
k=1
Under the hypothesis ui = si , P [e(k, i)|d(k)] = P[(δi−1 i > bi − (ti−1 i + ri−1 + σi−1 ))∧ ∧ (δi−2 i−1 ≤ bi−1 − (ti−2 i−1 + ri−2 + σi−2 ))∧ ∧ . . . ∧ (δk k+1 ≤ bk+1 − (tk k+1 + rk + sk ))]
(7)
We notice that P [e(0, i)] can be evaluated as in (7) with k = 0 and that P [d(1)] = P [δ0 1 > b1 − t0 1 ], hence we can recursively obtain P [d(i)]. Formula (7) can be analytically integrated, at least in theory, starting from the p.d.f. of variables δ and considering the sequence of the realizations of the same variables from δk k+1 to δi−1 i (see [4] for details). In practice, the integration requires to consider a number of cases exponential in the length of the path. Thus we resort to an -lookback approximation to evaluate P [d(i)], which consists in considering only the last delays before vi , that is, we reformulate equation (6) assuming that the first node of the path is vk¯ , k¯ = max{0, i − } and σk¯ = sk¯ . In this way we have only a fixed (possibly small) number of integration cases, regardless of the actual length of the path. The idea behind this approximation is that, since the delay resets whenever the vehicle arrives outside a time window, considering the whole delay history may
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
329
be redundant. Computational experiments in some particular cases show that a good approximation can be obtained even with small values of . Formulation (2)...(4) contains an exponential number of variables and we adopt a branch-and-price approach, that is, a branch-and-bound procedure using column generation to solve the linear relaxation: we initially consider a small subset of feasible paths, and a pricing subproblem is iteratively solved to determine further feasible paths that should be included in the formulation to get better solutions. The components we implemented for solving S-VRPHTW are essentially the ones of a standard branch-and-price for VRPTW (see e.g. [7]). The pricing problem asks for a feasible path with minimum reduced cost c¯(p) = c(p) − vi ∈p π ¯i , where π ¯i are the current dual variables associated with (3). This can be seen as an Elementary Shortest Path Problem with Resource Constraint (ESPPRC), where resources are time (to guarantee that arrival times fall into the time windows) and vehicle capacity, where defined. We follow the approach of [8] where ESPPRC is solved by a labeling algorithm. We use labels L(p) to store the features of partial feasible paths p from v0 to vi . Each label stores: c¯(p); si ; the total capacity used by p; a binary value Aj for each j = 1 . . . n (Aj is equal to 1 if node vj cannot be visited in any extension of p, because it has been already visited, or due to resource constraints); and D(p), the expected delay at the end of p. Notice that labels are referred to paths. In fact, when we extend a path p to p using arc (i, j), we have to compute c¯(p ) and D(p ) which depend on the node sequence (and not only on (i, j) and cumulated information up to vi , as e.g. in [6]). Of course, storing the labels of every partial path is not efficient, and, at each node, we do not store dominated labels, that is, labels having, for each field, a non-smaller value than an already stored one. This dominance criterion would be exact if we were able to store the actual delay cumulated on partial paths (the same extension of a path with greater delay and greater reduced cost would lead to a path with a greater reduced cost, since the expected number of disruptions would be greater). However, we cannot store the actual delays, so that storing their expected values ends up with a heuristic procedure for pricing (and hence for S-VRPHTW).
3
Computational results and conclusions
Preliminary tests have been conducted to determine the accuracy of the lookback approximation in estimating E(μp ). To this aim, we consider “regular” paths p of length n having, for any node vi , ri + ti i+1 = 2, ai = 2 i, bi = 2 i + 1 and exponentially distributed delay with average D. In this case,
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
330
Table 1 Expected number of disruptions on a regular path with 30 nodes. D 0.1 0.2 0.3 0.7
=1 0.01 0.98 3.52 11.83
=2 0.06 2.41 5.94 12.47
=3 0.20 3.98 7.08 12.22
=4 0.50 5.03 7.07 12.16
=5 0.98 5.31 6.73 12.18
=6 1.58 5.07 6.55 12.18
=7 2.20 4.72 6.56 12.18
=8 2.70 4.53 6.61 12.18
E(μp ) 2.31 4.62 6.60 12.18
Table 2 Average running times (in seconds). Group R1 C1 RC1 R2 C2 RC2
closed 12 over 12 9 over 9 8 over 8 10 over 11 6 over 8 6 over 8
det. 1.5 84.9 2.0 151.6 15.6 39.7
L1 1.6 92.4 2.2 331.7 33.6 96.0
L2 13.1 165.9 17.3 469.1 71.3 137.6
L3 8.5 169.2 12.6 399.2 73.2 148.5
we have analytically evaluated E(μp ), both exactly and with -lookback approximation for different values of (see [4]). Sample results with n = 30 are reported in Table 1. In general, with increasing , the approximation becomes better and better and the convergence to the exact value is faster for larger D. Notice that even small values of give good approximations. The branch-and-price based heuristic has been tested on instances from the Solomon’s benchmark [11]. We consider independent exponentially distributed delays associated to each arc (i, j), with average 0.1 · (ri + tij ), and we set P equal to twice the average distance of nodes from the depot. Our C++ implementation uses SCIP 3.0.2 libraries [1] and has been obtained by adapting the code for the basic VRP provided by Prof. Dr. Andreas Bley (Institut for Mathematics Department, University of Kassel). Tests have been run on a 3.4 GHz Intel i-5-3570k CPU with 4 GB RAM, using SoPlex 1.7.2 as LP solver. Tables 2 and 3 report a summary of the results for the 25-nodes instances (see [5] for details). We consider four cases: our approach using -lookback approximation, with ∈ {1, 2, 3} (results are reported in columns “L”), and the deterministic case with no delay (column “det.”), which can also be seen as a heuristic for the stochastic case. We say that an instance is closed if all the branch-and-bound nodes are closed within 30 minutes, for all the four cases: each row of Tables 2 and 3 reports, for each group, average results over the closed instances. Table 2 shows the group name, the number of closed instances, and average running times for the four cases. Just 5 instances over 56 are not closed, and running times are of the same order as
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
331
Table 3 Solution cost and simulated cost. Group R1 C1 RC1 R2 C2 RC2
Objective function det. L1 L2 L3 464.4 487.1 489.4 492.0 191.1 205.6 219.4 226.4 351.1 362.4 369.0 377.5 383.1 385.8 385.8 390.8 215.5 215.5 215.5 215.5 331.8 331.8 331.9 332.3
det. 534.3 282.0 524.3 400.6 228.4 333.3
Simulated Cost L1 L2 L3 506.0 506.0 494.5 261.7 247.7 240.4 433.5 411.5 406.1 389.5 389.5 390.8 220.8 220.7 220.8 332.3 332.3 332.3
Accuracy det. L1 15.0 3.9 47.5 27.3 49.3 19.6 4.6 1.0 6.0 2.3 0.4 0.4
(% gap) L2 L3 3.4 0.5 12.9 6.2 11.5 7.5 1.0 0.0 2.3 2.3 0.4 0.3
the deterministic case: the ratio is about 2 on average for = 1 and about 3 on average for both = 2 and = 3, meaning also that using = 3 does not remarkably take more time than = 2. Table 3 is devoted to solution costs. Columns “Objective function” report the heuristic value of the solution obtained in the four cases. We recall that this value is different from the expected cost (due to the -lookback approximation), and it is not the optimal one (due to heuristic pricing). A better estimate of the expected cost is obtained by extensively simulating the cost of the routes provided as a solution in the four cases, and average values are reported in columns “Simulated Cost”. The average per-cent gap between heuristic and simulated values is reported in columns “Accuracy”. As the effect of a better probability approximation, the trend of the objective function value is to increase with , while the simulated value decreases. In fact, = 3 is sufficiently large to make the objective function value very close to the simulated value (in most cases the difference is less then 4%). Preliminary results show that -lookback approximation and the heuristic pricing procedure based on extended labels represent a viable approach to the S-VRPHTW, providing meaningful results with running times that are comparable to the deterministic case. Several research issues remain open. We mention here the need for further testing and analysis to validate the results and assess the value of the stochastic approach. Also, the definition of lower bounds for the problem would allow a better assessment of the quality of the provided solution, and would open the way toward exact branch-and-price procedures for the S-VRPHTW.
References [1] Achterberg, T., SCIP: Solving constraint integer programs, Mathematical Programming Computation 1 (2009), pp. 1–41.
332
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
[2] Ando, N. and E. Taniguchi, Travel time reliability in vehicle routing and scheduling with time windows, Networks and Spatial Economics 6 (2006), pp. 293–311. [3] Andreatta, G., L. Capanna, L. De Giovanni, M. Monaci and L. Righi, Efficiency and robustness in a support platform for intelligent airport ground handling, Journal of Intelligent Transportation Systems 18 (2014), pp. 121–130. [4] Andreatta, G., M. Casula, C. De Francesco and L. De Giovanni, On the evaluation of the expected number of disruptions in VRP with hard time windows, Technical report, Dip. di Matematica, Universit` a di Padova. (2014), http://www.math.unipd.it/∼luigi/manuscripts/TR DM OR 01 2014.pdf. [5] Casula, M., “A branch-and-price algorithm for the VRP with time windows and stochastic travel times,” Master’s thesis, Dipartimento di Matematica, Universit`a di Padova (2014). [6] Christiansen, H. C. and J. Lysgaard, A branch-and-price algorithm for the capacitated vehicle routing problem with stochastic demands, Operations Research Letters 37 (2007), pp. 773–781. [7] Desrochers, M., J. Desrosiers and M. Solomon, A new optimization algorithm for the vehicle routing problem with time windows, Operations Research 40 (1992), pp. 342–354. [8] Feillet, D., P. Dejax, M. Gendreau and C. Gueguen, An exact algorithm for the elementary shortest path problem with resource constraints: Application to some vehicle routing problems, Networks 44 (2004), pp. 216–229. [9] Gendreau, M., G. Laporte and R. S´eguin, Stochastic vehicle routing, European Journal of Operational Research 88 (1996), pp. 3–12. [10] Li, X., P. Tian and S. Leung, Vehicle routing problems with time windows and stochastic travel and service times: models and algorithm, International Journal of Production Economics 125 (2010), pp. 137–145. [11] Solomon, M., Algorithms for the vehicle routing and scheduling problems with time window constraints, Operations Research 35 (1987), pp. 254–265. [12] Ta¸s, D., M. Gendreau, N. Dellaert, T. van Woensel and A. de Kok, Vehicle routing with soft time windows and stochastic travel times: A column generation and branch-and-price solution approach, European Journal of Operational Research (2013). [13] Toth, P. and D. Vigo, “The vehicle routing problem,” Society for Industrial and Applied Mathematics, 2002.
Electronic Notes in Discrete Mathematics 52 (2016) 325–332 www.elsevier.com/locate/endm
A branch-and-price based heuristic for the stochastic vehicle routing problem with hard time windows G. Andreatta 1, M. Casula 2, C. De Francesco 3, L. De Giovanni 4 Dipartimento di Matematica Universit` a degli Studi di Padova Padova, Italy
Abstract In the Vehicle Routing Problem with Hard Time Windows and Stochastic travel times, a disruption occurs if, due to stochastic events, a vehicle arrives too late at a customer. In this case a recourse action is required such that the service starts within the time window, and a relevant penalty cost is incurred. Despite the problem has been inspired by a real-life application in airport ground handling optimization, it has never been addressed before in literature, to the best of our knowledge. We discuss how the expected penalty cost can be evaluated and how this computation can be integrated in a branch-and-price procedure to obtain heuristic solutions. Preliminary tests on literature instances show the effectiveness of the approach. Keywords: Vehicle routing, Hard time windows, Stochastic travel times, Branch-and-price
1 2 3 4
Email: Email: Email: Email:
[email protected] [email protected] [email protected] [email protected]
http://dx.doi.org/10.1016/j.endm.2016.03.043 1571-0653/© 2016 Elsevier B.V. All rights reserved.
326
1
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
Introduction
The Vehicle Routing Problem (VRP) is a well known network optimization problem where a set of customers is located at the nodes of a graph and a set of node-disjoint routes have to be determined in order to visit all the customers. Routes have to start and end at a same node called depot and travel distances between nodes, as well as node service times, are given. The objective is to minimize the total cost, normally related to the total distance traveled (see [13], among others). One of the possible extensions of this problem is the VRP with time windows (VRPTW), where the service of each customer has to start within a specific time interval. Time windows may be soft (they can be violated at some penalty cost), or hard (early vehicles have to wait at no additional cost, while late vehicles are prohibited to serve). Real-life applications often require that uncertainty is taken into account so that stochastic variants of the VRP have been introduced (see e.g. [9] for a review). We consider the VRP with Hard Time Windows and Stochastic travel and service times (S-VRPHTW). Due to the uncertainty, a customer originally included in a route may be reached after the end of its time window: we call disruption this event, and we assume that a recourse action is taken such that the customer is served within its time window and a large fixed penalty cost is incurred. As a consequence, the objective is to minimize a weighted sum of the total distance traveled and the expected number of disruptions. The study of the S-VRPHTW has been inspired by a real application related to the optimization of ground handling operations in airport aprons, where the customers are the aircraft that need specific operations before/after departure/landing, and the vehicles are the specific Ground Service Equipment (GSE) required by the same operations (see [3]). If a GSE cannot serve an aircraft on time (for example because the previous aircraft is late, or the delay cumulated along the service route is sufficiently large, or because of a breakdown or other unforeseen events), a recourse GSE takes its place in the originally planned route, causing an additional relevant penalty cost. Among the stochastic variants of the VRP, the ones involving stochastic travel times are less studied and literature is concerned with soft time windows: we cite some heuristic approaches based on genetic algorithms or tabu search (see [2,10] among others), and an exact algorithm based on branch-and-price [12]. Also the method proposed in [6] for the VRP with stochastic demands may be adapted to soft time windows, under suitable restrictions on delay distributions. The S-VRPHTW is concerned with hard time windows and, to the best of our knowledge, it has not been yet considered in literature.
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
327
In this paper, we propose a first heuristic approach for the S-VRPHTW. Theoretical and practical issues related to the computation of the expected number of disruptions are addressed in Section 2, proposing an approximation and discussing how it can be integrated in a branch-and-price procedure to obtain heuristic solutions. Preliminary results on literature benchmarks show that the proposed approach is able to provide good solutions, trading off between computational efficiency and accuracy of the stochastic approximations, as outlined in Section 3, together with some lines for future research.
2
The evaluation of the expected number of disruptions and the branch-and-price procedure
Let G = (N, A) be a directed graph where N is the set of customers. For each arc (i, j) ∈ A, let tij be the travel time from customer i to customer j and cij the related cost. For each node v ∈ N, let rv be the service time and [av , bv ], with av ≤ bv the time window. The depot is represented by node v0 . Let p = v0 − v1 − . . . − vn be an elementary path on G, that is, a path without node repetitions. With relation to p, the starting time si of a node vi is recursively defined as follows: s0 = 0; si = max{ai , si−1 + ri−1 + ti−1 i } for every i = 1...n. Path p is feasible if vn = v0 and ai ≤ si ≤ bi , for every i = 0...n. We assume that travel and service times are subject to independent random delays δi j , for each (i, j) ∈ A. If P is the penalty cost associated to a disruption, the cost of p is defined as c(p) =
n−1
cvi vi+1 + P · E(μp )
(1)
i=0
where μp is a random variable representing the number of disruptions on p. We consider the following set-partitioning formulation for S-VRPHTW [7]: c(p)xp (2) min p∈Ω
s.t.
xp = 1
∀v∈N
(3)
∀p∈Ω
(4)
p∈Ω(v)
xp ∈ {0, 1}
where Ω is the set of all feasible paths, Ω(v) is the set of feasible paths containing v, and xp takes value one if p is included in the solution, 0 otherwise. For each path p we have to evaluate c(p), and hence E(μp ). If d(i) is the event
328
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
“a disruption occurs at the node vi of p”, then E(μp ) = ni=1 P [d(i)]. Clearly, P [d(i)] depends on the delay cumulated up to node vi . We notice that, because of the hard time windows, the recourse action is intended to absorb the effect of delays in case of disruption (for example, a recourse vehicle replaces the delayed one on the originally planned route). Let σi be the random variable denoting the random starting time for node vi . We have σ0 = 0 and, for i = 1...n mi = max{ai , σi−1 + ri−1 + ti−1 i + δi−1 i } if mi ≤ bi σi = (5) ui otherwise where ui is any value within the time window. We assume that the recourse action allows the service of customer i to start at the originally planned starting time, that is ui = si . By denoting with e(k, i) (for 0 ≤ k < i) the event “the first disruption after node vk occurs in vi (independently from the occurrence of other disruptions on the route from v0 to vk )”, we have P [d(i)] = P [e(0, i)] +
i−1
P [e(k, i)|d(k)] · P [d(k)]
(6)
k=1
Under the hypothesis ui = si , P [e(k, i)|d(k)] = P[(δi−1 i > bi − (ti−1 i + ri−1 + σi−1 ))∧ ∧ (δi−2 i−1 ≤ bi−1 − (ti−2 i−1 + ri−2 + σi−2 ))∧ ∧ . . . ∧ (δk k+1 ≤ bk+1 − (tk k+1 + rk + sk ))]
(7)
We notice that P [e(0, i)] can be evaluated as in (7) with k = 0 and that P [d(1)] = P [δ0 1 > b1 − t0 1 ], hence we can recursively obtain P [d(i)]. Formula (7) can be analytically integrated, at least in theory, starting from the p.d.f. of variables δ and considering the sequence of the realizations of the same variables from δk k+1 to δi−1 i (see [4] for details). In practice, the integration requires to consider a number of cases exponential in the length of the path. Thus we resort to an -lookback approximation to evaluate P [d(i)], which consists in considering only the last delays before vi , that is, we reformulate equation (6) assuming that the first node of the path is vk¯ , k¯ = max{0, i − } and σk¯ = sk¯ . In this way we have only a fixed (possibly small) number of integration cases, regardless of the actual length of the path. The idea behind this approximation is that, since the delay resets whenever the vehicle arrives outside a time window, considering the whole delay history may
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
329
be redundant. Computational experiments in some particular cases show that a good approximation can be obtained even with small values of . Formulation (2)...(4) contains an exponential number of variables and we adopt a branch-and-price approach, that is, a branch-and-bound procedure using column generation to solve the linear relaxation: we initially consider a small subset of feasible paths, and a pricing subproblem is iteratively solved to determine further feasible paths that should be included in the formulation to get better solutions. The components we implemented for solving S-VRPHTW are essentially the ones of a standard branch-and-price for VRPTW (see e.g. [7]). The pricing problem asks for a feasible path with minimum reduced cost c¯(p) = c(p) − vi ∈p π ¯i , where π ¯i are the current dual variables associated with (3). This can be seen as an Elementary Shortest Path Problem with Resource Constraint (ESPPRC), where resources are time (to guarantee that arrival times fall into the time windows) and vehicle capacity, where defined. We follow the approach of [8] where ESPPRC is solved by a labeling algorithm. We use labels L(p) to store the features of partial feasible paths p from v0 to vi . Each label stores: c¯(p); si ; the total capacity used by p; a binary value Aj for each j = 1 . . . n (Aj is equal to 1 if node vj cannot be visited in any extension of p, because it has been already visited, or due to resource constraints); and D(p), the expected delay at the end of p. Notice that labels are referred to paths. In fact, when we extend a path p to p using arc (i, j), we have to compute c¯(p ) and D(p ) which depend on the node sequence (and not only on (i, j) and cumulated information up to vi , as e.g. in [6]). Of course, storing the labels of every partial path is not efficient, and, at each node, we do not store dominated labels, that is, labels having, for each field, a non-smaller value than an already stored one. This dominance criterion would be exact if we were able to store the actual delay cumulated on partial paths (the same extension of a path with greater delay and greater reduced cost would lead to a path with a greater reduced cost, since the expected number of disruptions would be greater). However, we cannot store the actual delays, so that storing their expected values ends up with a heuristic procedure for pricing (and hence for S-VRPHTW).
3
Computational results and conclusions
Preliminary tests have been conducted to determine the accuracy of the lookback approximation in estimating E(μp ). To this aim, we consider “regular” paths p of length n having, for any node vi , ri + ti i+1 = 2, ai = 2 i, bi = 2 i + 1 and exponentially distributed delay with average D. In this case,
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
330
Table 1 Expected number of disruptions on a regular path with 30 nodes. D 0.1 0.2 0.3 0.7
=1 0.01 0.98 3.52 11.83
=2 0.06 2.41 5.94 12.47
=3 0.20 3.98 7.08 12.22
=4 0.50 5.03 7.07 12.16
=5 0.98 5.31 6.73 12.18
=6 1.58 5.07 6.55 12.18
=7 2.20 4.72 6.56 12.18
=8 2.70 4.53 6.61 12.18
E(μp ) 2.31 4.62 6.60 12.18
Table 2 Average running times (in seconds). Group R1 C1 RC1 R2 C2 RC2
closed 12 over 12 9 over 9 8 over 8 10 over 11 6 over 8 6 over 8
det. 1.5 84.9 2.0 151.6 15.6 39.7
L1 1.6 92.4 2.2 331.7 33.6 96.0
L2 13.1 165.9 17.3 469.1 71.3 137.6
L3 8.5 169.2 12.6 399.2 73.2 148.5
we have analytically evaluated E(μp ), both exactly and with -lookback approximation for different values of (see [4]). Sample results with n = 30 are reported in Table 1. In general, with increasing , the approximation becomes better and better and the convergence to the exact value is faster for larger D. Notice that even small values of give good approximations. The branch-and-price based heuristic has been tested on instances from the Solomon’s benchmark [11]. We consider independent exponentially distributed delays associated to each arc (i, j), with average 0.1 · (ri + tij ), and we set P equal to twice the average distance of nodes from the depot. Our C++ implementation uses SCIP 3.0.2 libraries [1] and has been obtained by adapting the code for the basic VRP provided by Prof. Dr. Andreas Bley (Institut for Mathematics Department, University of Kassel). Tests have been run on a 3.4 GHz Intel i-5-3570k CPU with 4 GB RAM, using SoPlex 1.7.2 as LP solver. Tables 2 and 3 report a summary of the results for the 25-nodes instances (see [5] for details). We consider four cases: our approach using -lookback approximation, with ∈ {1, 2, 3} (results are reported in columns “L”), and the deterministic case with no delay (column “det.”), which can also be seen as a heuristic for the stochastic case. We say that an instance is closed if all the branch-and-bound nodes are closed within 30 minutes, for all the four cases: each row of Tables 2 and 3 reports, for each group, average results over the closed instances. Table 2 shows the group name, the number of closed instances, and average running times for the four cases. Just 5 instances over 56 are not closed, and running times are of the same order as
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
331
Table 3 Solution cost and simulated cost. Group R1 C1 RC1 R2 C2 RC2
Objective function det. L1 L2 L3 464.4 487.1 489.4 492.0 191.1 205.6 219.4 226.4 351.1 362.4 369.0 377.5 383.1 385.8 385.8 390.8 215.5 215.5 215.5 215.5 331.8 331.8 331.9 332.3
det. 534.3 282.0 524.3 400.6 228.4 333.3
Simulated Cost L1 L2 L3 506.0 506.0 494.5 261.7 247.7 240.4 433.5 411.5 406.1 389.5 389.5 390.8 220.8 220.7 220.8 332.3 332.3 332.3
Accuracy det. L1 15.0 3.9 47.5 27.3 49.3 19.6 4.6 1.0 6.0 2.3 0.4 0.4
(% gap) L2 L3 3.4 0.5 12.9 6.2 11.5 7.5 1.0 0.0 2.3 2.3 0.4 0.3
the deterministic case: the ratio is about 2 on average for = 1 and about 3 on average for both = 2 and = 3, meaning also that using = 3 does not remarkably take more time than = 2. Table 3 is devoted to solution costs. Columns “Objective function” report the heuristic value of the solution obtained in the four cases. We recall that this value is different from the expected cost (due to the -lookback approximation), and it is not the optimal one (due to heuristic pricing). A better estimate of the expected cost is obtained by extensively simulating the cost of the routes provided as a solution in the four cases, and average values are reported in columns “Simulated Cost”. The average per-cent gap between heuristic and simulated values is reported in columns “Accuracy”. As the effect of a better probability approximation, the trend of the objective function value is to increase with , while the simulated value decreases. In fact, = 3 is sufficiently large to make the objective function value very close to the simulated value (in most cases the difference is less then 4%). Preliminary results show that -lookback approximation and the heuristic pricing procedure based on extended labels represent a viable approach to the S-VRPHTW, providing meaningful results with running times that are comparable to the deterministic case. Several research issues remain open. We mention here the need for further testing and analysis to validate the results and assess the value of the stochastic approach. Also, the definition of lower bounds for the problem would allow a better assessment of the quality of the provided solution, and would open the way toward exact branch-and-price procedures for the S-VRPHTW.
References [1] Achterberg, T., SCIP: Solving constraint integer programs, Mathematical Programming Computation 1 (2009), pp. 1–41.
332
G. Andreatta et al. / Electronic Notes in Discrete Mathematics 52 (2016) 325–332
[2] Ando, N. and E. Taniguchi, Travel time reliability in vehicle routing and scheduling with time windows, Networks and Spatial Economics 6 (2006), pp. 293–311. [3] Andreatta, G., L. Capanna, L. De Giovanni, M. Monaci and L. Righi, Efficiency and robustness in a support platform for intelligent airport ground handling, Journal of Intelligent Transportation Systems 18 (2014), pp. 121–130. [4] Andreatta, G., M. Casula, C. De Francesco and L. De Giovanni, On the evaluation of the expected number of disruptions in VRP with hard time windows, Technical report, Dip. di Matematica, Universit` a di Padova. (2014), http://www.math.unipd.it/∼luigi/manuscripts/TR DM OR 01 2014.pdf. [5] Casula, M., “A branch-and-price algorithm for the VRP with time windows and stochastic travel times,” Master’s thesis, Dipartimento di Matematica, Universit`a di Padova (2014). [6] Christiansen, H. C. and J. Lysgaard, A branch-and-price algorithm for the capacitated vehicle routing problem with stochastic demands, Operations Research Letters 37 (2007), pp. 773–781. [7] Desrochers, M., J. Desrosiers and M. Solomon, A new optimization algorithm for the vehicle routing problem with time windows, Operations Research 40 (1992), pp. 342–354. [8] Feillet, D., P. Dejax, M. Gendreau and C. Gueguen, An exact algorithm for the elementary shortest path problem with resource constraints: Application to some vehicle routing problems, Networks 44 (2004), pp. 216–229. [9] Gendreau, M., G. Laporte and R. S´eguin, Stochastic vehicle routing, European Journal of Operational Research 88 (1996), pp. 3–12. [10] Li, X., P. Tian and S. Leung, Vehicle routing problems with time windows and stochastic travel and service times: models and algorithm, International Journal of Production Economics 125 (2010), pp. 137–145. [11] Solomon, M., Algorithms for the vehicle routing and scheduling problems with time window constraints, Operations Research 35 (1987), pp. 254–265. [12] Ta¸s, D., M. Gendreau, N. Dellaert, T. van Woensel and A. de Kok, Vehicle routing with soft time windows and stochastic travel times: A column generation and branch-and-price solution approach, European Journal of Operational Research (2013). [13] Toth, P. and D. Vigo, “The vehicle routing problem,” Society for Industrial and Applied Mathematics, 2002.