A Survey on the Electric Vehicle Routing Problem: Variants and Solution Approaches

Abstract

In order to ensure high-quality and on-time delivery in logistic distribution processes, it is necessary to efficiently manage the delivery fleet. Nowadays, due to the new policies and regulations related to greenhouse gas emission in the transport sector, logistic companies are paying higher penalties for each emission gram of /km. With electric vehicle market penetration, many companies are evaluating the integration of electric vehicles in their fleet, as they do not have local greenhouse gas emissions, produce minimal noise, and are independent of the fluctuating oil price. The well-researched vehicle routing problem (VRP) is extended to the electric vehicle routing problem (E-VRP), which takes into account specific characteristics of electric vehicles. In this paper, a literature review on recent developments regarding the E-VRP is presented. The challenges that emerged with the integration of electric vehicles in the delivery processes are described, together with electric vehicle characteristics and recent energy consumption models. Several variants of the E-VRP and related problems are observed. To cope with the new routing challenges in E-VRP, efficient VRP heuristics and metaheuristics had to be adapted. An overview of the state-of-the-art procedures for solving the E-VRP and related problems is presented.

1. Introduction

The vehicle routing problem (VRP) is an NP-hard optimization problem that aims to determine a set of least-cost delivery routes from a depot to a set of geographically scattered customers, subject to side constraints [1]. The problem was first defined by Dantzig and Ramser [2] as the Truck Dispatching Problem. VRP is a generalization of the well-known traveling salesman problem (TSP), which aims to design one least-cost route to visit all the customers. The problem has applications in several real-life optimization problems, which has led to the definition of many problem variants over the years: limited vehicle load capacity (capacitated VRP, CVRP), customer time windows (VRP with time windows, VRPTW), multiple depots (multidepot VRP, MDVRP), pickup and delivery (VRP with pickup and delivery, VRPPD), time-dependent travel time (time-dependent VRP, TD-VRP), heterogeneous fleet (mixed fleet VRP, MFVRP), etc. [3, 4]. Due to the complexity of the problem, exact procedures are only capable of optimally solving small-sized problems: up to 360 customers for CVRP [5] and 50-100 customers for VRPTW [6]. Over the years, a vast number of heuristics, metaheuristics, and hybrid procedures were proposed for solving different VRP problems.

In the past decade, the European Union (EU) has announced many new actions and regulations related to greenhouse gas (GHG) emissions in the transport sector [7]. External factors and the rise of social and ecological awareness have prompted green initiatives in many companies. Conventional internal combustion engine vehicles (ICEVs), which are dependent on limited fossil fuels, severely pollute the environment, especially in congested urban areas. According to WEEA [8], the EU intends to decrease GHG emissions by 20% and 40% by 2020 and 2030, respectively. Sbihi and Eglese [9] introduced the research field of green logistics, which deals with the sustainability of delivery processes by taking into account environmental and social factors. With the electric vehicle (EV) market penetration, many logistic companies evaluated the use of the EVs in their vehicle fleet in order to decrease GHG emissions and, therefore, reduce the charges for every emission gram of /km. EVs have several advantages compared to ICEVs: (i) they do not have local GHG emissions; (ii) they produce minimal noise; (iii) they can be powered from renewable energy sources; and (iv) they are independent of the fluctuating oil price [10, 11]. There are two basic configurations of EVs: the battery electric vehicle (BEV), which is exclusively powered from batteries mounted inside the vehicle, and the hybrid electric vehicle (HEV), which can be powered from batteries inside the vehicle or by other energy sources, most commonly internal combustion engine. The plug-in HEV (PHEV) can be recharged by connecting the plug to the electric power source. In this paper mostly BEVs for logistic purposes are considered, where two main problems come to the fore: limited driving range and the need for additional recharging infrastructure.

Due to the limited battery capacity, the range that delivery BEVs can achieve with a fully charged battery is 160-240 km [12], which is much lower than the 480-650 km range of ICEVs [13]. To achieve a similar driving range as ICEVs, BEVs have to visit charging stations (CSs) more frequently. Today, there is still a lack of CSs in the road network infrastructure, and their locations and energy demand should be planned in future infrastructural plans. For an empty BEV to become operable again, battery energy has to be renewed at a CS. This can be performed in two ways: (i) by swapping empty batteries with fully charged ones at Battery Swapping Station (BSS) or (ii) by charging at CS [14, 15]. The former process can be performed in a time comparable to the refueling time of ICEVs. In the latter process BEVs recharge their batteries at CSs by plugging into the electric power source. The recharge time depends on the state of charge (SoC) when entering the CS, the desired SoC level when leaving CS, and the charging function.

1.1. Recent Literature Reviews and Scientific Contribution

Here we present, to our knowledge, the most recent literature reviews on the E-VRP and related problems.

Juan et al. [16] presented a review regarding the environmental, strategic, and operational challenges of EV integration in logistics and transportation activities. The authors performed a comprehensive analysis of environmental challenges including the transportation impact on the pollution and EVs’ possible contribution to the reduction of carbon emissions. Regarding the strategic challenges, the authors presented issues related to CSs: battery swap technology, different charging technologies, CS location problem, limited number of charges, and charging network; issues related to the EVs: mixed fleet of ICEVs and EVs, economic challenges when integrating EVs in a fleet, and routing constraints. The solution approaches for the E-VRP were presented in general, as a class of VRP solving procedures. Similarly, Margaritis et al. [17] presented practical and research challenges of EVs focusing on battery development, lack of charger compatibility, systematic energy management, the lack of optimization procedures that could minimize the EV routing and scheduling decisions, cooling/heating usage, financial sources, novel policies, measures for EV deployment in transport services, etc.

Montoya [18] researched several variants of the E-VRP: green VRP (GVRP), E-VRP with partial recharging and nonlinear charging functions, and the technician routing problem with a mixed fleet of ICEVs and EVs. For each problem, effective solving procedures were proposed: multispace sampling heuristic, iterated local search enhanced with heuristic concentration, and two-phase parallel metaheuristic based on solving a set of subproblems and extended set-covering formulation. The authors also formulated a fixed route vehicle charging problem (FRVCP) with and without time windows to optimize the charging decisions for a route with a fixed customer sequence. The authors did not focus as much on reviewing the E-VRP literature, especially not on the procedures for solving the problem. Compared to Montoya [18], this paper did not focus as much on the FRVCP and technician routing problem or on detailed procedures used to solve those problems.

The most recent survey on the E-VRP is presented by Pelletier et al. [19] and it includes technical background on EV types and batteries, EV market penetration, EV competitiveness and incentives, and an overview of the existing research regarding EVs in transportation science. The authors provided a comprehensive review on, at the time, the latest solving procedures for the E-VRP with mixed fleet and optimal paths and covered several key papers regarding partial recharges, hybrid vehicles, and different charging technologies.

In this paper, a survey on the E-VRP is presented, which includes approaches for solving the E-VRP and related problems that emerged with BEVs integration in the logistic processes. The focus is not on the economic and environmental challenges related to BEVs. Most of the latest literature reviews were published in 2016; hence, to the best of our knowledge, there is no published research that summarizes the state-of-the-art research in the E-VRP field. In this paper we outline the following contributions:(i)a review of the recent energy consumption models that could be used in BEV routing models;(ii)an updated literature review and a concise table summary of already reviewed E-VRP variants such as GVRP, mixed fleet, BSSs, partial recharges, and different charging technologies;(iii)a review of the additional emerged E-VRP variants, which include hybrid vehicles, CS siting, nonlinear charging function, dynamic traffic conditions and charging schedule optimization;(iv)a comprehensive analysis of operation research procedures in the E-VRP, which includes an overview of the procedures employed for solving various E-VRP variants, highlighting state-of-the-art procedures, and a concise table summary of the applied procedures.

1.2. Organization of This Paper

The remainder of this paper is organized as follows. In Section 2, the E-VRP is described and the literature review of the basic problem formulation is presented. In Section 3, basic characteristics and recent energy consumption models of BEVs are described together with the BEV’s application and evaluation in the delivery processes. In Section 4, variants of the E-VRP are presented with some related problems from the literature. In Section 5, approaches for solving E-VRP are presented, which include state-of-the-art exact, heuristic, metaheuristic, and hybrid procedures. The conclusion and future research directions are given in Section 6.

2. Electric Vehicle Routing Problem

With BEV penetration in logistic distribution processes, a problem of routing a fleet of BEVs has emerged: the E-VRP. The E-VRP aims to design least-cost BEV routes in order to serve a set of customers by taking into account often used constraints: vehicle load capacity, customer time windows, working hours, etc. [3, 20]. Additionally, BEVs have the limited driving range which directly corresponds to more frequent recharging events at CSs. CSs can be built at separate locations as public CSs or mounted at customers’ locations as private CSs. The time needed to travel to a CS and the recharging time are important aspects of fleet routing, especially if customer time windows are taken into account.

To the best of our knowledge, the first research regarding the routing of an electric fleet was published by Gonçalves et al. [21], with the authors observing VRPPD using a mixed fleet of EVs and ICEVs. Refueling of the vehicle is performed at the location where the need for the refueling occurred and the total refuel time is computed, based on the total distance traveled by the EV. Conrad and Figliozzi [22] formulated recharging VRP in which vehicles with limited driving range are allowed to refuel at customers’ locations during the route to up to 80% of the vehicle’s battery capacity. The authors described an application of the model and analyzed the impact of different driving ranges, fixed charging times, and time windows on EV routes and solution quality. The results indicated that customer time windows greatly limit route distance when recharging time is long and vehicle range is constrained. Erdoğan and Miller-Hooks [23] formulated GVRP in which a fleet of vehicles is powered by alternative fuels (alternative fuel vehicle, AFV): biodiesel, ethanol, hydrogen, methanol, natural gas, electricity, etc. AFVs can refuel at separately located stations, with fixed refueling time. The authors did not consider customer time windows and vehicle load capacity constraints. New problem-specific instances were developed and two heuristics for solving the problem were applied. The results showed that the limitation of driving range severely increased the number of refueling stations and the total traveled distance in the solution. A similar problem was researched by Omidvar and Tavakkoli-Moghaddam [24], in which the authors added vehicle load constraint, customer time windows, and congestion management, as the vehicle can stay at the customer’s location during the congestion hours. The authors minimized emission and pollution costs by applying commercial software for solving small instances and two metaheuristics for solving larger instances. Schneider et al. [25] published the first research on routing a BEV fleet by taking into account possible visits to CSs and charging time dependent on the SoC level when entering a CS. The problem was formulated as E-VRP with time windows (E-VRPTW) as they took into account load, battery, and time window constraints. The authors formulated E-VRPTW as the mixed integer linear program (MILP) on the complete directed graph , where customers are modeled as graph vertices and paths between customers are modeled as graph arcs. Here, the basic model formulation is given. Let be a set of geographically scattered customers who need to be served, and let be a set of CSs for BEVs. In order to allow multiple visits to the same CS, a virtual set of CSs is defined. Vertices 0 and denote the depot, and every route begins with vertex 0 and ends with vertex . Graph is defined as , where is set of arcs, . Depending on the real-life constraints, different arc values can be interpreted as distance , travel time , energy consumption , speed , cost , etc. The binary variable is equal to 1 if arc is traversed in the solution, and 0 otherwise. The whole MILP program for the E-VRPTW with equations for load, battery, time windows, flow, and subtour constraints is presented by Schneider et al. [25].

In the VRP, it is customary for the primary objective to minimize the total number of vehicles used (1) and then to minimize the total distance traveled (2) or some other objective functions [26]. Total vehicle number is a primary objective as generally greater savings can be achieved with fewer vehicles (vehicle fixed costs, labor cost, etc.). Such an objective is contradictory as with fewer vehicles total traveled distance increases and vice versa. By taking into account the high purchase cost of BEVs, such a hierarchical objective seems justifiable in BEV routing applications [25, 27].Objective functions can be complex with simultaneous minimization of vehicle number, total traveled distance [28], total travel times [29], total routing cost and planning horizon [11, 27, 30, 31], GHG emission [32, 33], energy consumption [34–36], etc. Total routing costs of BEVs usually consist of BEV acquisition cost, circulation tax, maintenance, costs related to the energy consumption (electric energy price), cost of battery pack renewal after its lifetime, labor costs, etc. Instead of a single-objective function, some authors use multiobjective function, i.e., fuel consumption and total driving time [37], fuel consumption and route cost [38], battery swapping and charge scheduling [39], etc. An overview of different objectives in E-VRP is presented in Table 1 in the column Objective.

Table 1 

Overview of the E-VRP variants and related problems.

3. Battery Electric Vehicles in Delivery Processes

The major problem that BEVs in delivery processes are facing is the limited driving range. Grunditz and Thiringer [101] analyzed over 40 globally available BEVs, which can be categorized into small, medium-large, high-performing, and sports cars. All of the BEV models utilize lithium-based batteries, especially lithium-ion [102] with battery capacity and distance varying within the ranges 12-90 kWh and 85-528 km, respectively. An average medium-sized personal BEV has a battery capacity of 30 kWh, which is enough to travel 250 km. In the delivery process, mostly light vans and freight BEVs are used, which have a shorter driving range (160-240 km) compared to the driving range of ICEVs (480-650 km) [13, 43, 103]. The reason is that the battery has lower specific energy (130 Wh/kg) than the fossil oil (1233 Wh/kg), and the amount of energy that can be stored in the battery is much lower than in the fossil oil. Batteries mounted in BEVs are mostly the main cause of high acquisition costs and technical limitations as the battery degrades over time, resulting in decreased maximal capacity. Pelletier et al. [104] concluded that the battery should be replaced after five to ten years or after 1,000 to 2,000 cycles with large SoC variations. The authors also described factors that influence such battery degradation: overcharging, overdischarging, high and low temperatures, high SoC during storage, large depth of discharge, etc., and they presented battery degradation models that can be used in goods distribution with BEVs.

