# The Hybrid Vehicle Routing Problem

## The Hybrid Vehicle Routing Problem

Transportation Research Part C 78 (2017) 1–12 Contents lists available at ScienceDirect Transportation Research Part C journal homepage: www.elsevie...

Transportation Research Part C 78 (2017) 1–12

Contents lists available at ScienceDirect

Transportation Research Part C journal homepage: www.elsevier.com/locate/trc

The Hybrid Vehicle Routing Problem Simona Mancini Politecnico di Torino, Corso Duca Degli Abruzzi 24, 10129 Torino, Italy

a r t i c l e

i n f o

Article history: Received 10 January 2016 Received in revised form 7 November 2016 Accepted 4 February 2017 Available online 7 March 2017 Keywords: Routing Hybrid vehicles Refueling station Matheuristic Mixed Integer Linear Programming

a b s t r a c t In this paper the Hybrid Vehicle Routing Problem (HVRP) is introduced and formalized. This problem is an extension of the classical VRP in which vehicles can work both electrically and with traditional fuel. The vehicle may change propulsion mode at any point of time. The unitary travel cost is much lower for distances covered in the electric mode. An electric battery has a limited capacity and may be recharged at a recharging station (RS). A limited number of RS are available. Once a battery has been completely discharged, the vehicle automatically shifts to traditional fuel propulsion mode. Furthermore, a maximum route duration is imposed according to contracts regulations established with the driver. In this paper, a Mixed Integer Linear Programming formulation is presented and a Large Neighborhood Search based Matheuristic is proposed. The algorithm starts from a feasible solution and consists into destroying, at each iteration, a small number of routes, letting unvaried the other ones, and reconstructing a new feasible solution running the model on only the subset of customers involved in the destroyed routes. This procedure allows to completely explore a large neighborhood within very short computational time. Computational tests that show the performance of the matheuristic are presented. The method has also been tested on a simplified version of the HVRP already presented in the literature, the Green Vehicle Routing Problem (GVRP), and competitive results have been obtained. Ó 2017 Elsevier Ltd. All rights reserved.

2

S. Mancini / Transportation Research Part C 78 (2017) 1–12

S. Mancini / Transportation Research Part C 78 (2017) 1–12

3

