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...

463KB Sizes 3 Downloads 323 Views

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.

1. Introduction In the last few years, the greenhouse effect has become a hot political topic worldwide. The European Commission has recently affirmed that the one-fifth of the EU’s total emissions of CO2 is due to road transport (Demir et al., 2014). Moreover, it has also remarked that, in the last two decades, emissions sensibly decreased in other sectors while they strongly increased in the transport one. Several governments across the world have adopted measures and regulations to reduce emissions due to vehicles emission. Such political decisions have had an important effect on the logistics industry. Many logistics companies have adopted Green Logistics projects, to reduce CO2 emissions. A reduction in pollution can be obtained in two different ways: through a better exploitation of current resources, and by using new, environmental friendly technologies. The first step in this process is to make a better use of available resources. This could be obtained by applying more efficient and sophisticated routing planning optimization methods and adopting smart distribution systems (Mancini, 2013b,a), which would help to decrease the traveling distance of vehicles and hence emissions. However, this generally results in a decline of emissions of only a few percent and the emission level of trucks and vans remains high. To strongly reduce emissions is necessary to exploit environmental friend technologies. One of the most adopted in logistics, specially in last-mile logistics, is the introduction of Electric Vehicles (EVs) in the delivery cargo fleet. This kind of vehicles, commonly adopted for urban E-mail address: [email protected] http://dx.doi.org/10.1016/j.trc.2017.02.004 0968-090X/Ó 2017 Elsevier Ltd. All rights reserved.

2

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