3.1. BEV Application

BEVs are more likely to be used on short distances and/or in urban areas where they are more effective than ICEVs due to the low driving speed, low noise production, frequent stops, and financial incentives. In cases when the average route length is short, such as the average FedEx route length in the USA, which is 68 km [12], BEVs can be applied directly and recharging can be performed on return to the depot. BEVs are already being applied in such occasions: DHL, UPS, FedEx, and Coca-Cola [105, 106] use BEVs mostly for last-mile deliveries as distances are shorter and vehicle loads are lower. Many companies are performing case studies of integrating BEVs in their delivery fleet. Lin et al. [60] were debating the use of BEVs in time-precise deliveries as the long recharging time at CS causes hard completion of an on-time delivery, which then significantly increases the overall routing costs. Davis and Figliozzi [10] and van Duin et al. [43] evaluated a wide range of scenarios to compare the routing costs of ICEVs and BEVs. Van Duin et al. [43] concluded that BEVs have the ability to efficiently perform urban freight transport and meanwhile to reduce the GHG emissions and noise nuisance. Davis and Figliozzi [10] reported that BEVs are not competitive if the solution to the same problem results in a higher number of BEVs than the number of ICEVs. For BEVs to be competitive, the authors pointed out a combination of several key elements: high daily distance (as much as maximum BEV driving range), low speeds and congestions, frequent customer stops, the reduction of a BEV’s purchase cost by tax incentives or technology development, and long planning horizon. On the other hand, Schiffer et al. [27] conducted a case study using freight BEVs for deliveries with a planning horizon of five years and compared the results to the delivery done by conventional freight trucks. Several characteristics of the performed case study are important: (i) delivery radius up to 190 km from the depot; (ii) CSs located at customers’ locations, which allows simultaneous charging while unloading goods; (iii) already existing strong current at CSs; and (iv) vehicles returning to the depot at least once a day. Overall, the authors concluded that there is no operational limitation when using BEVs compared to ICEVs, as the used number of vehicles and total traveled distance are competitive, and at the same time, the overall costs are lower with almost 25% less emission. A year later, Schiffer et al. [11] repeated the case study with more realistic costs when a strong current at the CS is not available and concluded that, with the higher CS investment costs, BEVs are no longer competitive. To fully assess the integration of BEVs in logistic processes, the authors pointed out three key elements: different network structures, future emission policies, and future technology development (battery capacity and charging infrastructure).

3.2. Energy Consumption

Due to the low specific energy, energy consumption should be precisely estimated in order to achieve a BEV’s maximal driving range and to reduce the overall routing costs. The energy consumption can be estimated by simulation models but due to the complexity of the E-VRP and unknown driving cycles in advance, mostly macroscopic models with several real-world approximations are applied in the BEV routing models. In the available literature, energy consumption is often estimated using longitudinal dynamics model (LDM). Here, the LDM of Asamer et al. [107] is presented. Force needed to accelerate and to overcome resistances (grade, rolling, and air) is given by (3), where is vehicle mass (mostly empty vehicle), acceleration, vehicle speed, gravitational constant, the inertia force of vehicle rotating parts (up to 5% of the total vehicle mass), road slope, rolling friction coefficient, air drag coefficient, air density, and vehicle frontal air surface. If , the vehicle is accelerating and power is needed for the movement of BEV (motor mode); otherwise, if , deceleration (braking) or driving downhill is occurring and energy is returned into the BEV’s battery as the electric engine has the ability to return the energy (recuperating mode). By process of recuperation, up to 15% of totally consumed energy can be returned [51, 108]. Electric power that comes from the battery is divided into the auxiliary power and mechanical power . Auxiliary power is spent on the electronic devices in the vehicle: heating, ventilation, light, etc., which can shorten the BEV’s range up to 30% [109]. Battery power can be computed by (4), where is the transmission coefficient between the electric motor and drivetrain, is the conversion ratio from chemical energy in the battery to electric energy, and is the conversion ratio from mechanical energy on wheels to chemical energy stored in the battery. Energy is returned into the battery only if the force is lower than zero and speed is higher than the experimentally determined value [107]. Energy consumption can be computed by the time integration of (4).

Goeke and Schneider [30] extended the energy consumption model by taking into account variable vehicle load mass when delivering goods without acceleration and braking processes. In the CVRP variants, as customers are being served, vehicle load is decreasing. The authors concluded that actual load strongly improves the solution quality, as a large number of solutions generated without taking into account load distribution tend to be infeasible due to the violation of battery capacity or time windows.

Genikomsakis and Mitrentsis [85] presented a more realistic electric engine model: the load-efficiency curve is approximated by piecewise function and the normalization factor is added to take into account the motor size. The authors used the recuperation energy factor dependent on the vehicle speed as follows: below there is no energy recuperation, beyond maximum energy is recuperated, and in between linear interpolation of energy recuperation is assumed. Therefore, regarding the energy recuperation the authors observed two cases: (i) when recuperated energy exceeds the consumption of the auxiliary devices and the excess of energy is stored into the battery; (ii) when recuperated energy is not sufficient to cover the consumption of the auxiliary devices and thus the energy is drawn from the battery. The authors compared their model on nine characteristic driving cycles (short/long, congested/uncongested, highway, etc.) to the FASTSim simulation tool [110], and this resulted in a relative energy consumption error of up to 4%. The developed energy consumption model was used to create a database of coefficients for several key characteristics: motor type, motor power, battery type, road type, road slope, road speed limit, etc.

Macrina et al. [99] developed an energy model for mixed GVRP with partial recharging and time windows. The authors modeled vehicle speed on arc with three phases : acceleration (), constant speed (), and deceleration (), resulting with either triangular function () or trapezoidal function (). The fuel consumption of ICEV is based on the similar LDM model, presented by (5) and (6), where is the load (cargo) weight, fuel-to-air mass ratio, heating value of typical diesel fuel, conversion factor, engine friction coefficient, radial engine speed, engine displacement, diesel engine efficiency, drivetrain efficiency, and travel spent in phase [40, 111]. For the consumption of BEV, the authors proposed a similar model to that presented by (7), where and are the efficiency of the electric engine in motor and recuperating mode. The authors compared the proposed energy consumption model to the basic model, where energy consumed is proportional to the distance traveled, and the energy consumption model of Goeke and Schneider [30], where acceleration and braking processes are not taken into account. The results showed that, compared to the proposed energy consumption model, basic energy consumption model and the model of Goeke and Schneider [30] produce an average relative error of 70% and 4%, respectively.

Basso et al. [92] developed an energy consumption model for two-stage E-VRP (2sEVRP) that includes detailed topography and speed profiles. The similar LDM model of Asamer et al. [107] is applied for computing the energy consumption on a road segment (link), but with time-dependent speed, varying mass, constant slope, and link distance. For each link, similar to Macrina et al. [99], the three or two characteristic phases of speed curve can be distinguished. In such a way, acceleration and breaking processes before and after the intersection are taken into account. As mass changes during the delivery, the energy consumption is expressed as a linear function of mass. Then, the shortest energy paths between all vertices pairs (customers, depot, CSs) without charging constraint are determined by applying the Bellman-Ford algorithm [112]. A nominal mass value between the mass of full and empty vehicle is used in the computation. The experiments showed an average energy estimation error of 2.28% and improved energy feasibility compared to some of the previous consumption models. The authors also reported that the use of average speed consumed at least 24% less energy.

Asamer et al. [107] analyzed the data collected from the BEV trips in order to determine the dependency between BEV energy share and average trip speed. Results showed that on lower speeds km/h, most of the energy is spent on the acceleration and auxiliary devices, while on the higher speeds km/h, approx. 70% of total energy is spent on overcoming the air drag force. The slope of the terrain in total energy consumption has the lowest share, and its value is even decreasing with the increase of average trip speed. The rolling resistance energy consumption share does not depend on the average trip speed.

Preis et al. [35] analyzed the energy optimal routing of BEV, where the routes were designed in different terrain slopes and battery capacity scenarios. The authors concluded that energy savings grow linearly with the maximum altitude difference and that decreasing the battery capacity of BEVs does not affect the recharging schedule but increases the total energy consumption. Fiori et al. [113] compared the power-based consumption model of BEV and ICEV [114]. The authors concluded that BEVs and ICEVs have different fuel/energy-optimized assignment. The faster routes increase the BEV’s energy consumption while the congested and low-speed arterial routes consume less energy. Masmoudi et al. [84] compared the realistic energy consumption model of Genikomsakis and Mitrentsis [85] to the constant consumption model (237.5 W/km) on instances for the dial-a-ride problem with EVs. The authors concluded that using the realistic model is more efficient with 0.14% difference between realistic and constant energy consumption from the best-known solutions (BKS). To emphasize the importance of the energy consumption model and energy minimization, Zhang et al. [36] compared the distance and energy minimization and concluded that the distance-minimizing objective consumes 16.44% more energy than the energy-minimizing objective.

In real-life conditions, speeds on the roads are time-dependent and can be described as the speed profile over the observed time period. Speed profile depends on the road type, driver behavior, traffic (accidents, recurrent congestions), weather conditions, etc. [115, 116]. A large number of parameters make it hard to predict a BEV’s energy consumption in different traffic scenarios. Therefore, some researchers apply data-driven approaches to predict the energy consumption of a BEV. De Cauwer et al. [117] developed a model for BEV energy consumption by applying multiple linear regression to real-world measured BEV data. The proposed multiple linear regression model (MLR) has seven features: distance, speed, energy consumed by auxiliary devices, positive elevation, negative elevation, temperature, and kinetic energy change per unit distance. Results showed that the prediction error of trip energy consumption is within 25%. De Cauwer et al. [115] applied a similar MLR model for predicting energy consumption on road segments, suited for BEV routing. The authors used a neural network based on the road-traffic and weather-related features to predict the speed profile of a road segment. The error prediction of the trip’s energy consumption is 12-14%. Lebeau et al. [52] used real-life data to model the energy consumed and recuperated on the trip through least square analysis. The used function depends on the trip duration, measured energy consumption, temperature, and correction parameter. The for the energy consumption model is 0.93 and 0.77 for the recuperated energy model. Wu et al. [118] presented an empirical and analytical method that can estimate a BEV’s instantaneous power in real time and overall trip energy consumption, with the average prediction error up to 15.6%. Fiori and Marzano [119] proposed a backward microscopic power-based LDM energy consumption model based on the known driving cycle. The acceleration, speed, and regenerative braking efficiency were modeled as functions of time. Based on the real driving BEV data, for each vehicle, six parameters were optimized: three parameters regarding the drivetrain and battery efficiency, the parameter for determining when the regenerative braking is occurring, and two traction power thresholds. The authors concluded that the proposed model is flexible enough to be applied to any kind of driving cycle and BEV type with given kinematic profile and characteristics.

Pelletier et al. [104] presented battery degradation models that could be used in BEV routing applications, as batteries are a crucial part of the economic BEV routing. Several lithium-ion battery degradation mechanisms with storage and operating conditions that affect battery lifespan were described. Barco et al. [42] applied a battery degradation model based on the temperature, SoC, and depth of discharge in their airport shuttle service with BEVs. It was observed that the consideration of the battery degradation model affects the charging patterns.

4. Variants of the Electric Vehicle Routing Problem

Many different VRP variants were researched over the years. With BEV appearance, researchers started to adapt them to the E-VRP context. Due to the specific characteristics of BEV routing, some new problem-specific variants emerged. Here, we present some of the most relevant E-VRP variants and related problems.

4.1. Energy Shortest Path Problem and the Electric Traveling Salesman Problem

The energy shortest path problem (E-SPP) and electric TSP (E-TSP) can be considered as two of the simplest forms of the E-VRP. In most of the VRP variants, graph arcs have positive weight values that can represent distance, travel time, cost, etc. To compute the shortest path between customers, the most commonly applied algorithms are Dijkstra, Bellman-Ford, , contraction hierarchies, etc. [108, 120]. By taking into account the recuperated energy of BEV, some arc weights could have a negative value, which makes most of the shortest path algorithms inapplicable. To overcome this problem, Artmeier et al. [108] formalized energy-efficient routing with rechargeable batteries as a special case of the constrained SPP (CSPP) in which the battery charge is limited with its capacity, and there is no recharging at CS. The authors adapted the Bellman-Ford shortest path algorithm with time complexity in order to solve the energy CSPP. To overcome negative graph edges, Eisner et al. [121] used Johnson’s shifting technique [122] to transform negative edge cost functions into nonnegative ones and applied Dijsktra’s algorithm with time complexity of [123]. Storandt [124] introduced CSPP for EVs with battery swapping, while Sweda and Klabjan [125] added recharging events at CSs with a differentiable charging function. Zündorf [50] solved the CSPP for BEV routing with battery constraints, different types of CSs, and nonlinear charging process. The author developed a charging function propagating algorithm in order to minimize travel time and applied CH and algorithms to solve the problem. Liao et al. [103] considered the BEV’s shortest travel path problem and presented a dynamic programming algorithm for solving the problem that runs in , where is the upper bound on the number of BSSs. More recently, Strehler et al. [126] observed recharging events in the graph vertices and arcs for the energy-efficient shortest routes of BEVs and HEVs, while Zhang et al. [127] dealt with the electric vehicle route planning with recharging problem (EVRC), which minimizes the overall travel and charging time and takes into account CSs’ locations, partial recharging, nonlinear charging functions, service time duration, and service frequency at CS.