I: Set of customers I0 : Set of customers and the depot I [ v 0 F: Set of refueling stations V: Set of nodes without the depot I [ F V 0 : Set of nodes I0 [ F r: Battery charge consumption rate (for km) (for electric propulsion) Q: Battery capacity q: Distance unitary penalty for using traditional fuel propulsion M: Maximum number of routes lj : Minimum distance from j to the nearest refueling stations or to the depot

4

S. Mancini / Transportation Research Part C 78 (2017) 1–12

Fig. 1. Optimal routes for VRP (a),GVRP (b) and HVRP (c).

The variables of the problem are now described.     

X ij : Binary variable equal to 1 if a vehicle travels from node i to node j and 0 otherwise Y j : Battery level upon arrival at node j T j : Arrival time at node j W j : Distance covered with traditional fuel to reach customer j U j : Distance covered with traditional fuel to reach the depot from customer j (defined only for j in I) The mathematical formulation of the HVRP is the following:

min

XX X dij X ij þ qðW j þ U j Þ i2V 0 j2V 0

j2V 0

ð1Þ

S. Mancini / Transportation Research Part C 78 (2017) 1–12

s.t.

X X ij ¼ 1 8 i 2 V

5

ð2Þ

j2V j–i

X X ij 6 1 8 i 2 F

ð3Þ

j2V 0 j–i

X X X ji ¼ X ij i2V 0 i–j

8 j 2 V0

ð4Þ

i2V 0 i–j

X X 0j 6 M

ð5Þ

j2V 0 j–0

X X j0 6 M

ð6Þ

j2V 0 j–0

T j P T i þ ðtij þ pj ÞX ij  T max ð1  X ij Þ 8 i 2 V 0 j 2 V and i – j

ð7Þ

0 6 T 0 6 T max

ð8Þ

t0j 6 T j 6 T max  ðt j0 þ pij Þ 8 j 2 V

ð9Þ

Y j 6 Y i  rðdij  W j Þ þ Q ð1  X ij Þ 8j 2 I i 2 V 0 and i – j Yj ¼ Q

8j 2 F [ v 0

ð10Þ ð11Þ

Y j P rðdj0 X j0  U j Þ 8j 2 I

ð12Þ

X ij 2 0; 1 8i 2 V 0 8j 2 V 0

ð13Þ

Y j P 0 8i 2 I

ð14Þ

The objective function is defined in (1). Constraint (2) implies that each customer is visited exactly once, while each refueling station may be visited only once, as stated in (3). Route continuity is guaranteed by constraint (4) while a maximum number of routes is imposed by constraints (5) and (6). The arrival time at each node is tracked by constraint (7). Constraints (8) and (9) make certain that each vehicle returns to the depot no later than a timelimit, T max . The battery charging level, upon arrival at each node, is ruled by constraint (10). The time and battery charge level tracking constraints, constraints (7) and (10), respectively, exclude the possibility of subtours. Constraint (11) set the value of battery charge equal to Q at the departure from the depot and resets it up to Q upon a visit to a recharging station. Constraint (12) implies that if the remaining battery charge is not enough to reach the depot, or a recharging station, the remaining distance is covered at a higher cost, by means of traditional fuel propulsion. Since the autonomy of traditional fuel tank is much larger than that of electric batteries, and since fuel stations are located frequently along the road network, no capacity restriction is considered for fuel tanks. Finally, the domain of the variables is specified by constraints (13) and (14). 5. A large neighborhood search based matheuristics for the HVRP Large Neighborhood search heuristics, (LNS), belong to the class of heuristics known as Very Large Scale Neighborhood search (VLSN) algorithms, as stated in Ahuja et al. (2002). All VLSN algorithms are based on the observation that searching a large neighborhood results in finding local optima of high quality, and hence a VLSN algorithm may return better solutions overall. However, searching a large neighborhood is very time consuming, hence various filtering techniques are used to limit the search. In VLSN algorithms, the search is usually restricted to a subset of the solutions belonging to the neighborhood which can be efficiently explored. Unlike what happens in other VLSN, the neighborhood is implicitly defined, in LNS, by the moves used to destroy and repair an incumbent solution. For a detailed survey on LNS applications to routing problems the reader may refer to Pisinger and Ropke (2010). The destroy operators may be defined in different ways. For routing problems, for instance, a destroy operator could consist in breaking k routes leaving the others unvaried, or in

6

S. Mancini / Transportation Research Part C 78 (2017) 1–12

removing a fixed percentage of the arcs in the current solution. A random (or randomized) component is used to select the arcs that have to be removed. The repair method rebuilds a feasible solution starting from the partially destroyed one. Generally, a greedy construction heuristic is used to rebuild the solution. This is a very fast but not always very accurate method, since only a sample solution is analyzed in the neighborhood. The innovative aspect of the LNS proposed in this paper concerns the possibility of addressing the whole neighborhood in reasonably short computational time. In fact, the large neighborhood search is exploited directly by the model. In this way, it is possible to obtain the local minimum with respect to the addressed neighborhood, which renders the intensification phase of the algorithm more powerful and precise. A detailed explanation of how the proposed matheuristic works is given in the following subsection. 5.1. Algorithm description The proposed Matheuristic, from now on called MH, works as follows. The procedure starts from a feasible solution. Two routes are destroyed at each iteration and the model is run again, with a short timelimit (i.e. 10 s) only on those customers who/that were previously involved in those routes, while the other routes are left unvaried. As the number of customer involved is small, the model is able to solve the problem to optimality (or almost to optimality) in a very short time. The procedure terminates when no further improvement can be found from destroying only two routes, or when a maximum number of iterations is reached. The idea is to apply a classical destroy operator while exploiting the model to efficiently explore a large neighborhood of solutions instead of reconstructing a new feasible solution with a greedy procedure, as it is normal practice in LNS algorithms. In this way, the reconstructed solution quality results to be much higher, because there is the possibility of exhaustively (or almost) exploring the neighborhood within a few seconds and of easily reaching the local minimum in the neighborhood, while in classical LNS algorithms, only a small part of the neighborhood is explored. The destroy operator implies a strong perturbation on the solution to ensure diversification in the search process. Let us introduce additional notations:  R ¼ r 1 ; r 2 ; . . . ; r n : set of routes in the current solution  Rbest ¼ r1 ; r2 ; . . . ; rn : set of routes in the best solution A pseudocode of the algorithm is reported hereafter: Algorithm 1. A Large neighborhood search matheuristics for the HVRP. set iter = 0 set Sb ¼ S0 for a ¼ 0; a 6 n; a þ þ do for b ¼ 0; b 6 n; b þ þ do if iter 6 MAXITER then set R ¼ Rbest destroy routes ra and rb run the model with 2 vehicles and considering as customers only customers that were involved in ra and rb obtain two new routes r 0a and r 0b if the sum of the cost of r 0a and r 0b 6 the sum of the cost of r a and r b then substitute ra and rb with r 0a and r 0b in R set Rbest ¼ R set a ¼ 0 and b ¼ 0 end if set iter = iter + 1 end if end for end for

6. Computational results In order to prove the efficiency of the proposed algorithm, it has been tested on a special case of the HVRP, the Green Vehicle Routing Problem (GVRP), for which benchmark results may be found in the literature. Computational tests have been carried out on a machine with a processor Intel i7 at 1.9 GHz with 8 GB of Ram. All the procedure has been coded in Xpressive and the model is solved with Xpress 7.8. The maximum number of iterations, MAXITER, has been fixed equal to 100 and the timelimit for each iteration of the LNS equal to 10 s.

7

S. Mancini / Transportation Research Part C 78 (2017) 1–12

6.1. A particular case: the GVRP The GVRP can be considered a special case of the HVRP in which the use of traditional fuel propulsion is not allowed, i.e. the fleet is composed of pure electric propulsion vehicles. This further constraints can be easily added to the HVRP model by forcing all the W j variables to zero. This problem has been introduced in Erdogan and Miller-Hooks (2012), where two constructive heuristics were proposed. The first one is a Modified Clarke and Wright Savings heuristic, while the second one is a Density-Based Clustering Algorithm. In the following, the heuristics are referred to as Erdogan1 and Erdogan2, respectively. In Schneider et al. (2014), the authors have presented a Variable Neighborhood search heuristic combined with a Tabu Search while a non-deterministic Simulated Annealing has been proposed in Felipe et al. (2014). Although these heuristics were originally proposed for two different extensions of the GVRP, they were also tested on the GVRP instances introduced by Erdogan and Miller-Hooks (2012); therefore it is possible to compare the MH results with those obtained by these heuristics. Computational tests have been carried out on two different set of instances: the first one contains small size instances (20 customers), while the second one is composed by large size instances with a range of customers varying from 100 to 500. 6.1.1. Small size instances The small size instances proposed in Erdogan and Miller-Hooks (2012) are divided in 4 sets:  S1: uniform customer distribution - 10 randomly generated instances of 20 uniformly distributed customers with 3 recharging station locations  S2: clustered customer distribution - 10 randomly generated instances of 20 clustered customers with 3 recharging station locations  S3: Impact of the spatial recharging station configuration - 10 instances, half selected from S1 and half from S2, each instance with 6 recharging stations randomly generated  S4: Impact of station density - 10 instances, half of which have been created from one instance of S1 and half from one instance of S2, while gradually increasing the number of recharging station from 2 to 10 in increments of 2 The initial solution used as the starting point for the MH is the best solution that can be obtained when the model is run with a timelimit of 5 s. Details of the results obtained with MH on the above described sets are reported in Tables 1–4. The tables are organized as follows. The first column report the name of the instance while columns 2 and 3 report the objective function obtained with Erdogan1 and Erdogan2. Columns 4 and 5 report the objective function and computational time (expressed in seconds) for Schneider while the same data, for Felipe and for the MH proposed in this paper, are reported in columns 6–7 and 8–9 respectively. The second last row reports the average values, while the percentage gap between MH and the other heuristics is reported in the last row. The gap is computed as ðSh  SMH Þ=SMH , where Sh represents the value of the objective function of the solution obtained with the hth heuristic and SMH the value of the one obtained by MH, where h varies in ½1; 4. Erdogan and Miller-Hooks (2012), do not provide explicit data about the computational time for each instances but they declare that they are of the order of seconds. The MH proposed in this paper has been executed on an Intel Core i7, 2.0 GHz, 4 GB RAM while Schneider on a Intel Core i5, 2.67 GHz, 4 GB RAM, Felipe on a Intel Core i5, 2.8 GHz, 8 GB RAM awhile both Erdogan1 and Erdogan2 were run on a Pentium4, 3.20 GHz, 2 GB RAM, but no computational times are reported by the authors. As shown in the resume reported in Table 5, MH clearly outperforms Erdogan1 and Erdogan2, obtaining better solutions for most instances and providing an average saving of around 12%, while it obtains comparable results with the best

Table 1 GVRP: Comparison of results on set S1. INSTANCE

Erdogan1

Erdogan2

Schneider

TIME

Felipe

TIME

MH

TIME

20c3sU1 20c3sU2 20c3sU3 20c3sU4 20c3sU5 20c3sU6 20c3sU7 20c3sU8 20c3sU9 20c3sU10

1818.35 1614.15 1969.64 1508.41 1752.73 1668.16 1730.45 1718.67 1714.43 1309.52

1797.51 1613.53 1964.57 1487.15 1752.73 1668.16 1730.45 1718.67 1714.43 1309.52

1797.49 1574.77 1704.48 1482 1689.37 1618.65 1713.66 1706.5 1708.81 1181.31

41 38 38 39 40 40 38 40 40 38

1805.41 1574.78 1704.48 1482 1689.37 1618.65 1713.67 1722.78 1708.82 1181.31

1.64 1.48 1.47 1.53 2.15 1.9 2.09 1.76 2.14 1.44

1797.49 1574.77 1704.48 1482 1689.37 1618.65 1713.9 1738.04 1708.81 1181.31

11 30 20 4 23 3 20 13 14 0.2

39.2

1620.127

1.76

1620.882

13.82

AVG

1680.451

1675.672

1617.704

3.68%

3.38%

0.20%

0.05%

8

S. Mancini / Transportation Research Part C 78 (2017) 1–12

Table 2 GVRP: Comparison of results on set S2. INSTANCE

Erdogan1

Erdogan2

Schneider

TIME

Felipe

TIME

MH

TIME

20c3sC1 20c3sC2 20c3sC3 20c3sC4 20c3sC5 20c3sC6 20c3sC7 20c3sC8 20c3sC9 20c3sC10

1300.62 1553.53 1083.12 1135.9 2190.68 2883.71 1701.4 3319.74 1811.05 2648.84

1300.62 1553.53 1083.12 1091.78 2190.68 2883.71 1701.4 3319.74 1811.05 2644.11

1173.57 1539.97 880.2 1059.35 2156.01 2758.17 1393.99 3139.72 1799.94 2583.42

37 35 15 32 36 43 11 37 36 27

1178.97 1539.97 880.2 1059.35 2156.01 2758.17 1393.99 3139.72 1799.94 2583.42

1.5 1.3 0.58 1.11 1.42 1.17 0.28 1.28 1.2 1.11

1178.97 1567.82 898.32 1109.73 2156.01 2758.17 1393.99 3139.72 1799.94 2583.42

20 21 6 13 13 30 1 55 24 11

AVG

1962.859

1957.974

1848.43

30.84

1848.97

1.10

1858.609

19.4

5.61%

5.35%

0.55%

0.52%

Table 3 GVRP: Comparison of results on set S3. INSTANCE

Erdogan1

Erdogan2

Schneider

TIME

Felipe

TIME

MH

TIME

S1_2i6s S1_4i6s S1_6i6s S1_8i6s S1_10i6s S2_2i6s S2_4i6s S2_6i6s S2_8i6s S2_10i6s

1614.15 1561.3 1616.2 1902.51 1309.52 1645.8 1505.06 3115.1 2722.55 1995.62

1614.15 1541.46 1616.2 1882.54 1309.52 1645.8 1505.06 3115.1 2722.55 1995.62

1578.12 1397.27 1560.49 1692.32 1173.48 1633.1 1532.96 2431.33 2158.35 1958.46

43 45 44 44 43 45 53 47 34 37

1578.12 1413.97 1571.3 1692.33 1173.48 1645.8 1505.07 2660.49 2175.66 1585.46

1.53 1.64 1.73 1.73 1.5 1.64 1.34 2.42 1.15 1.11

1578.12 1397.27 1597.31 1692.32 1173.48 1645.8 1544.98 2431.33 2175.66 1635.78

30 9 15 14 3 25 5 25 16 9

AVG

1898.781

1894.8

1711.588

43.38

1700.17

1.58

1687.205

15.1

12.54%

12.30%

1.45%

0.77%

Table 4 GVRP: Comparison of results on set S4. INSTANCE

Erdogan1

Erdogan2

Schneider

TIME

Felipe

TIME

MH

TIME

S1_4i2s S1_4i4s S1_4i6s S1_4i8s S1_4i10s S2_4i2s S2_4i4s S2_4i6s S2_4i8s S2_4i10s

1582.2 1580.52 1561.29 1561.29 1536.04 1135.89 1522.72 1786.21 1786.21 1783.63

1582.2 1580.52 1541.46 1561.29 1529.73 1117.32 1522.72 1730.47 1786.21 1729.51

1582.21 1460.09 1397.27 1397.27 1396.02 1059.35 1446.08 1434.14 1434.14 1434.13

38 41 45 49 51 31 36 41 45 47

1598.91 1483.19 1413.97 1397.27 1396.02 1059.35 1446.08 1434.14 1434.14 1434.13

1.56 1.58 1.64 1.62 1.67 1.12 1.31 1.37 1.39 1.4

1582.21 1460.09 1397.27 1397.27 1396.02 1059.35 1446.08 1434.14 1434.14 1434.13

14 40 9 19 50 6 8 5 25 40

42.36

1409.72

1.47

1404.07

21.6

AVG

1583.6

1568.143

1404.07

12.79%

11.69%

0.00%

0.40%

Table 5 Resume of results on the GVRP. SET

Erdogan1

Erdogan2

Schneider

Felipe

MH

S1 S2 S3 S4

1680.45 1962.86 1898.78 1583.60

1675.67 1957.97 1894.80 1568.14

1617.70 1848.43 1711.59 1404.07

1620.13 1848.97 1700.17 1409.72

1620.88 1858.61 1687.21 1404.07

12.79%

11.69%

0.00%

0.40%

9

S. Mancini / Transportation Research Part C 78 (2017) 1–12

heuristics in the literature (exactly the same results on average respect to Schneider and slightly better on average, 0.40% respect to Felipe) in terms of solutions quality, even if Felipe is much faster. 6.1.2. Large size instances The large size instances set proposed in Erdogan and Miller-Hooks (2012) is composed by 12 instances with a number of customers varying from 100 to 500. In order to address such larger problems, some changes to the algorithm have been adopted. The set of customers is divided in smaller clusters (in the range of 20–40 customer each). A feasible initial solution is obtained running the model, with a time limit of 5 s, on each cluster separately. Since, in such large instances, the number of routes created may be significantly large, the number of pairs of routes to be analyzed by the matheuristic may become huge, and most of the iterations may results useless. In fact, if we try to recombine customers belonging to two routes which are very far from each others, the probability to obtain an improvement is very limited. To overcome this issue, a proximity index, dra rb , is defined for each pair of routes ra and rb . For each pair of customers, i and j, a parameter aij is set equal to 1 if j is one of the h nearest customers to i and 0 otherwise. The proximity index dra rb is set equal to 1 if there is at least one pair of customers, i and j, for which aij = 1 and one of them is belonging to r a and the other one to rb , otherwise it is set equal to 0. Only pairs of routes having dra rb = 1 are considered for recombination by the algorithm. In this way only promising combinations are analyzed. The value of the parameter h has been set equal to 10. Computational results are reported in Table 6, which is organized as follows. For Erdogan are reported the objective function and the number of vehicles (computational times were not provided by the authors), while for Schneider, the two algorithms proposed by Felipe, and the MH proposed in this paper, are reported objective function, number of vehicles and computational times (expressed in minutes). The last two rows reports averaged values and percentage of improvement reached by the MH. As shown in the table, the MH strongly outperform Erdogan and obtain slightly better results (1.11% and 1.24%) than Felipe, while 1.75% worse than Schneider, in comparable computational time. To make a fairer comparison, it is important to remark that while MH and Felipe are performed only ones, Schneider reports the best results over 10 runs, therefore its computational times should be increased of a factor of 10. For this reason we can conclude that MH obtaines comparable results with the state of the art both in terms of efficiency and effectiveness. 6.2. The general case: The HVRP Since the HVRP has been introduced for the first time in this paper, no benchmark instances are available from the literature. Therefore, a new set of instances has been generated. This set, named SH, is composed of instances of 30, 50 and 75 customers. For each size of customers 3 instances have been constructed, with different spatial location of the refueling stations:  RS located within the customers area  RS located around the customers area  RS located far from the customers area Table 7 reports a resume of the characteristics of each instance. The battery capacity Q is fixed equal to 1000 for all the instances and the kilometric consumption rate is fixed equal to 10, which means that each vehicle has a maximum autonomy of 100 km. T max is expressed in minutes. The service time is equal to 10 min for the customers and to 60 min for the refueling

Table 6 Large size instances results. INSTANCE

111c_21 111c_22 111c_24 111c_26 111c_28 200c 250c 300c 350c 400c 450c 500c

Erodgan

Schneider

Felipe (48A)

Felipe(SA)

MH

OF

#v

OF

#v

TIME (min)

OF

#v

TIME (min)

OF

#v

TIME (min)

OF

#v

TIME (min)

5626.64 5610.57 5412.48 5408.38 5331.93 10413.59 11886.61 14229.92 16460.3 19099.04 21854.17 24517.08

20 20 20 20 20 35 41 49 57 67 75 84

4797.15 4802.16 4786.96 4778.62 4799.15 8963.46 10800.18 12594.77 14323.02 16850.21 18521.23 21170.9

17 17 17 17 17 35 39 46 51 61 68 76

21.76 23.56 21.9 25.12 24.17 76.65 120.9 182.23 232.03 305.12 525.52 356.01

4960.6 4914.2 4952.9 4934.11 4971.93 9276.63 11007.98 12869.17 14954.83 17351.92 19215.38 21636.59

18 18 18 18 18 32 39 46 54 61 68 76

21.74 21.23 21.27 21.25 21.29 76.41 120.67 184.27 227.3 302.03 514.68 352.16

5062.06 5029.17 5010.59 5092.06 5038.84 9206.28 10885.71 12827.35 14828.63 17327.27 19085.91 21475.71

18 18 18 18 18 32 38 45 52 60 67 75

12.35 12.3 12.13 12.07 11.99 73.84 108.1 302.04 332.38 384.68 456.26 154.45

4960.60 4914.20 5010.59 4778.62 4799.15 9473.27 10800.18 12869.17 14954.83 16850.21 19215.38 21989.00

18 18 18 17 17 33 39 46 54 61 68 78

20.32 21.33 22.02 20.19 22.44 79.94 108.92 208.02 280.08 355.19 566.06 304.11

156.05

10884.60

AVG

12154.226

10598.984

% IMP

100.00%

100.00%

159.58

10920.52 100.00%

157.03

10905.80 100.00%

167.39

10

S. Mancini / Transportation Research Part C 78 (2017) 1–12

Table 7 HVRP instances layout. NAME

RS LOCATION

CUSTOMERS

VEHICLES

RS

Tmax

HVRP30-1 HVRP30-2 HVRP30-3 HVRP50-1 HVRP50-2 HVRP50-3 HVRP75-1 HVRP75-2 HVRP75-3

WITHIN AROUND FAR WITHIN AROUND FAR WITHIN AROUND FAR

30 30 30 50 50 50 75 75 75

5 5 5 8 8 10 10 15 20

5 5 5 8 8 10 9 9 20

270 270 540 270 270 540 540 540 540

Table 8 UB and LB obtained by the model with a timelimit of 10,800 s. NAME

CUSTOMERS

UB

LB

GAP (%)

HVRP30-1 HVRP30-2 HVRP30-3 HVRP50-1 HVRP50-2 HVRP50-3 HVRP75-1 HVRP75-2 HVRP75-3

30 30 30 50 50 50 75 75 75

669.88 899.52 1628.84 504.02 731.56 374.53 141.1 327.6 3128.46

127.18 192.45 441.82 147.52 230.65 177.01 120.22 240.52 602.5

81.01 78.61 72.88 70.73 68.47 52.74 14.80 26.58 80.74

Table 9 UB and LB obtained by the matheuristic MH. NAME

MODEL

INIT SOL

MH

HVRP30-1 HVRP30-2 HVRP30-3 HVRP50-1 HVRP50-2 HVRP50-3 HVRP75-1 HVRP75-2 HVRP75-3

669.88 899.52 1628.84 504.02 731.56 374.53 141.1 327.6 3128.46

970.11 899.52 1940.13 1107.96 769.07 2612.23 283.63 408.32 3128.46

610.2 666.52 616.76 266.37 515.16 374.53 143.52 293.46 2221.21

AVG

933.95

1346.60

634.19

IMPROVEMENT

32.10%

52.90%

stations. The initial solution timelimit is fixed equal to 10, 30 and 60 s for instances with 30, 50 and 75 customers, respectively. The cost/km for traditional fuel covered distances is 3 time greater than the cost/km for electric propulsion. Results obtained by solving the model with the commercial solver Xpress 7.8, with a timelimit of 10,800 s, are reported in Table 8, while a comparison between the model and the matheuristic is reported in Table 9. As shown in the table, MH strongly outperforms the model results (32.1% of improvement) even starting from a very poor quality initial solution (52.9% of improvement with respect to the initial solution) which is clear a proof of the robustness of the method. The average computational time is around 100 s for the 30 customers instances, 200 s for the 50 customers ones and rises up to 1000 s on the 75 customers ones. In order to make a fair comparison, we run the model with a shorter time-limit comparable with the average computational time of the MH, and compare the best solution obtained by the model with those obtained by the MH. The results obtained by the model within the shorter time-limit, correspond to the initial solutions (with a timelimit of 5 s) reported in Table 9. This means that, if we run the MH we find a solution which is 52.9%, on average, better than those who would be obtained running the model for the same computational time. Therefore, MH can be considered a powerful tool, both in terms of efficiency and effectiveness, to address the HVRP.

S. Mancini / Transportation Research Part C 78 (2017) 1–12

11

7. Conclusions and future developments In this paper a new emerging Green Vehicle Routing Problem has been presented, the Hybrid VRP (HVRP), which is an extension of the Green VRP (GVRP) presented in the literature a few years ago, in which vehicles may operate with electric propulsion or a traditional fuel engine, paying a higher unitary cost for distances traveled in the latter mode, which can be seen as a penalty for introducing more pollution into the air other than a higher operational cost, given the fact that traditional fuel is much more expensive than electric propulsion. The main objection transportation companies make about electric vehicles usage is that, given their limited autonomy, visits to refueling stations must be planned along the routes, which implies long deviations from the original path with a consequent increase in traveled distances and route duration. This issue becomes more critical when the number of refueling stations along the road network is very limited. This problem may be overcome if hybrid vehicles are used, as they offer the possibility of covering some parts of the routes under traditional fuel propulsion in order to avoid long deviations to reach a refueling stations. A Large Neighborhood Search based Matheuristic (MH) for the HVRP has been proposed, in which neighborhoods are explored by means of the mathematical model, which is able to find the local minimum even for a large neighborhood, in a short computational time. This method is much more effective and efficient than classical Large Neighborhood Search approaches in which simple repair operators often produce lower quality solutions, and fail to reach the local minimum of the neighborhood, while more complex operators take very long computational times to reconstruct a high quality solution. The MH performances have been validated on the GVRP, for which benchmark results are available in the literature, and the method has obtained the same average performances as the best heuristic from the literature, and, in most of the cases, has reached the best known solution. When applied to the HVRP, MH shows high quality performances, and greatly improves results obtained with the model. The computational times are at least one order of magnitude shorter, even when starting from a very poor quality solution. This high robustness is an important point of strength of the algorithm. Future developments in this field could address the definition of extensions of this problem in which Traffic Limited Zones, in which it is forbidden to use traditional fuel propulsion are considered. Other developments could address the possibility of partial battery recharging and the usage of the so called regenerative breaking which allows EVs to recover a percentage of their battery charge while traveling on downhill roads. Moreover, multiple recharging technologies, with different recharging times, costs and availability could be considered, as in the work by Felipe et al. (2014) on the pure electric propulsion case. From a methodological point of view, the proposed innovative approach could be extended to other Vehicle Routing Problems.

References Ahuja, R., Ergun, O., Orlin, J., Punnen, A., 2002. A survey of very largescale neighborhood search techniques. Discr. Appl. Math. 123, 75–102. Bapna, R., Thakur, L., Nair, S., 2002. Infrastructure development for conversion to environmentally friendly fuel. Eur. J. Oper. Res. 1423, 480–496. Bektasß, T., Laporte, G., 2011. The pollution-routing problem. Transport. Res. Part B 45, 1232–1250. Bruglieri, M., Pezzella, F., Pisacane, O., Suraci, S., 2015. A variable neighborhood search branching for the electric vehicle routing problem with time windows. Electron. Notes Discr. Math. 7, 221–228. Conrad, R., Figliozzi, M., 2011. The recharging vehicle routing problem. In: Proceedings of the 2011 Industrial Engineering Research Conference. Dekker, R., Bloemhof, J., Mallidis, I., 2012. Operations research for green logistics an overview of aspects, issues, contributions and challenges. Eur. J. Oper. Res. 219, 671–679. Demir, E., Bektasß, T., Laporte, G., 2012. An adaptive large neighborhood search heuristic for the pollution-routing problem. Eur. J. Oper. Res. 223, 346–359. Demir, E., Bektasß, T., Laporte, G., 2014. The bi-objective pollution-routing problem. Eur. J. Oper. Res. 232, 464–478. Erdogan, S., Miller-Hooks, E., 2012. A green vehicle routing problem. Transp. Res. Part E 48, 100–114. Felipe, A., Ortuno, M., Righini, G., Tirado, G., 2014. A heuristic approach for the green vehicle routing problem with multiple technologies and partial recharges. Transp. Res. Part E: Logist. Transp. Rev. 71, 111–128. Franceschetti, A., Honhon, D., Van Woensel, T., Bektasß, T., Laporte, G., 2013. The time-dependent pollution-routing problem. Transp. Res. Part B 56, 265–293. Gonçalves, F., Cardoso, S., Relvas, S., Barbosa-Povoa, A., 2011. Optimization of a Distribution Network using Electric Vehicles: A VRP Problem. Technical Report. CEG-IST, UTL, Lisboa. He, F., Yin, Y., Zhou, J., 2015. Deploying public charging stations for electric vehicles on urban road networks. Transp. Res. Part C: Emerg. Technol., 227–240 Jung, J., Chow, J., Jayakrishnan, R., Park, J. Young, 2014. Stochastic dynamic itinerary interception refueling location problem with queue delay for electric taxi charging stations. Transp. Res. Part C: Emerg. Technol., 123–142 Kempton, W., Letendre, S., 1997. Electric vehicles as a new power source for electric utilities. Transp. Res. Part D: Transp. Environ., 157–175 Kuby, M., Lim, S., 2005. The flow-refueling location problem for alternative-fuel vehicles. Socio-Econ. Plann. Sci. 39, 125–145. Kuby, M., Lim, S., 2007. Location of alternative-fuel stations using the flow-refueling location model and dispersion of candidate sites on arcs. Netw. Spatial Econ. 7, 129–152. Lajunen, A., 2014. Energy consumption and cost-benefit analysis of hybrid and electric city buses. Transp. Res. Part C: Emerg. Technol., 1–15 Lin, Z., Ogden, J., Fan, Y., Chena, C., 2008. The fuel-travel-back approach to hydrogen station siting. Int. J. Hydrogen Energy 33 (12), 3096–3101. Mancini, S., 2013a. Logistics: Perspectives, Approaches and Challenges. Nova Publisher, New York, pp. 171–182 (Chapter Multi-echelon freight distribution systems: a smart and innovative tool for increasing logistic operations efficiency). Mancini, S., 2013b. Multi-echelon distribution systems in city logistics. Eur. Transp. 54 (2), 1–24. Mehrez, A., Stern, H., 1985. Optimal refueling strategies for a mixed-vehicle fleet. Nav. Res. Logist. Quart. 32, 315–328. Melkman, A., Stern, H., Mehrez, A., 1986. Optimal refueling sequence for a mixed fleet with limited refuelings. Nav. Res. Logist. Quart. 33, 759–762. Nicholas, M., Handy, S., Sperling, D., 2004. Using geographic information systems to evaluate siting and networks of hydrogen stations. Transp. Res. Rec. 1880, 126–134. Pearre, N., Kempton, W., Guensler, R., Elango, V., 2011. Electric vehicles: how much range is required for a days driving? Transp. Res. Part C: Emerg. Technol. 19, 1171–1184. Pisinger, D., Ropke, S., 2010. Handbook of Metaheuristics: International Series in Operations Research & Management Science. Springer, Berlin, pp. 1171– 1184 (Chapter Large Neighborhood Search. pp. 399–419). Schneider, M., Stenger, A., Goeke, D., 2014. The electric vehicle-routing problem with time windows and recharging stations. Transp. Sci. 48 (4), 500–520.

12

S. Mancini / Transportation Research Part C 78 (2017) 1–12

Tu, W., Li, Q., Fang, Z., Shaw, S., Zhou, B., Chang, X., 2016. Optimizing the locations of electric taxi charging stations: a spatial-temporal demand coverage approach. Transport. Res. Part C: Emerg. Technol. 65, 172–189. Upchurch, C., Kuby, M., Lim, S., 2007. A capacitated model for location of alternative-fuel stations. Geograph. Anal. 41 (1), 85–106. Waraich, R., Galus, M., Dobler, C., Balmer, M.G.A., Auxhausen, K., 2013. Plug-in hybrid electric vehicles and smart grids: Investigations based on a microsimulation. Transp. Res. Part C: Emerg. Technol., 74–86 Yamani, A., Hodgson, T., Martin-Vega, L., 1990. Single aircraft mid-air refueling using spherical distances. Oper. Res. 38 (5), 792–800. Yang, J., Sun, H., 2015. Battery swap station location-routing problem with capacitated electric vehicles. Comput. Oper. Res., 217–232 Yuan, Y., Mehrez, A., 1995. Refueling strategies to maximize the operational range of a nonidentical vehicle fleet. Eur. J. Oper. Res. 83 (1), 167–181.