Development and implementation of algorithms for vehicle routing during a no-notice evacuation

: This article develops and implements time-dependent shortest path and assignment algorithms for vehicle routing and assignment during a no-notice evacuation. A time-dependent shortest path algorithm with arc labeling is designed to improve computer storage space efficiency. Visual Basic for Application (VBA) is used to improve time efficiency of both shortest path and assignment algorithms. Performance of the algorithms is analyzed and compared through implementation in VBA and the General Algebraic Modeling System (GAMS). Since VBA seamlessly integrates with Microsoft Excel and enables efficient data manipulation, the algorithms implemented in VBA are more efficient and obtain optimal solutions faster than those implemented in the GAMS. The shortest path algorithm with arc labeling implemented in VBA may be used for real-time vehicle routing for large road networks during no-notice evacuations.


Introduction
In a generalized vehicle routing problem (VRP), multiple vehicles travel from origins (depots; Montoya-Torres, Franco, Isaza, Jiménez, & Herazo-Padilla, 2015) to destinations (demand points). One of the objectives of the VRP is to minimize total travel time or the latest arrival time for a fleet of vehicles while meeting demand at destinations. The VRP has a wide range of applications, e.g. bank or postal deliveries and school bus routing. During a no-notice evacuation of a metropolitan area, emergency response vehicles must be swiftly deployed to multiple demand points to manage traffic, including blocking ramps and lanes and diverting traffic. Once an emergency response vehicle arrives at an intended demand point, it is stationed at the point until the evacuation completes, which may take from a few hours to one day. Certain demand points may need multiple emergency response vehicles whereas others may need only one vehicle. Before a no-notice evacuation begins, emergency response vehicles are usually stationed at depots controlled by state or local government agencies, e.g. the state department of transportation. Since road conditions change rapidly, there is a strong need for real-time emergency vehicle routing during a no-notice evacuation.
The objective of this article is to investigate effective algorithms which may be implemented to identify optimal vehicle assignments and routes for the VRP in real time. The main contributions of this article are as follows: (a) Implementing and integrating time-dependent shortest path (TDSP) and assignment algorithms to minimize total travel time for emergency response vehicles during a no-notice evacuation of a metropolitan area or other emergency response and evacuation situations; and (b) Improving time and space efficiency of the algorithms for real-time vehicle assignment and routing in large road networks.
(i) The TDSP algorithm with arc labeling requires less computer memory compared to node labeling, and improves space efficiency. The TDSP with arc labeling and assignment algorithms are applied to road networks with up to 2,000 nodes, and are able to identify the optimal vehicle routes and assignments within eight minutes. (The largest road network studied previously had 1,000 nodes; Almoustafa, Hanafi, & Mladenović, 2013); and (ii) The algorithms are implemented in Visual Basic for Applications (VBA ® ). Input and output of the algorithms are stored and processed using Microsoft Excel ® . The seamless integration of VBA ® and Excel ® greatly improves time efficiency and enables real-time vehicle assignment and routing.
The VRP has many input parameters, including the number and origins of vehicles, the size and structure of a road network, the number and location of destinations, and vehicle travel times (dynamic or static). Numerous methods were introduced in recent years to solve the VRP; these methods stipulated various conditions for one or more parameters. One of the most common conditions was the upper limit for the size of road networks. Other conditions include static travel time and upper limit for the number of destinations. These methods become ineffective (solutions far from optimal) or inefficient (cannot identify a good or optimal solution within an acceptable amount of time) when stipulated conditions do not hold. To develop effective and efficient algorithms to solve VRPs with many parameter values remain a considerable challenge (Vidal et al., 2014).
Since exact methods that identify optimal vehicle assignments and routes are either ineffective or inefficient for generalized VRPs, heuristic methods including the genetic algorithm (Haghani & Jung, 2005), tabu search (Ichoua et al., 2003), branch and price (Almoustafa et al., 2013), and column generation algorithm (Spliet & Gabor, 2012) were studied. Haghani and Jung (2005) presented a genetic algorithm to solve a pick-up or delivery VRP with soft time windows. The study considered multiple vehicles with different capacities, real-time service requests, and dynamic travel times between destinations. Ichoua et al. (2003) conducted experiments to solve the VRP with time-dependent travel speeds, which satisfy the first-in-first-out (FIFO) property, using a parallel tabu search heuristic. Almoustafa et al. (2013) improved a branch-and-bound method to solve the asymmetric distance-constrained VRP suggested by Laporte et al. (1987). Chen et al. (2006) formulated a realtime TDVRP with time windows as a series of mixed integer programming models and developed a heuristic algorithm, which included route construction and improvement. Spliet and Gabor (2012) proposed a formulation of a time window asymmetric VRP and developed two variants of a column generation algorithm to solve the linear programming relaxation of this formulation. Kritzinger et al. (2012) applied variable neighborhood search algorithm to solve the TDVRP with time windows. Maden et al. (2010) proposed a heuristic algorithm for the VRP to minimize total travel time.
Road networks with different sizes were also studied. Laporte, Nobert, and Taillefer (1988) examined a class of asymmetrical multi-depot VRPs and location-routing problems for a network of 80 nodes. Haghani and Jung (2005) solved the TDVRP for networks with 30 demand nodes over 30 time intervals. Almoustafa et al. (2013) solved an asymmetric distance-constrained VRP for a network of 1,000 demand nodes. In summary, previous research predominately focused on developing heuristic methods for subsets of VRPs or TDVRPs. Effective and efficient algorithms which may be applied to generalized TDVRPs to obtain optimal vehicle assignments and routes were not available. Most algorithms and methods developed in previous research were not tested using real-world road networks and could not be validated for effectiveness or efficiency. The objective of this article is to provide effective and efficient algorithms for generalized TDVRPs in a no-notice evacuation. The rest of this article is organized as follows: Section 3 presents the problem and methodology. Section 4 validates the algorithms using two real-world and five simulated road networks. Section 5 concludes the article with future research directions.