As VRP is the generalization of the TSP, the E-VRP is closely related to the E-TSP, in which a set of customers has to be served by only one BEV. Roberti and Wen [64] formulated the E-TSP with time windows (E-TSPTW) as a compact MILP program with binary variables for recharging paths with intermediate stops. The authors observed full and partial recharge policies and applied three-phase heuristic for solving the problem. The authors solved multiple instances, among them the 13 E-VRPTW instances of Schneider et al. [25], which were solved optimally with only one vehicle. Doppstadt et al. [58] formulated the HEV-TSP with four working modes of HEV and provided new test instances. Liao et al. [103] introduced the EV touring problem, as the generalization of TSP, with the minimization of the total required time. Two scenarios with battery swap scheme were observed: the on-site station model, in which each city has a BSS; the off-site station model, in which BSS is located at an acceptable distance from the city. The authors proposed efficient polynomial time algorithms for the problem.

4.2. Heterogeneous or Mixed Vehicle Fleet

In today’s vehicle fleets, mostly ICEVs are present. Transition to an almost wholly electric fleet is a very challenging economic task. Therefore, most companies are gradually integrating BEVs into their existing ICEV fleet. Routing algorithms have to be upgraded and adapted for electric fleet characteristics, as on-time priority is harder to achieve.

Fleet size and mix VRP (FSM-VRP) was first introduced by Golden et al. [128], where the authors considered routing a fleet of vehicles with different acquisition costs and the routing cost is dependent on the vehicle type. Goeke and Schneider [30] formulated E-VRPTW and mixed fleet (E-VRPTWMF) with equal vehicle load capacities of BEVs and ICEVs, and equal BEV battery capacities. The authors compared the percentage share of BEVs to the total traveled distance obtained with three objective functions: total traveled distance, total costs with battery costs, and total costs without battery costs. The traveled distance on one battery pack was assumed to be 241,350 km (150,000 miles) with the payment of an additional $600 per kWh in case of replacement. The results showed that the BEV’s share in total traveled distance increased significantly when battery costs are not considered, while in the other two scenarios ICEVs performed most of the deliveries. Hiermann et al. [31] analyzed a similar problem, the electric fleet size and mix VRPTW and recharging stations (E-FSMFTW), with different vehicle load and battery capacities. The authors showed the overall positive influence of the heterogeneous vehicle fleet on the generalized total cost function. In most of the instances, three to four vehicle types are used in the final solution. A similar analysis on small-sized instances with different BEVs, ICEVs, and HEVs was undertaken by Lebeau et al. [52]. The authors defined seven groups of vehicle types that could be used for the delivery, from small vans and quadricycles through diesel-only or electric-only groups to a group of all vehicle types. The results showed the following aspects of routing: (i) the fleet with different vehicle types reduced the total routing costs; (ii) in the large van group, ICEVs outperformed BEVs; and (iii) HEVs showed a great application in deliveries solely made by trucks. Sassi et al. [47–49] also observed a heterogeneous fleet of ICEVs and BEVs with different load and battery capacities, different operating costs, time-dependent charges, and the compatibility of BEVs and chargers at CSs. The authors focused on the procedures for solving the problem, so they did not offer any comparison between the ICEV and BEV solutions. Hiermann et al. [96] introduced the Hybrid heterogeneous electric fleet routing problem with time windows and recharging stations (), in which the fleet consists of ICEVs, BEVs, and PHEVs. The authors compared the overall costs of the solutions obtained with a homogeneous fleet of ICEVs, BEVs, and PHEVs to the optimized solution with a mixed fleet. The gap for the ICEV fleet is the largest, with an average value of 60% in the extreme case when the fuel cost is the highest. In contrast, when the electric cost is the highest, the BEV fleet on average produces a 40% gap, while the PHEV fleet shows the lowest gap value, up to 25%. The authors pointed out that, in mixed solutions, BEVs are preferred for clustered instances, ICEVs for randomly distributed instances, and PHEVs for randomly clustered and distributed instances. Overall, the results showed that operational costs can be 7% lower in a mixed case compared to the homogeneous case.

4.3. Hybrid Vehicles

As compensation for the limited driving range of BEVs, HEVs, which have both an internal combustion engine and an electric engine, have been developed. Two main types of HEVs are present on the market: the series hybrid, in which only the electric motor drives the train and the internal combustion engine is used to recharge the battery pack, and the parallel hybrid, which uses both internal combustion engine and electric engine to drive the train, where the electric engine is more efficient in stop-and-go activities and the internal combustion engine is more efficient at high speeds. We focus on the PHEV as a version of a parallel hybrid in which batteries can be recharged by connecting a plug to the electric power source. The PHEVs have an option to decide during the route to run on either electric energy or fossil oil. This enables the visiting of customers far from the depot with almost no refuel/recharge during the route. As PHEVs have two engines, their load is heavier compared to the BEVs and ICEVs, and therefore they have a higher energy consumption rate. The time spent on the deliveries is shorter, which makes it easier to achieve time-precise deliveries and to reduce the costs of recharging, at the expense of higher costs due to fossil oil consumption.

One of the first papers that dealt with the PHEVs routing problem was published by Abdallah [41], where the author defined the problem as the plug-in hybrid electric VRPTW (PHEVRPTW). The objective of the proposed problem is to minimize the routing costs on the internal combustion engine while satisfying the demand and time window constraints. At each customer, the driver can either recharge the vehicle battery or go to the next open time window using an internal combustion engine when the electric energy has been depleted. Doppstadt et al. [58] defined the HEV-TSP in which the authors observed HEVs that are not plug-in and can only be charged while driving. Four working modes were observed, combustion-only mode, electric-only mode, charging mode, and boost mode, when the combined internal combustion engine and electric engine are used. As the authors assumed that there was not a charge left from the day before, the initial charging of the battery was set to zero. The authors presented the positive effects of using HEVs: (i) compared to the ICEV routes, the overall costs for the HEV routes in the tested instances reduced up to 13% and driving time increased up to 11%; (ii) savings depend highly on the depot location, and not on the limited battery capacity, which was an a priori assumption. Mancini [71] introduced hybrid VRP (HVRP), in which PHEVs can change propulsion mode at any time and the electric engine, rather than the internal combustion engine, is promoted. Vincent et al. [76] introduced a similar HVRP problem with PHEVs and mild-HEVs, which do not have an electric-only mode of propulsion. The results showed that PHEVs have a lower average cost per mile in all scenarios compared to the mild-HEVs and ICEVs. Hiermann et al. [96], in their E-FTW, used the electric motor of the PHEV as priority mode choice. When a time window of a PHEV route is violated, the mode could be changed at any time during the route by exchanging recharging time for the corresponding amount of fuel at no additional costs. The authors showed the positive impact of a PHEV fleet compared to ICEV and BEV fleets, especially for randomly clustered instances. PHEVs are usually represented with only 20% of the overall number of vehicles used in the solution due to their higher consumption and utility costs but, still, they constitute an important part of the fleet configuration due to their flexibility.

4.4. Partial Recharging

In the beginning, in most of the E-VRP problems, full recharge was considered when a BEV visited a CS [25]. This can be time-consuming because, depending on the SoC level, available charging technology, and battery capacity, the vehicle can charge from five minutes to eight hours [15]. Therefore, in real-life applications, partial recharging should be taken into account. The battery should be charged enough to complete the whole route or to surpass a fear that vehicle range will not be enough to perform designated tasks, the so-called range anxiety [75]. This particularly has an effect on customers with narrow time windows, where efficient charge scheduling can enable the feasibility of the route. On the economic side, significant savings can be achieved by applying partial recharging as a minimal amount of energy could be recharged during the day when electricity cost and energy network load are higher, and the rest of the energy could be replenished during the night [46, 47]. In some cases, it is natural to maintain the energy reserve. This can be done by SoC range limitation, i.e., [42, 47, 80]. Having an energy reserve seems even more important if energy consumption and range anxiety are taken into account because up to 30% of the consumed energy can be spent on BEV’s auxiliary devices [109]. Limiting SoC value also helps to preserve the battery as battery capacity decreases by overcharging and overdischarging [104].

Several papers have analyzed strategies of partial recharging and formulated the problem as E-VRPTW with partial recharging (E-VRPTWPR). Bruglieri et al. [29, 51] and Moghaddam [53] modeled the concept of partial charging in E-VRP. Felipe et al. [46], in GRVP with multiple technologies and partial recharges (GVRP-MTPR), and Keskin and Çatay [28], in E-VRPTWPR, ensured that after every change in the route configuration, at the previous CS, the BEV is charged only with the amount of energy sufficient to finish the segment of the route until the next CS or the depot. This results in the removal of certain CSs in the route, compared to full recharge strategy, and the arrival of a vehicle at the depot with an empty battery. The authors presented the positive impact of partial recharges on the total costs and energy savings. A similar procedure was applied by Sassi et al. [47–49] but with multiple charging technologies, different charging periods, and BEV chargers compatibility checks. Schiffer and Walther [73] compared full and partial recharge by solving small E-VRPTWPR test instances using commercial software. The authors concluded that, in some cases, a partial recharging strategy reduces the total traveled distance and number of visits to CSs. Montoya [18] and Montoya et al. [72] formulated the FRVCP for BEVs in which, for a fixed sequence of customers in the route, CS position and charging amount are optimized. The results on the newly derived test instances showed that good solutions tend to exploit partial recharges. A similar procedure was applied by Keskin and Çatay [79], Hiermann et al. [96], Froger et al. [68], and Schiffer and Walther [88, 89] to enhance the incumbent best solution by optimizing the charging decisions along the BEV route with a fixed customer sequence. Desaulniers et al. [57] presented the effects of partial recharging by optimally solving E-VRPTW instances containing up to 100 customers. In a case with a single recharge per route, the partial recharge reduced the routing costs by 0.97% and the number of vehicles by 2.25%, while in a case with multiple recharges per route these values are 1.91% and 3.80%, respectively, with a significant increase in the average number of recharges per route.

4.5. Different Charging Technologies

Today, multiple charging technologies are present: (i) slow, 3 kW (6-8 h); (ii) fast, 7-43 kW (1-2 h); and (iii) rapid, 50-250 kW (5-30 min) [15, 129]. To better control charging time in the E-VRP context, the selection of possible charging technology could also be optimized. This could make some customers who have narrow time windows more accessible by fast charging at previous CSs, or if the time windows are long, an economically better approach could be slow charging. Such a problem could be extended by taking into account CS working hours, time-dependent charging costs, the number of available chargers and their compatibility with BEVs, the power grid load, the charger power, etc. [46–49, 68, 80, 81, 87]. Felipe et al. [46] analyzed the effect of different charging technologies on the recharge cost. The authors concluded that none of the technologies dominated the others. The rapid and fast charging options provided slightly better results than the slow charging, but the best results were obtained when joint technologies were used as the most appropriate one could be chosen in each case. Çatay and Keskin [67] solved small-sized instances to present the insights of the quick charge option. The results showed that quick charging might reduce the fleet size and decrease the cost of energy needed to operate BEVs. Keskin and Çatay [79] formulated E-VRPTW and fast charging (E-VRPTW-FC) with partial recharges and three different charging technologies. The authors provided MILP formulation and minimized total recharging costs while operating a minimum number of vehicles. The authors compared two mixed integer programming (MIP) models: (i) model 1, in which binary variables are used as decision variables to choose which charging technology to use; (ii) model 2, based on the model of Schneider et al. [25] and Keskin and Çatay [28], in which the authors copied each station vertex three times depending on the charging technology. The results showed that in all cases model 1 produced better results with shorter computation time. Testing on larger instances showed that the fast charging option is more beneficial when customers have narrow time windows and that multiple charging technologies improved the overall results (in 28 of 29 instances), while in the instances with larger time windows the average improvement was 0.17%.

4.6. Nonlinear Charging Function

In most of the E-VRP related literature, either linear or constant charging time is considered. Most of the BEVs have lithium-ion batteries installed, which are often charged in constant-current constant-voltage (CC-CV) phases: first by constant current until approx. 80% of the SoC value and then by a constant voltage. In the CC phase, SoC increases linearly, and in the CV phase, the current drops exponentially and SoC increases nonlinearly in time, which prolongs charging time [102, 104]. In only a few reviewed papers the author considered nonlinear charging process and solved the E-VRP either by linearization per segments or by estimating charging time by a data-driven approach [50, 68, 72, 87, 95]. Zündorf [50] developed a charging function propagation algorithm for determining an EV battery-constrained route by taking into account piecewise linear and concave functions of the recharging process. Montoya et al. [72] formulated the E-VRP with nonlinear charging functions (E-VRP-NL) and presented the MILP formulation. The authors fitted piecewise linear functions for 11, 22, and 44 kW CSs to the real measured data. The results showed that neglecting the nonlinear charge can lead to infeasible or overly expensive solutions: 12% of the routes in good solutions recharged the battery in the nonlinear part, after 80% of the SoC value. Froger et al. [95] proposed two new formulations for the E-VRP-NL: (i) the improved MILP version of Montoya et al. [72] by arc-based tracking of the SoC value; (ii) a path-based model without CS replication, where a sequence of vertices (a path) between a pair of customers or depots is determined. The second one yielded much better results with shorter computing time. Froger et al. [68] explored similar formulations for the E-VRP-NL with capacitated CSs (E-VRP-NL-C): (i) the CS replication-based formulation with flow and event-based formulations for CS capacity constraints; (ii) recharging path-based formulation.