transport in many large cities across the world (Lajunen, 2014), are becoming more and more popular for last mile deliveries. EVs are not only much more cleaner respect to the conventional vehicles, but also more energy efficient. In fact, the energy source for EVs are renewable sources, such as solar and wind, and therefore, more sustainable. Additionally, EVs are good for the Electric Grid for two reasons. First, they may balance the excess supply of energy by renewable resources, such as wind, during night time, because most EVs are recharged at night. Second, in a smart grid, EVs can act as a power source during daytime peak hours, in which more power is requested (Kempton and Letendre, 1997). A Vehicle to Grid (V2G) system allows the charged EV batteries to give back their energy to the grid. A study on the impact of electric vehicles on smart grids have been proposed in Waraich et al. (2013). However, due to their limited range, they require frequent visits to recharging stations along their route. In fact, latest technologies allows for maximum 100 miles of autonomy, while standard daily routes are generally much longer, as stated in the study proposed by Pearre et al. (2011). Furthermore, long recharging times are requested for a fully complete recharge and the impact of these times on route planning can be extremely significant, especially when tight delivery time windows and service times are considered. A solution to overcome this issue is represented by the battery swapping strategy, according to which, is possible to reduce recharging times by replacing the vehicles batteries with fully charged ones (Yang and Sun, 2015). Although this approach is much faster, it is expensive because the number of battery to be purchased by the fleet holder is much higher than the number of vehicles, and a planning of fully charged battery repositioning at recharging stations must be carried out, with additional operative costs. Recharging times may be reduced introducing fast recharging technologies (Felipe et al., 2014), but these new technologies are very expensive. Another way to reduce loss of time, due to deviations from the original route path to visit recharging stations, is to optimize their location and improve their dissemination along the road network. Many studies have been conducted on optimal recharging stations location, among which (He et al., 2015; Tu et al., 2016; Jung et al., 2014). This approach may be adopted for public transport problems, in which the municipality cover the set up costs for new recharging infrastructures, but cannot be applied by private companies. For all the above discussed reasons, the management of battery recharging is still an open issue and the diffusion of electric vehicles is still slowed by their limited autonomy. This limitation could be overcome using hybrid vehicles, which can work both with traditional fuel and electric propulsion. The electric battery could be exploited on short routes (which can be performed without recharging the battery) or in zones in which recharging stations may be easily reached, while traditional fuel propulsion could be used in cases in which a visit to a recharging station implies a long deviation, or when the electric battery would allow almost the whole route to be covered. In fact, in such cases, it would be more advantageous to cover a few kilometers with a traditional fuel engine than to plan a battery recharging stop. The use of hybrid vehicles could constitute a fair compromise between economic interest and environmental issues. The main contribution of this paper is twofold. First, a new, and more realistic, extension of the Green VRP, the Hybrid VRP, has been introduced, in which vehicles may use both electric and traditional propulsion, and a mathematical formulation has been presented. Second, an efficient Large Neighborhood Search Matheuristic (MH) is proposed. To validate the efficiency and the effectiveness of the proposed approach it has been tested also on the Green VRP for which several heuristics have been presented in literature. Results obtained with MH are competitive with the best methods in literature. The paper is organized as follows. In Section 2 a literature review is reported. Section 3 deals with the problem description, while the mathematical formulation of the problem is given in Section 4. The Matheuristic is presented in Section 5, while computational results are reported in Section 6. Finally, Section 7 is devoted to conclusions and future developments. 2. Literature review In the last few years, increasing attention has been paid to Green Logistics, which involves the integration of environmental aspects in logistics. Many papers concerning Operations Research applications to Green Logistics have been proposed in the literature. A wide set of issues has been addressed, such as intermodal transportation, mode choice models, fleet choice and exploitation, smart distribution systems and fuel choice. For a complete survey on this subject, the readers may refer to Dekker et al. (2012). One emerging research area concerns pollution emission minimization. In Bektasß and Laporte (2011), the authors introduced the Pollution-Routing Problem (PRP), an extension of the classical Vehicle Routing Problem with Time Windows, which consists in routing a number of vehicles to serve a set of customers, and determining their speed on each route segment in order to minimize a function that includes fuel, emission and driving costs. The same problem has been addressed in Demir et al. (2012) where an Adaptive Large Neighborhood Search based heuristics approach is proposed. A time-dependent version of the PRP has been addressed in Franceschetti et al. (2013), while (Demir et al., 2014) introduced the bi-level pollution routing problem. In the nineties and the first decade of the 21st century, refineries focused on removing lead additives from gasoline, in order to preserve air quality. Biofuels based on organic waste can easily be mixed with standard gasoline. However, this technology involves the necessity of adapting engines, which can be quite rather expensive. Electric vehicles are environmentally friendly, since their engines emit almost no emissions. However, due to their limited autonomy, they are more popular for the movements of in-city goods than for medium-long range freight transport. In order to compensate for their short range, a dense power re-supply network would need to be set-up, possibly in conjunction with the possibility of changing batteries. Unfortunately, the present re-supply network is still very limited. Several works, related to power recharging and to location optimization of supply stations, have been conducted (Nicholas et al., 2004; Kuby and Lim, 2005, 2007; Upchurch et al.,

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

3

2007; Lin et al., 2008; Bapna et al., 2002), while some earlier works which were focused on military applications and considered issues related to the limited capacity of fuel tanks, may be found in the literature (Mehrez and Stern, 1985; Melkman et al., 1986; Yamani et al., 1990; Yuan and Mehrez, 1995). Only a few papers have actually dealt with VRPs pertaining alternative fuel vehicles. In Gonçalves et al. (2011), the authors dealt with a VRP concerning Pickup and Delivery (VRPPD) with a mixed fleet which consisted of both electric and traditional fueled vehicles. The objective was to minimize the total costs, which consist of vehicle related fixed and variable costs. They considered time and capacity constraints and simply added an extra time for recharging the electric vehicles batteries, where needed. However, they did not explicitly consider recharging stations location in their model. Therefore, the problem resulted in a mixed fleet VRPPD with an additional distancedependent time variable. In Conrad and Figliozzi (2011) the recharging vehicle routing problem (RVRP) was introduced, in which vehicles were allowed to recharge directly at customer locations, adding a time penalty to the route duration. Erdogan and Miller-Hooks (2012) were the first to combine a VRP with the possibility of refueling a vehicle at a station along the route. They considered a limited number of refueling infrastructures located along the network and a fixed recharging time (independent of the remaining level of battery charge). The proposed Green Vehicle Routing Problem (GVRP), considers a maximum route duration and fuel constraint. Fuel is consumed with a given rate per traveled distance and is tank is totally replenished at recharging stations. Schneider et al. (2014), present the Electric Vehicle Routing Problem with Time Windows and Recharging Stations (E-VRPTW), which can be seen as an extension of the GVRP in which time windows are considered. An extension of the E-VRPTW, where partial recharging are allowed, has been studied in Bruglieri et al. (2015). Finally, multiple recharging technologies, characterized by different recharging times and costs, have been introduced in Felipe et al. (2014). 3. Problem description The Hybrid Vehicle Routing Problem (HVRP) consists in visiting a set of customers, I, starting from a depot v 0 . The available fleet is composed of M identical vehicles. Since multiple trips are not allowed, at most M routes may be scheduled. Each vehicle must start and end its route at the depot. Each node j is characterized by a service time pj . Each pair of customers, ði; jÞ, is associated with a travel time t ij and a distance dij . The travel speeds are assumed to be constant over a link. In addition, no limit is set on the number of stops that can be made for refueling. Recharging stations may be used several times by the fleet. In order to allow that, dummy copies of each station are created, so we can assume that each dummy copy must be visited at most once following the approach adopted in the formulations of the Green VRP and of the Electric VRP with Time Windows proposed by Erdogan and Miller-Hooks (2012) and Schneider et al. (2014) respectively. The number of dummy copies to consider is a critical issue, in fact, if too many copies are added, the dimension of the network grows up too much and the model becomes inefficient, while if too few copies are added it may yield to suboptimal solutions. Schneider et al. (2014), have carried out an analysis on the number of copies necessary, and, in this work, the same number of copies has been considered. When refueling is undertaken, it is assumed that the tank is filled up to capacity. The routes start and end at the depot and must contain visits to refueling stations, if they are necessary. Each route may be covered using only the electric engine, or a portion of it can be covered using traditional fuel propulsion. The distance covered by a vehicle in traditional fuel mode, over the link connecting to customer j is defined as W j . A distance unitary additional penalty cost, q, is added when traditional fuel propulsion is used in order to dissuade this mode and to promote the use of the electric engine. A limit on the duration of the routes, T max , is imposed. An example on how the optimal routes of a VRP, GVRP and HVRP may differ (from each other) is reported in Fig. 1. More in detail, given its limited autonomy, the electric vehicle must perform a very long deviation from the optimal VRP route, in order to reach a refueling station, while the hybrid vehicle may cover the same route as in the VRP, but must cover a small portion of the path with traditional fuel propulsion, thus paying an extra cost which is still much lower than the cost of reaching the refueling station. 4. Mathematical formulation The additional notations used in formulating the HVRP are defined hereafter.          

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.