Problem definition
Let A, B, Γ represent sets of vehicles, depots, and demand points, respectively, in a time-dependent road network. There are total |A| vehicles stationed at |B| depots at the beginning of a planning period, which is a time period during which vehicles must be dispatched to demand points to meet the demand. Some or all of the |A| vehicles need to be dispatched to |Γ| demand points, each of which requires d γ vehicles, where γ represents a demand point, ∈ Γ. Let c i,j,t be the cost (time) it requires for a vehicle α, ∈ A, to travel from node i at time t and to node j. i, j ∈ V, where V is the node set in the road network. B, Γ ⊂ V. Let (i, j) represent an arc that originates from node i and points at node j, (i, j) ∈ E, where E is the set of arcs in the road network. c i,j,t = ∞ if (i, j) ∉ E. When i = β, β is a depot and ∈ B, c ,j,t = ∞, ∀j, if α is not ready to travel from β at time t.
In the TDVRP, the first objective is to identify the earliest time for α to arrive at a demand point γ from its depot β. α may travel from β at any time after α is ready for travel. Suppose α may travel from β at time t 1 and arrive at γ at time t 2 . Alternatively, α may travel from β at time t 3 and arrive at γ at time t 4 . According to the FIFO property, t 2 < t 4 if t 1 < t 3 . Therefore, α should travel from β as soon as α is ready for travel. Let t β be the time at which α becomes ready for travel at β. c ,j,t = ∞ when t < t β , http://dx.doi.org/10.1080/23311835.2016.1163767 ∀j. c ,j,t << ∞ when t ≥ t β and ( , j) ∈ E. The first objective of the TDVRP is therefore to identify the TDSP for α to travel from β at t β to γ. Equation (1) is a model whose optimal solution is the TDSP between β and γ when α travels from β at t β .
Subject to: The second objective of the TDVRP is to minimize total travel time. Let s , represent the optimal value of Equation (1), i.e. the minimum time for α to travel from β and arrive at γ. Equation (2) models an assignment problem (AP) that determines which vehicles are dispatched to a demand point to meet the demand. Note that the objective of Equation (2) is not to minimize the summation of arrival times or the latest arrival time. Since Equation (1) identifies the earliest arrival times, Equation (2) intends to minimize total travel time. In many VRPs, less travel time implies less uncertainty and more reliable vehicle routing and assignment (Chen et al., 2013). Both Equations (1) and (2) are pure integer programming problems. If (2) is a balanced transportation problem and both constrains may be changed to equality constraints. If Subject to:

Solution steps
Optimization software packages may be used to solve models in Equations (1) and (2). For example, the General Algebraic Modeling System (GAMS; GAMS Development Corporation, 2013) is a high-level modeling system for mathematical programming and optimization. Equations (1) and (2) may be described in algebraic statements in GAMS input files and solvers may be used to find optimal solutions and values. For medium to large road networks, however, this approach is ineffective or inefficient due to the limitation of computer memory and the inconvenience of extremely long Y , = 1 travels to 0 otherwise computation time. For real-time transportation planning, time-and space-efficient algorithms must be developed to solve the models in Equations (1) and (2). Figure 1 shows three steps to identify optimal solutions and values to the TDVRP.