4.7. BEV Routing and CS Location

Due to the currently low BEV market share, the number of CSs installed in the road infrastructure is also relatively low. Therefore, great potential lies in the simultaneous decision-making of CS locations and BEV routes. The classic location routing problem (LRP) consists of determining the locations of the depots and vehicle routes supplying customers from these depots [130]. A modification of the LRP that deals with CS facilities is formulated as electric LRP (E-LRP) [11, 27, 73, 88].

Schiffer et al. [11, 27] presented a case study for solving ELRPTWPR. It was assumed that CSs can be located at customers’ locations and that service time can be used for recharging. Overall, results indicated the viability of combining CSs siting and BEV routing for specific cases when a delivery range is not far from the depot. Schiffer and Walther [73] compared ELRPTWPR and E-VRPTWPR solutions on down-scaled test instances and concluded that the ELRPTWPR gives a better solution in all instances. The reason comes from the fact that BEV can be charged while serving, which reduces the overall route time, as there is no travel time wasted on traveling to and from CS. Hof et al. [69] and Yang and Sun [56] addressed the problem of E-LRP but with BSSs. Results indicated that decreasing the construction cost of BSS led to the expected increase in their number in the solution. Schiffer and Walther [88] formulated the LRP with intraroute facilities (LRPIF) as a more general problem, in which intraroute facilities can be used for refueling, loading, or unloading goods. Sun et al. [131] proposed a location model for CSs without routing based on the travel demands of urban residents. For short-distance travelers, slow CSs are utilized, while for long-distance travelers, the rapid CSs are considered. The CSs’ locations are determined to maximize coverage and flow according to the concept of vehicle refueling. In the case study, the authors pointed out two critical aspects: the BEV driving range and the budget constraint.

4.8. Battery Swap

Instead of charging at CS, at specially designed BSS, empty or nearly empty batteries can be replaced with fully charged ones [14]. The main advantages of such procedure are the time in which it can be performed and the ability to recharge when energy network load and electricity costs are lower, i.e., during the night. A whole replacement procedure could last less than ten minutes, which is competitive to the refueling time of ICEVs and much faster than one of the fastest charging BEV technologies. The drawbacks of such procedure are the nonstandardized batteries and their installation in BEVs, which makes it hard to swap empty batteries with fully charged ones. Adler and Mirchandani [44] observed the E-VRP with swappable batteries, already determined BSS’s locations, a fixed number of batteries per BSS, full recharge of four hours, and a fixed swapping time of two minutes. Full battery recharge is considered so it is possible that when the vehicle arrives at the station, there is no fully charged battery available and the vehicle has to wait. First, the total routing cost is minimized and then the battery reservations are made, so that the vehicle could avoid BSSs without available batteries. Yang and Sun [56] and Hof et al. [69] simultaneously determined the locations of BSSs and the vehicle routing plan and provided different metaheuristics for solving the problem. Yang and Sun [56] formulated the problem as BSS-EV-LRP and compared the solutions of the BSS-EV-LRP instances to the corresponding best-known CVRP solutions. The total routing costs for 150 km and 300 km driving ranges increased to 7% and 3.5%, respectively, and in some cases, there was no gap between BSS-EV-LRP and CVRP solutions. Hof et al. [69] on the same problem reported that if construction costs of BSSs are equal to zero, the BSSs would be constructed in 83% of total available candidate sites.

4.9. Two-Echelon Routing Problem

Breunig et al. [93] proposed the electric two-echelon VRP (E2EVRP), in which goods are transported in two echelons: (i) in the first echelon, goods are transported by conventional freight vehicles from the depot to the satellite facilities; (ii) in the second echelon, goods are transported from the satellite facilities to the customers by light BEVs. Two vehicle types are observed in the problem: ICEVs with higher load capacity located at the depot and BEVs with lower load capacity located at satellite facilities. BEVs are used for the last-mile deliveries due to their lower pollution, noise, and size. The authors formulated the problem on the multigraph and applied exact procedures on smaller instances to solve the problem optimally and a metaheuristic procedure on larger newly developed instances. The results showed that the increase in CS density decreased the detour costs following the expression . Also, the vehicle range has a great impact on the solution quality: a range below 70 km produced infeasible solutions while a range higher than 150 km decreased recharging visits to almost zero. Jie et al. [97] analyzed a similar problem, the two-echelon capacitated E-VRP with BSS (2E-EVRP-BSS), in which, in both echelons, BEVs are used. The BEVs in the first echelon have higher load and battery capacities than the BEVs in the second echelon. The authors reported that the battery driving range is the most important aspect of routing and that it slightly depends on the number of BSSs used.

4.10. Charging Schedule

Many companies that use BEVs prefer charging the vehicles at their own facilities in order to charge the vehicles between the delivery routes and during specific periods of the day. In such occasions, there are usually a limited number of chargers at the depot, typically fewer than the fleet size; therefore, the efficient charging schedule at the depot has to be determined.

Pelletier et al. [87] considered the electric freight vehicles charge scheduling problem (EFV-CSP) at the depot for the BEVs that deliver goods to a set of customers over multiple days. The authors included multiple charging technologies at the depot, realistic charging process (piecewise linearization of nonlinear charging function), time-dependent charging costs, grid power restriction, battery degradation costs (cyclic and calendar aging), and facility-related demand (FRD) charges representing the maximal demand registered over the billing period. The fixed vehicle routes are known in advance and the number of charging events in the depot is limited to avoid impractical solutions of constantly moving a vehicle from one charger to another. The authors presented the following aspects of charge scheduling tests conducted in summer (higher electricity costs) and winter (lower electricity costs) periods: (i) model tries to keep the SoC lower when battery degradation costs are included; (ii) in summer, vehicles are rarely charged in peak hours, which results in more vehicle charging simultaneously or the use of fast chargers that retrieve more power from the grid in nonpeak hours but incur higher FRD charges; (iii) to avoid cycling the battery in high SoC, it is preferable to split the long routes into smaller ones; (iv) fast chargers are heavily used in a high BEV utilization context; (v) grid power restriction increases overall energy costs, especially in summer months, and leads to infeasible solutions, limiting the number of vehicles simultaneously charging; and (vi) total costs are always lower with larger batteries, as smaller batteries require larger discharge cycles.

Barco et al. [34, 42] also analyzed charge scheduling for assigned routes in an airport shuttle scenario. The authors defined the set of charging actions (charge profile) over the determined programming horizon and minimized overall operating costs. Sassi et al. [47, 49] also dealt with determining the charging schedule of BEVs at the depot and included time-dependent charging costs and chargers compatibility checks with BEVs. Adler and Mirchandani [44] provided an online routing model of BEVs by taking into account battery reservations to minimize the average delay of all vehicles by occasionally detouring them. Wen et al. [65] observed the service time of CSs, meaning that a CS can be visited only in some specific time period, usually working hours. The authors proposed a battery reservation model in order to charge the reserved battery at the defined time period. Sweda et al. [75] dealt with the possibility that a CS might not be available at some point in time and rerouting should be performed. The problem is formulated as the adaptive routing and recharging policies for EVs. First, optimal a priori routing and recharging policies were determined and then heuristic procedures were applied for the adaptive routing.

4.11. Dynamic Traffic Conditions

Most of the E-VRP research considers static conditions on the road network. The traffic states change recurrently, depending on the time of the day, day of the week, and season, or nonrecurrently when a traffic incident occurs, such as an accident [132–134]. TD-VRP routes a fleet of vehicles by taking into account variable travel time on the road network [135, 136]. Shao et al. [74] observed BEVs in such time-dependent context in their E-VRP with charging time and variable travel time (EVRP-CTVTT). The recharging time is fixed to 30 minutes and batteries are always charged to full capacity. The authors discretized one day into two-minute intervals and applied the dynamic Dijkstra algorithm to find the shortest travel time path when weights in the graph are not constant. The authors presented a real-life problem and solved it by applying a genetic algorithm with a running time of three hours. Omidvar and Tavakkoli-Moghaddam [24] presented a model for routing AFVs that aims to avoid routing during congestion hours, when the pollution costs are high, by waiting at a customer’s locations. These costs were computed based on the road speed profile. Mirmohammadi et al. [62] presented MILP formulation for the periodic green heterogeneous VRP with time-dependent urban traffic and time windows, which is green in terms of the emission minimization and not the utilization of AFVs. The planning horizon is divided into several periods, during which static traffic conditions are considered.

4.12. Related Problems

4.12.1. Route Scheduling and Other Electric Variants

Many researchers are dealing with the scheduling of bus/taxi routes that are fixed or could be slightly altered. Li et al. [82] addressed the mixed bus fleet management problem (MBFM) composed of electric, diesel, compressed natural gas and hybrid-diesel buses. The authors maximized the total benefit of replacement of old vehicles with new ones under budget constraints while optimizing the route assignment for each bus during the planning period. Two routing procedures were developed to solve the recharging problem: single-period routing and routing across multiple periods of a day. Analysis of the results on the case study of Hong-Kong showed a cost-effective scheme of fleet configuration in the following order: mixed fleet, electric, natural gas, diesel, and then hybrid-diesel buses. Wen et al. [65] also dealt with the scheduling of electric buses (electric vehicle scheduling problem, E-VSP) in order to efficiently serve a set of timetabled bus trips. The linear partial and full recharge were considered, with the minimization of bus numbers and the total traveled distance. Lu et al. [83] addressed the problem of optimal scheduling of a taxi fleet with mixed BEVs and ICEVs to service advance reservations. The authors presented a multilayer taxi-flow time-space network in order to minimize the total operating fleet cost.

Bruglieri et al. [77] formulated the MILP for a one-way electric car-sharing problem as the electric vehicle relocation problem (E-VReP). Relocation of BEVs is performed by workers who come with bicycles to the pickup points, put the bicycle into the BEV’s trunk, and drive the BEV to one of the delivery points. The authors balance a trade-off among the customers’ satisfaction, the workers’ workload balance, and the car-sharing provider’s objective.

Masmoudi et al. [84] formulated the dial-a-ride problem with EVs and BSSs (DARP-EV) for customers with special needs and disabilities. The authors observed a fleet of heterogeneous BEVs with equal battery capacity but different available resources for the customers: handicapped person’s seat, stretcher, wheelchair, or accompanying person’s seat. The authors adapted the energy consumption model of Genikomsakis and Mitrentsis [85] with constant speed, acceleration, road slope, and varying mass while loading/unloading passengers.

Paz et al. [86] addressed the multidepot electric vehicle location routing problem with time windows (MDEVLRPTW) with three variants: (i) with BSS (MDEVLRPTW-BS); (ii) with partial recharging (MDEVLRPTW-PR) at CSs that could also be located at customers’ locations; and (iii) with BSS and partial recharging (MDEVLRPTW-BSPR) where conventional charging is considered at customers’ locations and battery swap is performed at a separately located BSSs. Three MIP models were presented for each case. The authors determined the number and location of BSSs and/or CSs (without investment costs), the number and location of the depots, the number of vehicles, and the customers’ sequence in route. The proposed model allows minimizing the number of vehicles given a fixed number of CSs or BSSs, and vice versa. For the CS, linear recharge is considered, while the battery swap time is set to 10% of the corresponding linear recharge time. Results were presented by solving small instances with a commercial solver. The results showed that BSPR and PR variants provide better results than the BS variant and that BSPR tends to act as PR in most of the instances.

Schiffer and Walther [89] extended the ELRPTWPR model of Schiffer and Walther [73] with uncertain customer patterns related to customers spatial distribution, demands, and service time windows. The authors formulated robust ELRPTWPR (RELRPTWPR) as a MILP program. Five scenarios were observed, which represent customer profile for five working days. The authors used set-based uncertainty representation, which is covered within scenarios that consider customer patterns in daily deliveries. To achieve a computationally tractable model of robust counterpart formulation, the adversarial approach was used, which uses a finite set of scenarios. The robust model was compared to the deterministic model (selected customer pattern), median model (average customer pattern), and worst-case model (worst-case customer pattern), where first the E-LRP was solved, and then the E-VRP. The deterministic models produced infeasible solutions in 35-63% of instances and incurred, on average, 5% higher costs than the robust approach. Also, the robust approach showed the most homogeneous configuration of CSs.

Keskin et al. [80] extended E-VRPTW by considering waiting times at the CSs and proposed a solution method for solving small-sized instances. The authors assume single charger at CS, Poisson arrivals, exponential distribution of service times, first-in-first-out strategy, and penalties for late arrivals. The planning horizon is split into morning, afternoon, evening, and night period. The routing decisions are determined based on time-dependent waiting times at CSs.

Kullman et al. [81] introduced E-VRP with public-private recharge strategy (E-VRP-PP), where demand at public CSs is unknown. The queuing at CS is modeled as with first-in-first-out strategy, and represents the number of identical chargers at CS. The problem is modeled as a Markov decision process, and an approximate dynamic programming solution is proposed. The authors extended the FRVCP formulation of Montoya et al. [72] with time-dependent waiting times and discrete charging decisions. Dynamic and static routing policies are proposed, where the static policies are incorporated as base policies for dynamic routing. Using dual bound on the optimal policy, the authors concluded that proposed routing policies are within 4.7%. The authors pointed out that the use of a public-private strategy outperformed the use of just private CSs, as overall savings were 20% higher.

4.12.2. Green and Pollution Routing Problem