TDSP algorithms
Real-time TDVRPs remain a great challenge due to time and space complexities. In the TDVRP, characteristics of the road network change with time; a shortest path computed from a snapshot of the road network may not be the shortest path at a different time. The TDSP algorithm is developed to find the shortest paths between depots and demand points in real-time. Three assumptions related to the TDSP algorithm are: (a) The road network satisfies the FIFO principle, which specifies that if two vehicles take the same route from the same depot to the same demand point, the vehicle leaving the depot early always arrives at the demand point early. According to the FIFO principle, a vehicle should leave its depot or other nodes whenever it is ready. Waiting at any node is never beneficial because a vehicle that leaves later always arrives later; (b) The planning period is "discretized" into sufficiently small time intervals, δ's, = 1, 2, 3, … and ∈ Δ where Δ is the set of time intervals over the planning period; and (c) Travel time between two nodes connected by an arc depends on the time at which a vehicle leaves the beginning node of the arc.
The TDSP algorithm with node labeling was first developed by Dreyfus (1969) and further validated by Kaufman and Smith (1993). The TDSP algorithm with node labeling listed below finds the earliest arrival time of a vehicle at a demand point given that the vehicle is stationed at a depot when travel begins. The algorithm is executed for each pair of demand point and depot for each vehicle. Vehicles stationed at the same depot may be ready to travel from the depot at different times because different vehicles may require different preparation times, e.g. time for a driver to arrive at the depot and get ready may vary. The earliest arrival times of vehicles, which are stationed at the same depot, at the same demand point may be different and need to be calculated separately using the TDSP algorithm.

TDSP algorithm with node labeling
(1) Assign to every node i, i ∈ V, in a transportation network a value representing the arrival time of a vehicle α, ∈ A, at the node. For a depot node β where α is stationed, set the value to a finite positive integer number representing the time at which α is ready to travel from β. Set the value to infinity for all other nodes i's, i ∈ V and i ≠ β; (2) Mark β visited. Mark all other i's unvisited. Set β as the current node; (3) Calculate the arrival time of α at each unvisited neighbor j, j ∈ V, of the current node i, i ∈ V. A node j is a neighbor of i if there is an arc begins at i and points at j, i.e. (i, j) ∈ E. The arrival Step 1: Develop and implement the TDSP algorithm to identify the shortest paths between pairs of depots and demand points Step 2: Apply the assignment algorithm to output of Step 1 to minimize total travel time Step 3: Validate and present optimal solutions and values to the TDVRP time at j is the summation of the value set for i and the travel time c i,j, between i and j. The travel time c i,j, is obtained from a three-dimensional matrix, |V| × |V| × |Δ|, which stores time-dependent travel times. For example, if i = 1, j = 2, and the value set for node 1 is 15, the travel time between nodes 1 and 2, c 1,2,15 , is a component in the matrix identified by node 1, node 2, and δ = 15. δ = 15 indicates α begins traveling from node 1 at time 15; (4) For each unvisited neighbor j of the current node i, compare the arrival time at j calculated in Step 3 and the value set for j. Set the value for j as the smaller one between the two; (5) Identify the unvisited node with the smallest value. Mark the node visited. Set the node as the current node. If the node is the desired demand point γ, ∈ Γ, stop. The value set for γ is the earliest arrival time of α traveling from β at γ. Otherwise go to Step 3.
The TDSP algorithm with node labeling is implemented in VBA ® and the GAMS. The algorithm can only solve TDVRPs for small road networks with less than 500 nodes when there are more than 1,000 time intervals using a Windows 7 ×64 PC, Intel i7-3770 CPU @3.40 GHZ, and 16.0 GB RAM. The main reason for size limitation on the road network is the large storage space required by node labeling.
Each component in the matrix is travel time from one node to the other during a time interval. These travel times are often obtained through field observations (Rakha, El-Shawarby, Arafeh, & Dion, 2006) and real-time monitoring of road network conditions. |V| is the size of the road network, i.e. the number of nodes. |Δ| is the number of time intervals. For example, if the road network size increases by tenfold, the storage space for the three-dimensional matrix increases by 100 times.
The degree of a node in a network is the number of arcs connected to the node. Most real-world road networks have a mean degree between two and four (Barabasi, 2002;Jeong, 2003). On average, each node is connected to two to four arcs. If arc labeling is used, the TDSP algorithm manipulates a two-dimensional matrix, |E| × |Δ|, where |E| is the number of arcs in the road network. Each component in the two-dimensional matrix is travel time along an arc during a time interval. Travel times in the two-dimensional matrix are the same as those in the three-dimensional matrix, but are organized in a different format that reduces storage space requirement. Since |E| ≤ 4 × |V| according to the mean degree of a road network, storage space requirement for arc labeling is at most 4 × |V| × |Δ|, which is much less than |V| × |V| × |Δ| required for node labeling for large road networks. The TDSP algorithm with arc labeling described below is used to identify the earliest arrival times of vehicles at demand points. This is a significant improvement for the TDSP algorithm and greatly increases space efficiency of the algorithm.