The green vehicle routing problem focuses on the reduction of routing pollution on the environment. The key idea is to promote the use of sustainable energy sources and minimize overall emissions. Regarding the emission of BEV, Álvarez Fernández [137] indicated that emission savings of BEVs can be significant, up to 80% in some scenarios compared to ICEVs. The authors presented a model for prediction of a BEV’s GHG emissions linked to a route for each country, according to empirical measurements and prediction of routes’ energy consumption.

Erdoğan and Miller-Hooks [23] first defined the GVRP with AFVs and full refuel, where also BEVs were included. Therefore, the E-VRP can be observed as a special case of the GVRP. The authors provided MILP formulation of the problem with replication of CSs vertices. Koç and Karaoglan [59] improved the MILP formulation of Erdoğan and Miller-Hooks [23] by removing dummy vertices that represent multiple visits to refueling stations. Instead, a binary variable related to traveling between two vertices by using the station in between is introduced (arc-duplicating). This formulation improved the average gap to the optimal integer solution on all scenarios compared to the model of Erdoğan and Miller-Hooks [23] by 25%. The authors proposed a combination of exact branch-and-cut procedure and simulated annealing with several moving strategies to solve the problem. Leggieri and Haouari [70] formulated the problem as a nonlinear MIP, which was linearized to derive MILP formulation by applying a reformulation-linearization technique. The authors used similar arc-duplicating to that of Koç and Karaoglan [59] but defined the binary variable as the sum of binary variables for the direct arc between the customers and binary variables for the same direct arc but with stations inserted. The proposed formulation without additional constraints was able to solve to optimality 97.5% of test instances, which substantially outperformed the formulations of Erdoğan and Miller-Hooks [23] and Koç and Karaoglan [59]. In addition, the same formulation with prepossessing consistently outperformed the exact branch-and-cut procedure of Koç and Karaoglan [59]. Andelmin and Bartolini [66] took a different approach for the exact solving of GVRP based on a set partitioning formulation in which columns correspond to feasible routes. The problem is modeled using a multigraph in which vertices correspond to customers and each arc represents a possible sequence of stations visited when traveling between customers. The authors added weak subset row inequalities, subset row inequalities, and k-path cuts, with the latter reported as the crucial part of the algorithm efficiency. The results showed tight lower bounds and high efficiency of the algorithm as instances with up to 110 customers were solved to optimality. Bruglieri et al. [94] presented a path-based approach for solving the GVRP, based on two phases. In the first phase, a set of feasible paths is generated, removing from it all those paths that are dominated. In the second phase, through MILP formulation, paths are combined to create routes for the final solution. The authors compared their results to previously observed research papers of Erdoğan and Miller-Hooks [23], Koç and Karaoglan [59], Leggieri and Haouari [70], and Andelmin and Bartolini [66]. On small-sized instances, the proposed algorithm produced optimal solutions with the lowest running time while on the larger instances it produced 2.3% worse solutions compared to Andelmin and Bartolini [66], but in much less time.

Poonthalir and Nadarajan [38] formulated biobjective fuel efficient GVRP (F-GVRP) with the minimization of route costs and fuel consumption using goal programming. The authors analyzed the impact of varying speed on fuel consumption where speed is modeled using triangular distribution. To efficiently solve the problem, the authors applied a population-based metaheuristic. The tests were performed by constraining the route cost of the solution to the BKS of GVRP and then the fuel minimization was performed. The authors pointed out that the use of varying speed instead of constant speed reduced the fuel consumption on average by 20.48%, and in half of the instances it reduced the number of vehicles. The lower number of vehicles is a consequence of the lesser amount of time spent at the refueling stations as less fuel is consumed, so the vehicle can visit more customers. The authors pointed out that if the route cost constraint was relaxed, most of the instances might have better route costs than the BKSs of GVRP. The authors also analyzed the impact of speed interval size on the fuel consumption and concluded that fuel consumption and average speed increased with larger interval size.

Normasari et al. [100] formulated capacitated GVRP (CGVRP) as a MILP program and solved the problem by applying a simulated annealing heuristic. Compared to the GVRP, the results indicated that the number of tours and total costs for the CGVRP increased by 38.8% and 32.5%, respectively.

Macrina et al. [33] investigated the green mixed fleet VRP with partial recharges and time windows (GMFVRP-PRTW), in which polluting emissions are modeled through the functions of both traveled distance and vehicle load. The problem was formulated as a MILP and solved by a commercial solver for small instances and by an iterative local search metaheuristic for large instances. The overall pollution emissions are maintained below the designated threshold value, to simulate emission limits of freight vehicles. Macrina et al. [99] presented MILP formulation for a similar GMFVRP-PRTW problem but added more realistic consumption model presented in Section 3.2, where the authors considered different charging technologies (single technology per CS), cost of energy recharged at the depot, and the fuel cost. A large neighborhood search metaheuristic is applied to solve the problem.

Koyuncu and Yavuz [98] compared the node- and arc-duplicating MILP formulation of mixed fleet GVRP (MGVRP) by taking into account the following: (i) a mixed fleet of AFVs and ICEVs; (ii) load and time window constraints; (iii) internal station visits (at customer) and external station visit (at separately located stations); (iv) fixed full, variable full, or variable partial refuel at the station; (v) site-dependent refueling rate; and (vi) the visit of two or more consecutive stations being prohibited. In the node-duplicating formulations [23, 25], each refueling vertex is replicated at most times, where is the number of customers. In the arc-duplicating formulation, two subsets of arcs are contained: arcs between customers (direct arcs) and refueling arcs which contain inserted station between customers [70]. In the preprocessing step for both formulations, unfeasible arcs are removed and tighter lower and upper bounds for auxiliary variables and number of vehicles are presented. On the GVRP test instances of Erdoğan and Miller-Hooks [23], the authors reported that the node- and arc-duplicating formulations yielded an average optimality gap of % and %, respectively. On the GVRP instances of Yavuz [138], significant conclusions could not be drawn.

As BEVs have no local emissions, the E-VRP is closely related to the minimization of GHG emissions, where a problem-specific GVRP variant called the pollution routing problem (PRP) was introduced [32, 40]. Bektaş and Laporte [40] formulated the PRP as a nonlinear MIP where the main objective is the optimization of vehicle speed and GHG emissions. Demir et al. [32] concluded that speed optimization improves the PRP solution and minimizes fuel consumption and driver costs.

It can be noted that several papers are dealing with fuel consumption and pollution emissions and not AFVs, but the term green VRP is used to indicate a more ecologically aware routing problem [62, 139].

4.12.3. More General Formulation of the E-VRP

The E-VRP can be formulated more generally as the VRP with intermediate stops (VRPIS), in which the vehicle visits intermediate or intraroute facilities to replenish/unload the goods or to refuel [55]. The LRPIF extends the problem in terms of the facilities location optimization [88]. The service time depends on the load and fuel level at the arrival at the facility. Schneider et al. [55] presented the E-VRP with recharging facilities (EVRPRF) as a special case of VRPIS, where BEVs have a fixed recharging time and limited vehicle load capacity without time windows constraints.

4.12.4. Multiobjective Variants

Most of the researchers are using a single-objective function that represents travel distance, total costs, energy consumption, total time, etc. The multiobjective variants are still relatively scarcely applied. Generally, the solution of the multiobjective problem is not the optimal solution for all of the objectives, but rather it is satisfactory in those terms.

Demir et al. [37] introduced the biobjective PRP with the minimization of fuel consumption and total driving time. The fuel consumption is based on a model similar to the one presented by (6) but with constant speed. An adaptive large neighborhood search was used to generate nondominated/Pareto optimal solutions while an additional hybrid procedure was used to obtain the solution to the biobjective model. The experiments showed that there is no need to prolong driving time to achieve a significant reduction in fuel consumption.

The already mentioned F-GVRP of Poonthalir and Nadarajan [38] (Section 4.12.2) used biobjective minimization of route costs and fuel consumption. Wang et al. [91] considered multiobjective function to determine optimal BEV routes between origin-destination pairs. Three objectives are minimized: travel time, energy consumption, and charging cost. Travel time includes driving time, queuing time, and charging time, while charging cost includes electricity cost, service cost, and parking fee. The multiobjective model is transformed into a single objective by applying a fuzzy programming approach for objectives’ membership functions and fuzzy preference relations to obtain weighting coefficients of each objective. A genetic algorithm is designed to obtain an optimal solution. The influences of various driver trade-offs (weighting conditions) among different objectives are explored using Pareto curves. It was observed that driver trade-offs have an effect on the travel times, and almost no effect on the energy consumption and charging costs.

Amiri et al. [39] optimized the location of the battery swapping in the BSSs network and the charging schedule of depleted batteries. The multiobjective model minimizes three objectives: battery charging and power loss costs, deviation from the nominal voltages, and network capacity releasing. The nondominated sorting genetic algorithm is applied to solve the problem. The results on the test case showed that an optimal charging schedule in terms of cost and network constraints can be achieved. Besides, the electricity load profile is more leveled.

Bruglieri et al. [77] focused on the E-VReP with three objectives: minimization of workers used to relocate the BEVs (service provider objective), maximization of relocation requests served (users’ objective), and minimization of the longest routes’ duration (workload balance). The problem is solved by approximating the Pareto optimal front based on a two-phase approach: first, feasible solutions are generated by randomized search heuristics, and in the second phase, to find a set of nondominated solutions, the MILP formulation is solved by epsilon constraint programming.

4.13. Summary of the E-VRP Variants

Table 1 presents a summary of the E-VRP variants and related problems in the available literature. In total, 79 papers were presented in the table. The problem characteristics are observed in the following order: vehicle load (cargo) capacity (48), linear charging/refueling function at station (34), customer time windows (34), partial recharging strategy (30), fixed refueling/recharging time (19), different charging technologies at CSs (17), heterogeneous (mixed) fleet (16), energy consumption model (15), BSS , nonlinear charging function , dynamic traffic conditions , location routing , capacitated stations , and HEVs . We can point out that only a few papers are dealing with (i) simultaneous CS siting and BEV routing (E-LRP) where great savings can be achieved in companies whose business plans include CS investments, as E-LRP generally gives better routing configuration than the E-VRP [73, 88]; (ii) nonlinear charging function, negligence of which can lead to infeasible solutions and additional penalty costs [68, 72]; (iii) HEVs (PHEVs) that have an advantage in several routing scenarios but are being researched only recently due to the higher complexity of the problem as HEVs can change their propulsion mode at any time during the route [96]; (iv) dynamic traffic conditions, which significantly influence BEVs’ energy consumption [74]; (v) different charging technologies at CSs, which could reduce the overall costs by selecting the best possible charger option in each occasion [46, 79]; (vi) CS related challenges: capacitated CS [68, 75, 80], CS reservations [44], and public-private recharge strategy [81]; and (vii) robust E-VRP variant resistant to uncertainties in travel time, service time, and demand [89].

5. Problem Solving Methods

Since VRP is a well-researched problem, a large number of procedures for solving the problem have been proposed. Due to the NP-hardness of the problem and a large number of customers in real-life problems, most of the procedures used in real-life applications are heuristics, metaheuristics, and hybrid combinations. For small-sized problems, a great number of exact procedures have been proposed. Many of the VRP procedures in the available literature are with an adaptation applicable for solving the E-VRP. An overview of the applied procedures for solving the E-VRP and related problems is presented in Table 2.

Table 2 

Overview of the applied procedures.

5.1. Feasibility

When creating or modifying a VRP solution, two types of search procedures regarding the feasibility of the solution can be observed: allowing only feasible solutions or acceptance of infeasible solutions. The infeasible solution means that some customers are served without satisfying all of the problem constraints. Related, feasible procedures are searching in the space of feasible solutions, while infeasible procedures allow searching the space of infeasible solutions, which broadens the search. In the E-VRP, the feasibility mostly refers to vehicle load, customer time window, and battery (energy) constraints. In infeasible procedures, objective function is often defined with penalization coefficients, which are updated during the search process. At the beginning of the search, infeasible solutions are allowed in order to search a larger solution space. As the search process comes to an end, penalties for infeasible solutions increase. An example of such an objective function is given by (8) where is the total traveled distance in the solution ; penalties coefficients are , , and ; is a diversification penalty; and , , and are values of constraint violation, respectively, load capacity, time windows, and battery capacity [25]. An important part of such a function is the efficient computation of constraint violations [25, 31, 96, 140].Several researchers in the E-VRP field are allowing infeasible solutions during the search [11, 25, 30, 31, 69, 84, 88, 96] while some others are only allowing feasible solutions [23, 28, 46–49, 52, 68, 72, 79]. For example, in the E-VRPTWPR of Keskin and Çatay [28], the customer is inserted in the solution if it satisfies the load and partial time window constraints. If the solution is energy infeasible the CS is inserted to make the solution completely feasible. If the former one was not able to produce the energy feasible solution, the procedure returns to the latest feasible solution. Similarly, Montoya et al. [63] and Froger et al. [68] address the fuel infeasible route by solving a CSPP on the auxiliary graph. A list of papers allowing infeasible or feasible-only solutions is presented in Table 2 in column NF.

5.2. Exact Procedures