TDSP algorithm with arc labeling
(1) Assign to every arc (i, j), (i, j) ∈ E, in a transportation network a value representing the arrival time of a vehicle α, ∈ A, at j, i, j ∈ V. For (i, j) in which i = β, a depot node where α is stationed, set the value to a finite positive integer number representing the time at which α arrives at j. The arrival time at j is the summation of time at which α is ready to travel from β and travel time from β to j. The travel time from β to j, c ,j, , is obtained from a two-dimensional matrix, |E| × |Δ|, which stores time-dependent travel times. For example, if β = 1, j = 2, and the time at which α is ready to travel from β is δ = 15, the travel time between nodes 1 and 2, c 1,2,15 , is a component in the matrix identified by arc 1, 2 and δ = 15. Set the value to infinity for all other arcs (i, j), (i, j) ∈ E and i ≠ β; (2) Mark all (i, j) unvisited; (3) Identify the unvisited (i, j) with the smallest value. Mark (i, j) visited. Set the destination node, i.e. the second node j, in (i, j) as the current node. If j is the desired demand point γ, ∈ Γ, stop. The value set for (i, ) is the earliest arrival time of α traveling from β at γ; (4) For each unvisited (i, j) whose i is the current node, compare the arrival time at j and the value set for (i, j). time from i to j, c i,j, , is obtained from the matrix |E| × |Δ|. δ is the value set for the arc marked as visited in Step 3; (5) For each unvisited (i, j) whose i is the current node, set its value as the smaller one between the arrival time at j calculated in Step 4 and the value set for (i, j). Go to Step 3.

Assignment algorithm
After calculating all the shortest paths between depots and demand points, an assignment algorithm needs to be implemented to determine which vehicles from a depot will travel to a demand point to meet the demand. There are several variants of APs (Pentico, 2007), e.g. generalized AP, bottleneck AP, quadratic AP, and semi-AP. A semi-AP is a modified classic AP which differs in the disproportionality between demand points and depots. The TDVRP is a semi-AP except that a vehicle is not required to be assigned to a demand point. The Hungarian method (Kuhn, 1955) is used as the assignment algorithm to solve the TDVRP.

Case study of a small road network
A road network in Midwest of the United States of America ( Figure 2) is analyzed to validate the TDSP and assignment algorithms. In this case study, a no-notice evacuation prompts the need for transportation agencies to deploy vehicles and personnel to specific locations for traffic control. The network has 346 nodes (diamonds in Figure 2), out of which six are depots (red diamonds) and 15 are demand points (blue diamonds). The depots are locations where a state or federal transportation agency, e.g. the state department of transportation, dispatches service vehicles and equipment. The demand points are locations where the agency identifies a need for traffic control during a no-notice evacuation event. Overall, there are 654 roads (arcs in Figure 2) connecting all the nodes; thus the All computation results are obtained using a Windows 7 ×64 PC, as described previously. Table 1 summarizes run times.
Algorithms implemented in VBA ® and the GAMS provide the same optimal solutions and values. Table 1 shows, however, that there is a substantial advantage in using VBA ® to identify the optimal solution and value for real-time TDVRPs. The GAMS employs a suite of solvers, e.g. CPLEX, and algorithms, e.g. simplex and branch-and-cut algorithms, to solve linear programming and mixed integer programming problems, but lack efficiency in computer memory management. To compute the results in Table 1, both the TDSP algorithm with arc labeling and assignment algorithm are coded in the GAMS to identify the optimal solution and value. It takes even more computation time when models in Equations (1) and (2) are directly used in the GAMS to find the optimal solution and value. The TDSP and assignment algorithms implemented in VBA ® directly manipulate data and are customized to solve the TDVRPs; they are more efficient in terms of run time and computer storage space. Table 1 does not show a substantial difference in run times between node labeling and arc labeling for the TDSP algorithm. The difference between these two methods needs to be further investigated using larger road networks. Figure 3 shows a road network in the City of San Francisco (Brinkhoff, 2002), which includes 174,956 nodes and 223,001 arcs with a mean degree of 2.55 (= 223,001 × 2 174,956 ). A portion of the network (area inside the red box at the lower left corner of Figure 3), including 492 nodes and 559 arcs that connect the nodes, is selected to validate the TDSP and assignment algorithms. The selected portion has a mean degree of 2.27 (= 559 × 2 492 ). Six depots and 15 demand points are randomly and uniformly selected among the 492 nodes. Total 92 vehicles are randomly and uniformly stationed at six depots. Each time interval is one minute. The GAMS requires more than 10 h computing the optimal solution to the TDVRP for this network and sometimes stops unexpectedly due to insufficient computer memory. Run times for VBA ® are summarized in Table 2, which does not include run times for the GAMS since they are too large and the authors choose not to complete the runs in the GAMS. The authors choose five different planning periods, each of which specifies a time period during which vehicles must be dispatched to demand points to meet the demand. Run time increases as the planning period increases, indicating more computation time is needed for a larger planning period.