Exact procedures are able to find an optimal solution for a smaller number of customers: up to 360 customers for CVRP and up to 100 customers for VRPTW [5, 6]. Some of the most used exact algorithms are branch-bound-cut-and-price, dynamic programming (DP), MIP, set partitioning (SP), etc. [1]. In the E-VRP, many researchers formulate the problem as a MIP and solve the small instances with commercially available software. Koç and Karaoglan [59] applied the branch-and-cut algorithm coupled with simulated annealing and four improvement strategies to solve the GVRP. The authors reported that 22 out of 40 instances were solved to optimality. Hiermann et al. [31] solved the heterogeneous E-VRP by applying a branch-and-price algorithm on small instances up to 15 customers where the pricing problem was solved using a bidirectional labeling algorithm. The algorithm was able to find an optimal solution for all of the instances in a couple of minutes, with only three instances having computation time longer than five minutes. Desaulniers et al. [57] presented branch-price-and-cut algorithms for four E-VRPTW variants depending on the full or partial recharge and on the single or multiple recharges per route. On the 696 instances with 25, 50, and 100 customers, 98%, 90%, and 27%, respectively, were solved optimally. Schiffer et al. [11, 27], Schiffer and Walther [88], and Hiermann et al. [96] applied DP on resource CSPP (RCSPP) to find the optimal facility route configuration in order to enhance the solutions produced by metaheuristic, while Pourazarm et al. [54] applied DP to find the fastest BEV routes in the single and multivehicle environment. Similarly, Montoya et al. [63] used pulse algorithm to determine the optimal placement of refueling facilities in a route with fixed customers. Similar to the FRVCP formulation of Montoya [18], Keskin and Çatay [79] in their E-VRPTW-FC formulated the MILP program to optimize charging related decisions in BEV route with the fixed customer positions. For every iteration, the authors tried to solve this problem exactly by applying CPLEX solver. Liao et al. [103] applied a branch-and-bound algorithm to solve the EV touring problem. Montoya et al. [63, 72] and Hiermann et al. [96] used set partitioning on the pool of routes to select the best subset of routes guaranteeing no overlapping in customers assignments. Froger et al. [68] used a branch-and-bound algorithm to solve the set partitioning model and to discard the selection of routes that are infeasible or for which the total time is underestimated.

5.3. Heuristic Procedures

Heuristic procedures seek to solve the problem based on the specific knowledge of the problem, usually suboptimal or close enough to a satisfactory solution. In the research field of VRP, heuristic procedures can be split into constructive and improvement heuristics.

5.3.1. Constructive Heuristics

Constructive heuristics are often used to generate an initial solution by either serial or parallel route construction. Solutions are constructed in a greedy way, which often produces solutions of the VRP that are 10-15% far from an optimal solution [4]. In the E-VRP, constructive heuristics are modified and adapted to the BEV characteristics and feasibility checks. Some well-known general constructive heuristics used in E-VRP are presented.

The modified savings method of Clarke and Wright (MCWS) [141] with the insertion of CSs is used to generate an initial solution of several E-VRP and related problems [23, 32, 52, 69]. The algorithm starts with the creation of back-and-forth routes for each customer. Next, CSs are inserted in a least-cost manner, in either feasible or infeasible routes. Then, the routes that increase the savings the most are merged together and the procedure is repeated until all customers are served.

The sweep algorithm [142] inserts customers in the active route in a circular manner, resulting in efficient space division. The customers are sorted based on the value of a polar angle between the depot and randomly chosen point. Then, customers are inserted in the active route at the position causing minimal increase until a violation of route constraint, when usually a new route is opened [25, 56, 97].

Nearest neighbor heuristic (NNH) starts from the depot, and in each iteration the customer with the least-cost increase from the previously selected customers is added to the route. The route is terminated when some/any constraints are violated and a new route is opened [76].

Route-first cluster-second [143] approach constructs a giant infeasible TSP tour, usually with NNH, which is then split into several routes. The split procedure can be solved as the SPP on the acyclic graph. The procedure generates energy infeasible solutions, so the CS placement and charging schedule need to be determined to make the solution energy feasible [63, 68, 72].

The cluster-first route-second [144] approach divides the customers into clusters and then each cluster is solved as a subrouting problem. One example used in GVRP is the density based clustering algorithm (DBCA) proposed by Erdoğan and Miller-Hooks [23], which determines customer clusters based on their density, by ensuring that each customer has at least some predefined number of customers in the radial neighborhood. In the second step to obtain the routes, the MCWS is performed.

Many researchers are adopting some form of the insertion heuristics, which have some characteristics of already reviewed heuristics. The general idea of the insertion heuristics is to iteratively insert or add customers in available routes. In each step of the algorithm, an unserved customer together with a route and its position in the route are determined to least increase the objective function. This selection can be deterministic but is often stochastic to select at random one from -best inserts to diversify the solution, especially in population metaheuristics. In most feasibility ensured cases when the energy constraint is violated, a CS is added to the route to make it energy feasible, and when the other constraints are violated such as vehicle load capacity or time windows, the new route is opened, and the procedure is repeated. Examples of adaptations to the E-VRP are push forward insertion heuristics (PIH) [41, 145], constructive-k-PseudoGreedy [46], [28, 30, 31, 33, 36] modified insertion heuristics, and charging routing heuristic (CRH) [47–49], in which first the charging scheme at the depot is determined, and then the least-cost feasible customer insertions are performed.

5.3.2. Improvement Heuristics

Improvement heuristics or local search (LS) explores the neighborhood of the incumbent solution, searching for a better solution. The neighborhood is explored by applying perturbation moves based on the composite neighborhood operators. The local search stops when no improving solution can be found in the neighborhood of the incumbent solution, which is then called local optima. A great variety of researchers have performed LS procedures to intensify the search, coupled together with perturbation moves to escape the local optima. Often, the perturbation moves are similar to the neighborhood operators used in LS phase. Most of the classical VRP neighborhood operators [3, 4] are used in the E-VRP, but some additional neighborhood operators have been developed. Depending on whether the operators perform on only one route or between the routes, they can be divided into intraroute and interroute operators. Here, we present some of the most used LS neighborhood operators in the E-VRP literature:(i)2-Opt [146] replaces two arcs with the two new ones, with a possibility of route direction reversal.(ii)2- [147], unlike the 2-Opt operator, avoids the reversal of a route direction.(iii)4-Opt [84, 148] replaces four consecutive arcs with four new ones.(iv)Or-Opt [149] replaces three arcs with three new ones such that a sequence of three vertices is relocated.(v)Exchange (Swap) [150] exchanges the two vertices in the solution.(vi)FacilityInsertion [55] and StationInsertion [11] operators insert refueling or loading facility at best position in between two consecutive facility visits where fuel or load violation occurred. The similar operator InsertRemoveIf [31] inserts CS into the part of the route where the violation occurred to make the route energy feasible and repair violations, but it also removes redundant CS visits. Schneider et al. [55] also proposed the facility-related removal and exchange operators FacilityReplacement and Exchange.(vii)RechargeRelocation [46] removes all visits to the CSs and then tries to insert CSs at positions leading to a better solution. Similarly, the GlobalChargingImprovement [72] operator removes all CS visits in the route, and if the route is energy feasible, it stops; otherwise, it solves FRVCP to optimize CSs positions and charging amount.(viii)Relocate [150] removes one vertex from the solution and inserts it into another position in the solution. Felipe et al. [46] proposed a modification of the relocation operator called Reinsertion, , for solving GVRP-MTPR. Binary parameters and denote whether the first-improve or best-improve strategy is used and whether the values of savings are updated or not. The RemoveTwoInsertOne [84] operator removes two vertices at random and inserts them one by one into another vehicle at best possible positions.(ix)CrossExchange [151] removes consecutive vertices from route to route , and consecutive vertices from route to route [84]. The RemoveSequence [84] operator is similar as it removes consecutive vertices from only one route and inserts it into another one.(x)Resize [31], RelocateAndResize [31], -OptWithModeChange [58], and VehicleSwap [96] change the vehicle type assignment or the propulsion mode in the heterogeneous E-VRP and HEV routing problems.(xi)StationInRe [25] inserts and removes CSs from the solution in so that if an arc to the CS is already in the solution, it is removed; otherwise it is inserted.

The order of the operators affects both the solution quality and the execution time. Often, there is a question whether to make a first better or the best improvement move in the LS phase [152]. Some E-VRP researchers perform one or the other, and some combine the two approaches, i.e., best move of first 100 [11] or first 50 [31].

Neighborhood operators explore only the space of the immediate vicinity of the current solution, which often leads to local optima. In order to search a larger solution space, Shaw [153] proposed the large neighborhood search (LNS) heuristic, which is based on destroying and repairing the solution. Heuristic efficiency depends on the implemented destroy and, especially, repair operators. The main drawback of the heuristic is the repeated use of the same destroy and repair operators, which can lead to local optima. A modification of LNS by Christiaens and Vanden Berghe [154] is giving state-of-the-art results on a wide range of VRP variants. In the E-VRP, LNS is successfully implemented by Sassi et al. [47–49] and Zhang et al. [36], but most often it still needs some other guiding heuristic to escape the local optima. Macrina et al. [99] applied an LNS scheme to solve GMFVRP-PRTW. A list of used destroy and repair operators in the E-VRP is presented as a part of the adaptive LNS version in the appendix.

5.3.3. Route or Move Evaluation

Operator evaluation is an important part that highly affects execution time and has to be efficiently implemented. For example, in VRPTW, for most of the basic neighborhood operators, the load capacity check and time window check can be computed in constant time , if for each customer endogenous variables, i.e., remaining load capacity, latest arrival time, etc., are tracked. In the E-VRP, there is also a constraint regarding the vehicle battery capacity (battery feasibility), which also influences the time window feasibility. Therefore, time window evaluation is more complex. For example, when the customer is trying to be inserted into the route, the load and energy check can be computed in constant time for most of the basic neighborhood operators. The endogenous variables for vertices between the insertion position and next CS have to be recomputed in order to know the charging time at next CS. The complexity of such approach is where is the number of customers between the position of the insertion and next CS. If there is no CS in the second part of the route, no recalculation is needed. By storing additional information about the recharging position and its amount and by the use of general concatenation operators presented by Vidal et al. [155], Hiermann et al. [31] reported that all of the basic neighborhood operators can be evaluated in . Goeke and Schneider [30], instead of storing additional information about recharging, used a surrogate cost function whose evaluation is less demanding. In the search process, the authors store best moves in the composite neighborhood and evaluated either exact or surrogate cost, and at the end if there is any move with the surrogate cost, its exact cost is evaluated and the best move is chosen. Hiermann et al. [96] presented a different approach for E-FTW in which customer sequence without CSs is determined, and then in the route evaluation part time-consuming CSs insertions and charging schedule are performed. Whenever the route needs to be evaluated, the heuristic finds the routing cost by applying a labeling algorithm to determine CSs placement and a greedy algorithm to determine charging amount and running mode of PHEV. As such an evaluation is time-consuming, the authors performed the following: (i) use neighborhood restrictions to consider only promising arcs for each customer vertex; (ii) keep track of move evaluations in cache memory; and (iii) lower bounds to first evaluate the move as ICEV in , and then if it improves the objective function, perform the route evaluation.

5.4. Hybrid Metaheuristic Procedures

Many researchers employ metaheuristics to continue the exploration after the first local optima occurrence. Metaheuristics can be defined as heuristics guiding other heuristics and can be divided into neighborhood-oriented metaheuristics and population metaheuristics. Most often, multiple metaheuristics and heuristics are combined together and adapted to the problem; therefore, the term hybrid is used in such occasions.

5.4.1. Neighborhood-Oriented Metaheuristics

Neighborhood-oriented heuristics iteratively explore the neighborhood of the incumbent solution. Here, we present the most used neighborhood heuristics in the E-VRP literature.

Simulated annealing (SA) [156] controls the search by the temperature parameter, which mimics the physical cooling process of the material. Most often, it is used as a control heuristic for accepting solutions produced by other heuristics. The probability of accepting solution , which is worse than the current best solution , is computed as (Boltzmann function [157]), where is the temperature parameter. The higher the temperature value, the higher is the probability of accepting a deteriorating move to escape the local optima. The procedure starts with high temperatures at the beginning and then dynamically lowers the temperature as the search process comes to an end. In the E-VRP, most often initial temperature is selected such that a % worse solution than the initial solution is accepted with a probability of [25, 28, 30, 46].

Tabu search (TS) [158, 159] is a well-known metaheuristic, which uses memory to forbid the exploration of the space that has already been explored. It escapes the local optima by accepting a worse solution if it is the best in the explored space. Arcs that are deleted from the solution are stored in the tabu-list, which prohibits the reinsertion of the deleted arcs into a specific part of the solution for a designated number of iterations, tabu-tenure. In some implementation cases, when the new best solution is found, some of the moves and arcs are lifted and removed from the tabu-list, aspire criterion [25]. In the E-VRP, TS with neighborhood operators is used to intensify the search [25, 35, 48].

Variable neighborhood search (VNS) [160] is a metaheuristic that systematically changes the neighborhood each time there is no improvement in LS phase. The changes are based on predefined neighborhood structures. In each step, the new neighborhood operator is used, which increases the vicinity from the incumbent solution. The selection of neighborhood operators and customers can be deterministic (variable neighborhood decent (VND)) or probabilistic. Schneider et al. [25] successfully adapted it in the E-VRPTW for diversification of the solution based on the cyclic exchange moves, resulting in 15 different neighborhood structures. The VNS in the E-VRP was also combined with exact algorithms like branch-and-price [51], evolutionary algorithms [84], deterministic LS operators [64, 72], and adaptive mechanisms (AVNS) [11, 27, 55, 69].

Iterated local search (ILS) [161] is based on successively repeating LS on the incumbent solution. When LS ends up in local optima, perturbation move is applied to escape from local optima. The perturbation move can be scaled dynamically to overcome local optima. The effectiveness of the procedure highly depends on the local search procedures used. Montoya et al. [72] applied ILS together with VND for solving E-VRP-NL with perturbation performed by splitting perturbed giant TSP tour. Macrina et al. [33] used typical ILS with some of the LS operators used in the perturbation phase to solve GMFVRP-PRTW. Zhang et al. [36] used ILS and LNS in an ant colony framework for solving GVRP. The LNS is used as a perturbation mechanism while the LS is used to further improve the solution. Sassi et al. [47, 49] also used a multistart combination of ILS and LNS but without the LS phase. Sassi et al. [48] applied similar ILS paradigm, but they called it iterated tabu search (ITS), as they used TS for the intensification and LNS for the perturbation. A similar ITS procedure is applied by Doppstadt et al. [58] for the HEV-TSP.

Adaptive large neighborhood search (ALNS) [162, 163] is the extension of LNS heuristic in which, throughout the search process, different procedures for destroying and repairing the solution are selected in an adaptive manner based on their past performance. The idea is to have versatile destroy and repair operators in order to escape the local optima. For each destroy and repair operator, the adaptive weight and score are increased in the following order (scores ): (i) if a new best solution is found; (ii) if the solution is better than the previous one; (iii) if the solution is worse than the previous one but it is accepted due to the acceptance criteria; and (iv) no increase if the solution is not accepted. The score of worse solution is greater than the score of better solution , to reward unimproved solutions that diversify the search [28]. For each iteration, the weights are updated based on the previous weight values and new score values. As a great number of the researchers used ALNS for solving the E-VRP and related problems (Table 2), a vast number of destroy and repair operators were applied. Some of them originated from Pisinger and Ropke [162] but also some new problem-specific ones regarding the CS have been developed [28, 30–32, 79, 88]. The description of the applied destroy and repair operators in E-VRP is presented in the appendix. Schiffer and Walther [89] applied a parallel ALNS version of Schiffer and Walther [88] to solve the RELRPTWPR. The OpenMP was used to solve multiple scenarios in parallel, where the customer destroy and repair operators, as well as the LS and DP, are performed in parallel. The change of facility configuration, adaptive evaluations, and stopping criterion are not parallelized.

SA is one of the most used metaheuristics for acceptance criteria, in order to escape the local optima, but some others can also be applied. Tiwari et al. [164] proposed the Cauchy function (CA) for the acceptance of the solution as it gives more opportunities to escape from local optima. The record-to-record (RR) principle [165], which adds negative value to the objective function of the best solution so far, can also be used to escape the local optima. Several papers in the E-VRP field applied such heuristics for solution acceptance [48, 49, 84].

5.4.2. Population Metaheuristics

Population metaheuristics are based on the natural selection to evolve a population and survival of the fittest. They have been widely applied in the VRP field: genetic algorithm [166], path relinking [167], scatter search [168], ant colony [169], bee colony [170], particle swarm optimization [171, 172], etc.

Genetic algorithm (GA) consists of a pool of individuals (VRP solutions) called a population, which goes through the process of evolution. Such evolution process includes evaluation of individuals, selection of parents, crossover, mutation, and a replacement of the old population by a new one. The important aspects of the algorithm are the representation of the VRP solution, efficient management of population diversity, and fitness function to compare the individuals in the population. For further algorithm description and application, we refer to Potvin [173]. The ant colony (AC) algorithm is based on ant behavior when they search for food. As they search, they leave pheromones, which accumulate over the most used path, leading to the optimal path between the food and the nest. AC algorithm consists of three procedures: pheromone initialization, probability adjustment, and pheromone update. Particle swarm optimization (PSO) mimics the behavior of organisms such as fish schooling or bird flocking. Physical movement of the individual (particle) in the swarm is tracked and then well-balanced moves are applied to adapt the individual to the global and local exploration. Each particle has a position in the swarm and corresponding velocity in order to optimize the problem following personal best solution and global best solution [172].

Barco et al. [34, 42] applied differential evolution to solve the E-VRP and charge scheduling. The initial individual represents the valid route assignment and charging profile. Next, mutation based on OR and XOR operators is performed on three randomly selected individuals. The crossover is a random 2-point while the selection chooses between the target individual and new individual, which is generated by the crossover. Alesiani and Maslekar [45] applied GA for solving a single BEV routing problem. The individual is represented as a two-dimensional binary array, where the rows represent the routes and the columns represent CSs. The initial population contains feasible individuals, individuals without CSs, and individuals containing battery CSPP routes. The authors select parents at random or based on the probability and apply 2-point crossover. Mutation is applied at random to remove, insert, or exchange CS, while the LS is performed to improve the solution. Shao et al. [74] applied GA without detailed description for solving real-life EVRP-CTVTT instances. The authors encoded the individual with all visited vertices including customers, depots, and CSs. The crossover operator tries to create a new individual that satisfies all the constraints. Also, elitism is applied to preserve the elite individuals. The authors reported the average computational time of three hours. A similar algorithm with the LS after the crossover and mutation is applied for solving E-VRP by Shao et al. [90] to obtain an optimal operation scheme including routes, charging plan, and driving paths. Masliakova [61] applied the GA to find the optimal placement of CSs for electric bus routes. The individual is represented as the sequence of CSs and bus stops visited, and a crossover is implemented as the swap between the subroutes. Mutation moves the CS to the nearest bus stop or creates a new bus stop with CS. For the fitness function, the author used an AC approach, which ranks the individuals in the population. Joo and Lim [78] used the AC approach to determine the energy consumption and speed in energy-efficient BEV routing. Each ant selects the next vertex with a greedy stochastic search, which combines heuristic information about energy consumption, speed, and pheromones. The authors concluded that AC strategy is suitable for routing BEVs in terms of energy efficiency. An efficient version of evolutionary algorithm and VNS (EVO-VNS) is applied by Masmoudi et al. [84] for solving the DARP-EV. Each population is divided into groups in such a way that each group contains good and bad solutions. Further on, solutions from each group are selected in the VNS phase. If the best solution in the current group is not improved for three iterations, the crossover between the selected solutions and the rest of solutions in the group is performed. The rest of the individuals in the population are generated by applying destroy-repair procedure with regret heuristic. Poonthalir and Nadarajan [38] applied particle swarm optimization with greedy mutation operator and time-varying acceleration coefficient (TVa-PSOGMO) to solve the F-GVRP. Each particle is encoded as the sequence of customers who are converted to feasible routes by inserting depots and refueling stations. Time-varying inertia and acceleration coefficients are employed to better manage diversification and intensification exploration in the multiobjective environment. Additionally, to prevent premature convergence, a greedy mutation operator (GMO) is performed. Hiermann et al. [96] used hybrid GA (HGA) with LS to obtain good assignment and routing decisions for E-FTW. The authors represented the solution as a set of vehicle routes and used a biased fitness function. During the search, the authors maintained two populations of feasible and infeasible individuals. Parents are selected using binary tournament scheme from those two populations. The recombination is performed on the parent routes combined in a giant single tour, with an OX crossover and a split algorithm to obtain a complete solution. The mutation phase is performed as LNS with destroy and repair operators similar to ones of Hiermann et al. [31]. In the end, the solution is improved by the LS procedure.

5.5. Summary of the Applied Procedures

The overview of the applied procedures in the E-VRP field is presented in Table 2, where 78 papers are presented. Searching in the space of infeasible solutions is applied in 17 papers, while the LS for the intensification is applied in 25 papers. Most of the researchers do not use any special criteria for accepting the solution, but 17 papers use SA as acceptance criteria, three RR and one CA. For the initial solution, mostly some form of the universal insertion heuristic is applied (18). Then follow MCWS (11), NNH , random , split procedure , and sweep . Regarding the applied heuristics for solving the problem, the order is the following: ALNS (14), GA (10), LNS , TS and ITS , ILS , AC , VNS , AVNS , VND , PSO , and differential evolution . From all of the exact procedures, the is the most applied one (12) and then there is the branch-bound-cut-and-price . Forty-three papers use commercially available software like CPLEX (32), Gurobi , MATLAB , and XPRESS for solving either the whole problem or some subproblem during the search, mostly charging decisions.

We can note that a vast number of researchers use hybrid heuristics procedures with only a few focusing on the exact procedures [31, 57]. It can be noted that most of the procedures that are giving state-of-the-art results are using exact procedures at some point during the search, usually for (i) optimal determination of charging visits and charging amount [79, 88, 96] and (ii) set partitioning over the set of routes at the end of the search [65, 72, 95, 96].

The ALNS is one of the most applied procedures and it is giving the BKSs on various problem variants. Among them, we can point out the following: ALNS of Hiermann et al. [31] (HPH) and ALNS of Goeke and Schneider [30] (GS) for heterogeneous fleet, ALNS of Keskin and Çatay [28] (KC1) and ANLS of Keskin and Çatay [79] (KC2) for partial recharges and different charging technologies, and ALNS of Schiffer and Walther [88] (SW) originally designed for the LRPIF (ELRPTWPR) which is also giving BKSs on the E-VRPTW and E-VRPTWPR. Schiffer and Walther [88] compared their algorithm on the E-VRPTW instances to the HPH, GS, KC1, and VNS-TS of Schneider et al. [25] (SSG). The authors found nine BKSs with an average gap of 0.27% to the previous BKSs and average running time of 3.77 minutes. On the E-VRPTWPR, compared to the KC1, the authors found 41 new BKSs. The proposed ALNS of SW is efficiently parallelized to solve the RELRPTWPR [89].

Most of the population-based metaheuristics are not efficiently implemented so they are not giving high-quality solutions. But, the HGA of Hiermann et al. [96] (HHPV) gives some of the BKSs on the E-VRPTW, E-VRPTWPR, and E-FTW. The authors compared their results on the E-VRPTW and E-VRPTWPR to the SSG, HPH, GS, and KC1. On the E-VRPTW instances, the HHPV found 11 BKSs with 1.36% increase in the vehicle number, 0.46% distance gap from the previous BKSs, and an average running time of 11.59 minutes. On the E-VRPTWPR, the authors found 28 new BKSs. By comparing the results of HHPV and SW on E-VRPTW, we can conclude that in 20 instances SW outperformed HHPV, while in the opposite case the HHPV was better in 17 instances.

6. Conclusion and Future Research Directions

With the increase in popularity of electric-powered vehicles, many logistics companies have started to integrate electric vehicles into their delivery fleets. Routing a fleet of electric vehicles for delivering goods was formulated as the electric vehicle routing problem. Besides the often used load (cargo) capacity and time window constraints, E-VRP routing models have to account for the limited driving range of BEVs, which directly corresponds to the more frequent recharging needs at CS.

In this paper, a literature review on recent developments in the E-VRP field is presented. We consider BEVs to be vehicles powered only from batteries mounted inside the vehicle. A general overview of the BEV’s characteristics for goods delivery includes the driving range, battery capacity, application, and case studies. As energy consumption estimation is an important part of BEV routing, we summarized the recent research on theoretical and data-driven energy consumption models.

Due to the BEVs specific characteristics, new E-VRP variants have emerged: a mixed fleet of electric and conventional vehicles, partial recharging, simultaneous CS siting and BEV routing, nonlinear charging, different charging technologies, battery swap technology, hybrid vehicles, green routing, etc. The development of efficient heuristics was necessary to find optimal or near-optimal solutions to the new routing problems. We reviewed the state-of-the-art exact, heuristic, and hybrid procedures applied for solving various E-VRP variants. The adaptive large neighborhood search [88] and hybrid genetic algorithm [96] produced high-quality solutions on various E-VRP variants. Most of the state-of-the-art procedures use exact algorithms during the search to determine the optimal placement of CSs and the charging amount.

From the literature review, it can be noted that the electric vehicle routing research community has grown rapidly in the last few years and many problem variants have already been explored. Nevertheless, we can highlight possible future directions as follows.

By our observation, there is a lack of papers regarding case studies and application cases, where actual E-VRP models could be evaluated, and some meaningful insights could be drawn. Several energy consumption models were reviewed, but only few are predicting realistic energy consumption at the road segment level in the road network. Only recently, have researchers started to integrate a nonlinear charging process, CS location problem, and hybrid electric vehicles into the E-VRP models. A few papers addressed the problem of CS capacity, as a limited number of BEVs can charge simultaneously at a CS. Several variations were observed: waiting times [68, 80, 81], battery reservations [44], and adaptive routing with uncertainties in CS availability [75, 81]. Furthermore, the following real-life characteristics have not been sufficiently studied in the E-VRP field: dynamic traffic conditions [74], uncertainties in demand, travel time, time windows, service time and charging time [75, 81, 89], compatibility of BEVs with chargers in CSs [48, 49], and recharging at public or private CSs [81].

For future research regarding the solution procedures, we can highlight the following. Although there have been a couple of exact procedures developed for GVRP, only a few exact procedures have been proposed for the E-VRP and its extensions. A great number of researchers applied a population metaheuristic to solve the problem, but only a handful of them produced high-quality solutions in reasonable computation time. To improve the computation time, parallelized procedures could be used [18, 89]. Furthermore, more general procedures could be investigated as most of the state-of-the-art procedures are using complex problem-specific heuristics.

Appendix

ALNS Destroy and Repair Operators

Here, we present an overview of the used destroy and repair operators in ALNS metaheuristic for solving EVRP variants and related problems. We divide the destroy operators in terms of whether they include customer removal (C), station/facility removal (S), or route removal (R).