Node labeling and arc labeling in VBA ® for a large road network
Arc labeling takes more time to compute the optimal solution than node labeling when the planning period is less than 1,440 min. The difference in data processing time between arc labeling and node labeling is relatively insignificant when the planning period is short. Once the TDSP is identified, the optimal solution, i.e. a sequence of nodes on the shortest path, must be determined by backtracking from the destination to origin. Arc labeling is expected to take more time in backtracking since it must find the two end nodes of a visited arc. The advantage of arc labeling becomes evident when the number of nodes or planning period increases. In Table 2, when the planning period is 1, 440 min (24 h), node labeling could not compute the optimal solution due to insufficient memory;

Run times of TDSP and assignment algorithms for simulated road networks
The road network in Section 4.1 has a mean degree of 3.78 whereas the two road networks in Section 4.2 have a mean degree of 2.55 and 2.27, respectively. To further validate the algorithms implemented in VBA ® , simulated road networks are generated for a given degree distribution. A random network, one of the most studied networks, is described using a degree distribution where n is the total number of nodes, d is node degree, and p is the probability that an arc between two nodes exists. The mean degree d = (n − 1)p. Random networks were used to describe the structure of the US highway system (Barabasi, 2002;Jeong, 2003). Five random road networks with total 500, 1,000, 1,200, 1,500, and 2,000 nodes are generated using a mean degree of 2.5.
The TDSP and assignment algorithms in VBA ® are applied to the five simulated networks and their run times are summarized in Table 3, which reveals several important observations: (a) Run times increase as the total number of nodes increases; (b) Run times increase in most cases as the planning period increases. This is consistent with results in Tables 1 and 2; (c) If node labeling can compute an optimal solution to the TDVRP, it takes less time than arc labeling. This is due to the additional time required for arc labeling to backtrack from the destination to origin and identify the optimal solution; (d) Node labeling cannot compute an optimal solution when size of the matrix that stores travel times between nodes becomes too large. For instance, when total number of nodes is 1,500 and 2,000, node labeling cannot compute an optimal solution for any planning period; and  (e) Arc labeling can efficiently compute optimal solutions to the TDVRP for large road networks over a long planning period. The largest run time in Table 3 is 479 s, which is about 8 min and acceptable for real-time TDVRPs.

Conclusions and future research
This paper develops a new TDSP algorithm with arc labeling and integrates it with the assignment algorithm to identify the optimal solution and value for the TDVRP in real time. The experiments show that the algorithms implemented in VBA ® may be used for real-time vehicle routing in large road networks. The case study presented in Section 4.1 implements the algorithms in a no-notice evacuation for a metropolitan area with acceptable performance and does not require above average computing power. Both the TDSP and assignment algorithms are validated using a large road network and five simulated networks. The TDSP with arc labeling may be used for real-time TDVRP over a long planning period. Future research will investigate the performance of the TDSP and assignment algorithms in a variety of real-world road networks with various degree distributions, and explore certain heuristic algorithms, e.g. genetic algorithm, to solve the TDVRP. In addition, other programming languages, e.g. C++, Java, and Julie, will be used to further improve the efficiency of the algorithms.

Funding
The authors received no direct funding for this research.