Destroy operators:(i)Add, Swap, Drop, SwapPerfect, and SwapPerfectOut (S) [88, 174] operators change the configuration of CS in the solution. As the name says, Add, Swap, and Drop operators add, remove, and exchange the CSs in the solution. SwapPerfect and SwapPerfectOut [88] close arbitrarily or no longer needed CSs and iteratively open new ones with the best insertion strategy.(ii)BatteryBestUse (C) [97] operator detects the route with the highest charge level at the return to the depot and the customers between the depot and last-used BSS are removed from the solution.(iii)BSSCostBased (C) [97] operator removes customers that have high traveling costs of arcs connecting to the BSS. First, the random BSS is selected, and the customers with the highest cost are removed from the solution.(iv)BSSCustomerRoute and BSSProprotionRoute (R) [97] operators aim to remove routes that visit more BSSs and fewer customers (BSSCustomerRoute) and routes that have a relatively large cost of visiting BSSs (BSSProprotionRoute).(v)Cluster (C,S) [30, 163] operator removes clusters of geographically close vertices. First, a random route is selected and divided into two clusters. One of the two clusters is selected at random, and all the vertices in the cluster are removed. If more vertices need to be removed, first the random vertex from the removed cluster is chosen. Then, the route closest to the chosen vertex is selected, and the procedure is repeated.(vi)ExpensiveStation (S) [79] operator removes stations that incur the biggest costs.(vii)FullCharge (S) [28] operator removes stations at which BEVs are fully charged.(viii)HistoricalKnowledge (C,S) [32, 163] operator remembers the best position cost for each vertex during the search and then the vertex that has the highest difference from its best remembered position is removed from the solution.(ix)InefficientRouteAndNeighbour (R) [31] operator removes one route based on its average cost per unit transferred and then the next distance closest route is selected for removal.(x)LeastUsedStation (S) [79] operator aims to reduce frequent visits to CSs. This operator ranks the CSs in the nondecreasing order of charge amount and removes the top of them.(xi)Neighborhood (C,S) [32] operator removes the vertices that have the highest shares in routes average distance.(xii)NodeNeighborhood (C,S) [32] operator randomly selects a vertex in the solution and removes it and vertices in rectangular area around it.(xiii)Random (C,S) [32] operator randomly removes vertices.(xiv)RandomRoute and GreedyRoute (R) [28, 32, 174] operators randomly or greedily remove routes. GreedyRoute operator removes routes that have the lowest number of customers.(xv)Related, RandomAndRelated, WorstAndRelated [31, 162], and TimeRelated [65] (C) operators remove a set of customers that are in some sense related. The first vertex is selected at random or based on the detour cost, and the next vertices are selected deterministically or by a roulette wheel selection based on the relatedness measure. The relatedness measure of Pisinger and Ropke [162] takes into account only the distance between the vertices, while Hiermann et al. [31] took into account distance and fractions of time window violations between the vertices.(xvi)Remove customer options (C,S) [28]: (i) only customer removal; (ii) removal of customer with preceding CS; and (iii) removal of customer with succeeding CS.(xvii)RequestGraph (C) [163] operator creates the undirected graph where each customer is a vertex and arc weight corresponds to the number of times the arc is used in currently known best solutions. Then, similar to the Shaw procedure, related customers are removed from the solution.(xviii)Shaw (C) [153] operator tries to remove customers that are similar to each other by taking into account geographical distance, demand difference, earliest start time difference, and assigned route difference. The first customer is selected at random while the other customers, instead of a deterministic way, can be chosen in a probabilistic way from most similar customers to the current one [30]. Demir et al. [32] and Keskin and Çatay [28] use ProximityBased, TimeBased, and DemandBased as a special case of Shaw removal, where either only distance, time, or demand is taken into consideration. ProximityBased operator is used by Hiermann et al. [96] in LNS but it is called a similar operator.(xix)SinglePoint, TwoPoint, and Binary (C,S) [56] operators remove vertices in the BSS service zones. First, the route is selected at random. In the SinglePoint version, the first position in the route is selected at random, and then all the customers between the selected position and the depot are removed. The Binary version is similar to SinglePoint, but the first position is selected as the middle point in the route. In the TwoPoint version, two positions are selected at random, and customers between these positions are removed.(xx)StationBased (C,S) [56] operator removes all vertices connected to a randomly selected BSS.(xxi)StationVicinty (C,S) [30] operator removes customers and CS in the radial vicinity of the selected CS.(xxii)Target (R) [31] operator selects a vertex with the highest contribution to the total distance of its assigned route. Then, the routes that are closest to the selected vertex are removed from the solution.(xxiii)Worst (C,S) [163] operator removes the solution vertices that have high costs where cost is computed as routing cost between preceding and succeeding vertex. Applied examples are WorstDistance and WorstTime [28, 32] and the operator of Goeke and Schneider [30] where a more complex objective function is used. Yang and Sun [56] add a feasibility measure to the worst selection method in order to prefer the removal of vertices that improve the feasibility of a solution.(xxiv)WorstChargeUsage (S) [28] operator removes stations at which BEV arrives with a relatively high charge level.(xxv)Zone (C,S) [28, 32] operator is based on the removal of vertices in predefined zones. The whole configuration region is divided into smaller zones. The operator at random selects zones and removes customers and stations from it.

Repair operators:(i)Greedy (C,S) [162] operator iteratively inserts removed vertices at the best possible position in the solution. Modification of Greedy operators are (i) greedy randomized adaptive search (GRASP) [175] insertion, which instead of selecting the best possible insertion selects random insertion from a pool of best insertions (Goeke and Schneider [30] applied GRASP only on customer vertices); (ii) GreedyWithNoise [32] operator adds noise to the cost function in order to diversify the search (similar to RR acceptance strategy); (iii) AdvanceGreedy [56] operator adds feasibility value to the cost function; (iv) SequentialNode [31] operator inserts vertices based on the removal order followed by restricted candidate list of five best options; and (v) SequentialPerturbed [88] operator extends the SequentialNode operator by perturbing insertion costs. For the station insertions, Keskin and Çatay [28] detect the first arc in the route that is energy infeasible and then try one of the following three operators to insert the station: Greedy inserts the best station between the detected arc and the previous customer; GreedyWithCopmarison compares the Greedy insertions on the detected arc and on the previous arc and chooses the one that increases the cost at least; and the Best operator compares the Greedy insertions on arcs between the detected arc and the previous CS and chooses the best one.(ii)Regret (C,S) [162] computes the -regret values for each removed vertex as a difference between the insertion cost in the best route and -best route. This so-called regret value quantifies the amount of increase in routing cost if a vertex is inserted in -best route and not the first best route, meaning that these vertices have a lower number of cost-effective insertions. In EVRP mostly regret-2 and regret-3 heuristics are used. Similar to the greedy case, noise or feasibility values can be added to the cost function: RegretWithNoise [32], AdvanceRegret-k [56].(iii)SemiParallelConstruction and SemiParallelInsertion (C,S) [31] operators work as the initial construction algorithms and either construct new routes or insert customers in the existing routes.(iv)TimeBased and ZoneBased (C,S) [28, 32] operators insert vertices in the best route positions that minimally increase the overall route time. In ZoneBased case, the vertices can be inserted only within a certain zone.

Nomenclature

2E-EVRP-BSS:	Two-echelon capacitated electric vehicle routing problem with battery swapping stations
2sEVRP:	Two-stage electric vehicle routing problem
AC:	Ant colony algorithm
AFV:	Alternative fuel vehicle
ALNS:	Adaptive large neighborhood search
BEV:	Battery electric vehicle
BKS:	Best-known solution
BSS:	Battery swap station
BSS-EV-LRP:	Electric vehicles battery swap stations location routing problem
CA:	Cauchy function
CGVRP:	Capacitated green vehicle routing problem
CRH:	Charge routing heuristic
CS:	Charging station
CSPP:	Constrained shortest path problem
CVRP:	Capacitated vehicle routing problem
DARP-EV:	Dial-a-ride problem with electric vehicles and battery swapping stations
DBCA:	Density based clustering algorithm
DP:	Dynamic programming
E-FSMFTW:	Electric fleet size and mix vehicle routing problem with time windows and recharging stations
E-LRP:	Electric location routing problem
ELRPTWPR:	Electric location routing problem with time windows and partial recharges
E-SPP:	Energy shortest path problem
E-TSP:	Electric traveling salesman problem
E-TSPTW:	Electric traveling salesman problem with time windows
E-VReP:	Electric vehicle relocation problem
E-VRP:	Electric vehicle routing problem
E-VRP-NL:	Electric vehicle routing problem with nonlinear charging functions
E-VRP-NL-C:	Electric vehicle routing problem with nonlinear charging functions and capacitated CSs
E-VRP-PP:	Electric vehicle routing problem with public-private recharging strategy
E-VRPTW:	Electric vehicle routing problem with time windows
E-VRPTW-FC:	Electric vehicle routing problem with time windows and fast charging
E-VRPTW-SF/MF/SP/MP:	Electric vehicle routing problem with time windows, S (single recharge), M (multiple recharges), F (full recharge), P (partial recharge)
E-VRPTWMF:	Electric vehicle routing problem with time windows and mixed fleet
E-VRPTWPR:	Electric vehicle routing problem with time windows with partial recharging
E-VSP:	Electric vehicle scheduling problem
E2EVRP:	Electric two-echelon vehicle routing problem
EFV-CSP:	Electric freight vehicles charge scheduling problem
EV:	Electric vehicle
EVFSMVRPTW:	Electric vehicle fleet size and mix vehicle routing problem with time windows
EVO-VNS:	Evolutionary variable neighborhood search
EVRC:	Electric vehicle route planning with recharging
EVRP-CTVTT:	Electric vehicle routing problem with charging time and variable travel time
EVRPRF:	Electric vehicle routing problem with recharging facilities
F-GVRP:	Fuel efficient green vehicle routing problem
FRD:	Facility-related demand
FRVCP:	Fixed route vehicle charging problem
FSM-VRP/MFVRP:	Fleet size and mix vehicle routing problem
FSMVRPTW:	Fleet size and mix vehicle routing problem with time windows
FSMVRPTW-EV:	Fleet size and mix vehicle routing problem with time windows for electric vehicles
GA:	Genetic algorithm
GHG:	Greenhouse gas
GMFVRP-PRTW:	Green mixed fleet vehicle routing problem with partial battery recharging and time windows
GMO:	Greedy mutation operator
GRASP:	Greedy randomized adaptive search procedure
GVRP:	Green vehicle routing problem
GVRP-MTPR:	Green vehicle routing problem with multiple technologies and partial recharges
H2E-FTW:	Hybrid heterogeneous electric fleet routing problem with time windows and recharging stations
HDARP:	Heterogeneous dial-a-ride problem
HEV:	Hybrid electric vehicle
HEV-TSP:	Hybrid electric vehicle traveling salesman problem
HEVRP-TDMF:	Heterogenous electric vehicle routing problem with time-dependent charging costs and a mixed fleet
HGA:	Hybrid genetic algorithm
HVRP:	Hybrid vehicle routing problem
ICEV:	Internal combustion engine vehicle
ILS:	Iterated local search
ITS:	Iterated tabu search
LDM:	Longitudinal dynamics model
LRP:	Location routing problem
LRPIF:	Location routing problem with intraroute facilities
LS:	Local search
MBFM:	Mixed bus fleet management problem
MCWS:	Modified Clark and Wright savings method
MDEVLRPTW (BS/PR/BSPR):	Multidepot electric vehicle location routing problem with time windows (battery swapping/partial recharging)
MDVRP:	Multidepot vehicle routing problem
MDVRPI:	Multidepot vehicle routing problem with inter-depot routes
MGVRP:	Mixed fleet green vehicle routing problem
MILP:	Mixed integer linear program
MIP:	Mixed integer program
MLR:	Multiple linear regression
MTFSP:	Mixed taxi fleet scheduling problem
NNH:	Nearest neighbor heuristic
PHEV:	Plug-in hybrid electric vehicles
PHEVRPTW:	Plug-in hybrid electric vehicle routing problem with time windows
PIH:	Push forward heuristics
PRP:	Pollution routing problem
PSO:	Particle swarm optimization
RCSPP:	Resource constrained shortest path algorithm
RELRPTWPR:	Robust electric location routing problem with time windows and partial recharging
RR:	Record-to-record procedure
SA:	Simulated annealing
SoC:	State of charge
TD-VRP:	Time-dependent vehicle routing problem
TS:	Tabu search
TSP:	Traveling salesman problem
TSPTW:	Traveling salesman problem with time windows
TVa-PSOGMO:	Particle swarm optimization with greedy mutation operator and time-varying acceleration coefficient
VND:	Variable neighborhood decent
VNS:	Variable neighborhood search
VRP:	Vehicle routing problem
VRP-HFCC:	Vehicle routing problem with a mixed fleet of conventional and heterogenous electric vehicles including new constraints
VRP-MFHEV:	Vehicle routing problem with mixed fleet of conventional and heterogenous electric vehicles
VRPIRF:	Vehicle routing problem with intermediate replenishment facilities
VRPIS:	Vehicle routing problem with intermediate stops
VRPPD:	Vehicle routing problem with pickup and delivery
VRPTW:	Vehicle routing problem with time windows.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this article.

Acknowledgments

The research has been supported by the European Regional Development Fund under grant KK.01.2.1.01.0120 and has also been partially supported by other two sources: Croatian Science Foundation under project IP-2018-01-8323 and project KK.01.1.1.01.0009 (DATACROSS) within the activities of the Centre of Research Excellence for Data Science and Cooperative Systems supported by the Ministry of Science and Education of the Republic of Croatia.