CVRP
Lộc
2024-01-24
In this section here, we will learn about the concepts of Vehicle Routing Problem in Supply Chain and how to implement in R.
The dataset you just need 3 main columns: longitude, latitude and demand.
Input:
Data:
#Call packages:
pacman::p_load(sf,
rio,
here,
janitor,
tidyverse,
dplyr,
magrittr,
ggplot2,
purrr,
lubridate,
mice,
plotly)
optimize<-read.csv(r"(C:\\Users\\locca\\Documents\\Xuân Lộc\\VILAS\\Final project\\Optimize_df.csv)")
#New manufacter:
new_manufacter= data.frame(Lat =c(21.12256201),
Long = c(105.9150683))
route<-data.frame("Manufacter",
new_manufacter[,c("Long","Lat")],
0)
colnames(route)<-c("Customers",
"Longitude",
"Latitude",
"Total.transactions")
route<-rbind(route,
optimize%>% select(Customers,
Longitude,
Latitude,
Total.transactions))
colnames(route)[4]<-"demand"
route$node<-1:nrow(route)There are functions to solve the VRP that I learned from Github manufactering by Argaadya author.
I tested all so you confidently copy and paste to your script. Quick explanation such as:
giant_path function is for calculate the least distance between multiple points based on the Hamiltonian Cycle algorithm.
assign_vehicle function is for separating the output from giant_path to groups segment appropriate to the conditions such as: the limited capacity of vehicle, the limited fuel,…
fitness is for computation base on the Genetic algorithms which is stochastic search algorithms inspired by the basic principles of biological evolution and natural selection.
Functions:
# Giant_path function for the Hamiltonian Cycle:
giant_path <- function(distance, demand){
visited_spot <- NULL
vehicle_num <- 1
post <- 1:ncol(distance)
names(demand) <- 1:length(demand)
# Randomly select initial spot
initial_spot <- sample(2:length(demand), 1)
while (any(demand != 0)) {
available_spot <- which(demand != 0)
# Calculate the distance to unvisited spot
initial_dist <- distance[ c(available_spot), initial_spot]
initial_dist <- initial_dist[ which(names(initial_dist) != initial_spot)]
visited_spot <- c(visited_spot, initial_spot)
demand[ initial_spot ] <- 0
if ( length(initial_dist)>1) {
initial_spot <- which(initial_dist == min(initial_dist)) %>% names() %>% as.numeric()
} else {
initial_spot <- which(demand != 0)
}
}
visited_spot <- c(1, visited_spot, 1)
names(visited_spot) <- NULL
total_distance <- embed(visited_spot, 2)[ , 2:1] %>% distance[.] %>% sum()
result <- list(route = visited_spot,
total_distance = total_distance)
return(result)
}
# Assign Vehicle function to assign the route:
assign_vehicle <- function(x, demand, capacity, distance){
vehicle_load <- capacity
visited_spot <- NULL
vehicle_num <- 1
for (i in x) {
initial_spot <- i
if (vehicle_load >= demand[i]) {
visited_spot <- c(visited_spot, initial_spot)
vehicle_load <- vehicle_load - demand[ initial_spot ]
} else {
# Go back to depot
vehicle_num <- vehicle_num + 1
vehicle_load <- capacity
visited_spot <- c(visited_spot, 1)
# Revisit the spot
visited_spot <- c(visited_spot, initial_spot)
vehicle_load <- vehicle_load - demand[ initial_spot ]
}
}
visited_spot <- c(visited_spot)
total_distance <- embed(visited_spot, 2)[ , 2:1] %>% distance[.] %>% sum()
return(list(distance = total_distance,
route = visited_spot,
num_vehicle = vehicle_num))
}
#Fitness function for ga method:
fitness <- function(x, capacity, demand, distance, ...){
vehicle_load <- capacity
visited_spot <- 1
vehicle_num <- 1
for (i in x) {
initial_spot <- i
if (vehicle_load > demand[initial_spot]) {
# Go to the spot
visited_spot <- c(visited_spot, initial_spot)
vehicle_load <- vehicle_load - demand[ initial_spot ]
} else {
# Go back to depot
vehicle_load <- capacity
visited_spot <- c(visited_spot, 1)
vehicle_num <- vehicle_num + 1
# Go to the spot
visited_spot <- c(visited_spot, initial_spot)
vehicle_load <- vehicle_load - demand[ initial_spot ]
}
}
visited_spot <- c(visited_spot, 1)
total_distance <- embed(visited_spot, 2)[ , 2:1] %>% distance[.] %>% sum()
return(-total_distance)
}
fitness_explain <- function(x, capacity, demand, distance, ...){
vehicle_load <- capacity
visited_spot <- 1
vehicle_num <- 1
total_demand <- NULL
for (i in x) {
initial_spot <- i
if (vehicle_load > demand[initial_spot]) {
# Go to the spot
visited_spot <- c(visited_spot, initial_spot)
vehicle_load <- vehicle_load - demand[ initial_spot ]
} else {
total_demand <- c(total_demand, 6000 - vehicle_load)
# Go back to depot
vehicle_load <- capacity
visited_spot <- c(visited_spot, 1)
vehicle_num <- vehicle_num + 1
# Go to the spot
visited_spot <- c(visited_spot, initial_spot)
vehicle_load <- vehicle_load - demand[ initial_spot ]
}
}
total_demand <- c(total_demand, 6000 - vehicle_load)
visited_spot <- c(visited_spot, 1)
total_distance <- embed(visited_spot, 2)[ , 2:1] %>% distance[.] %>% sum()
result <- list(route = visited_spot,
total_distance = total_distance,
vehicle_num = vehicle_num,
total_demand = total_demand)
return(result)
}Output:
Nearest Neighbour Algorithm:
Let say the limited capacity of vehicle in one run is 3000 transactions.
The output show that the output in list form:
$total_distance is 18.39416 meaning that if we don’t need to change our travel route and just go to all the points
$distance is 22.77717 higher that 18.39416 because we added the assumption that the limited capacity of vehicle is 3000 transactions.
$num_vehicle is 5 meaning 5 is suitable route number.
#Create distance matrix for route:
route_distance <- dist(route[ ,c("Longitude","Latitude")],
upper = T,
diag = T) %>%
as.matrix()
#
best_route <- list(total_distance = Inf)
for (i in 1:10^4) {
try_route <- giant_path(demand = route$demand,
distance = route_distance)
if (best_route$total_distance > try_route$total_distance) {
best_route <- try_route
}
}
best_route## $route
## [1] 1 18 36 20 34 37 42 50 24 22 8 35 23 41 3 2 40 10 31 9 11 49 19 28 38
## [26] 25 4 45 33 21 5 47 17 15 46 29 32 13 27 16 7 6 44 14 26 43 12 48 39 51
## [51] 30 1
##
## $total_distance
## [1] 18.39416nearest_result <- assign_vehicle(x = best_route$route,
demand = route$demand,
capacity = 3000,
distance = route_distance)
nearest_result## $distance
## [1] 22.77717
##
## $route
## [1] 1 18 36 20 34 37 42 50 24 22 1 8 35 23 41 3 2 40 10 31 1 9 11 49 19
## [26] 28 38 25 4 45 33 21 5 1 47 17 15 46 29 32 13 27 16 1 7 6 44 14 26 43
## [51] 12 48 39 51 30 1
##
## $num_vehicle
## [1] 5Now we will plot the output. I will present by 2 way:
ggplot package but I suggest if there are little numbers of points. You can see that the plot is not clear and coherent.
leaflet package is the best approach if you don’t have to show in a static version.
#ggplot:
library(ggrepel)
nearest_result$route %>%
as.data.frame() %>%
rename(from = ".") %>%
mutate(to = lead(from)) %>%
slice(- nrow(.)) %>%
left_join(route,
by = c("from" = "node")) %>%
rename(from_x1 = Longitude,
from_x2 = Latitude) %>%
left_join(route,
by = c("to" = "node")) %>%
rename(to_x1 = Longitude,
to_x2 = Latitude) %>%
rename(demand = demand.y) %>%
select(-demand.x) %>%
mutate(unit = c( rep("A", 10),
rep("B",10),
rep("C", 13),
rep("D", 10),
rep("E",12))) %>%
ggplot() +
geom_point(data = route,
aes(Longitude,
Latitude),
show.legend = F) +
geom_curve(aes(x = from_x1,
xend = to_x1,
y = from_x2,
yend = to_x2,
color = unit),
curvature = -0.05,
size = 0.5,
alpha = 0.7,
arrow = arrow(type = "closed",
angle = 30,
length = unit(3, "mm")) ) +
geom_point(data = route[1,],
aes(Longitude,
Latitude),
color = "red",
size = 3) +
geom_label_repel(data = route,
aes(Longitude,
Latitude,
label = node),
alpha = 0.7,
max.overlaps = Inf) +
scale_color_manual(values = c("firebrick",
"orange",
"dodgerblue",
"green3",
"purple"))+
theme_minimal() +
theme(panel.grid = element_blank()) +
theme(legend.position = "top") +
labs(x = quote(x[1]),
y = quote(x[2]),
title = "Vehicle Routing Problem with Nearest Neighbour",
subtitle = "Vehicle Capacity = 3000",
color = "Unit") ## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
Also, I use osrm package to calculate the closest value of distance in real life.
#leaflet:
df<-nearest_result$route %>%
as.data.frame() %>%
rename(from = ".") %>%
mutate(to = lead(from)) %>%
slice(- nrow(.)) %>%
left_join(route,
by = c("from" = "node")) %>%
mutate(unit = c( rep("A", 10),
rep("B",10),
rep("C", 13),
rep("D", 10),
rep("E", 12)))
# Calculate the shortest round trip between the facility locations
library(osrm)## Warning: package 'osrm' was built under R version 4.2.3## Data: (c) OpenStreetMap contributors, ODbL 1.0 - http://www.openstreetmap.org/copyright## Routing: OSRM - http://project-osrm.org/## Trip A:
dfA<-df %>% filter(unit == "A")
tripA <- osrmTrip(dfA[,c("Longitude","Latitude")],
overview = "full")
dfB<-df %>% filter(unit == "B")
tripB <- osrmTrip(dfB[,c("Longitude","Latitude")],
overview = "full")
dfC<-df %>% filter(unit == "C")
tripC <- osrmTrip(dfC[,c("Longitude","Latitude")],
overview = "full")
dfD<-df %>% filter(unit == "D")
tripD <- osrmTrip(dfD[,c("Longitude","Latitude")],
overview = "full")
dfE<-df %>% filter(unit == "E")
tripE <- osrmTrip(dfE[,c("Longitude","Latitude")],
overview = "full")
##Create matrix for summary:
summary = matrix(c(as.numeric(tripA[[1]]$summary),
as.numeric(tripB[[1]]$summary),
as.numeric(tripC[[1]]$summary),
as.numeric(tripD[[1]]$summary),
as.numeric(tripE[[1]]$summary)),
ncol = 2,
nrow = 5,
byrow = T)
colnames(summary) = c("Duration","Distance")
rownames(summary) = c("Trip A","Trip B","Trip C","Trip D","Trip E")
library(DT)
datatable(summary)#Create labels:
labels<- paste0("<strong> Customers </strong> ",
optimize$Customers, "<br/> ",
"<strong> Longitude: </strong> ",
optimize$Longitude, "<br/> ",
"<strong> Latitude: </strong> ",
optimize$Latitude, "<br/> ",
"<strong> Total Transactions </strong> ",
optimize$Total.transactions, "<br/> ") %>%
lapply(htmltools::HTML)
#Plot map with leaflet:
library(leaflet)
leaflet(data = optimize) %>%
addTiles() %>%
addMarkers(lng = optimize$Longitude,
lat = optimize$Latitude,
label = ~labels) %>%
addPolylines(data = tripA[[1]]$trip,
label = "Trip A",
color = "firebrick") %>%
addPolylines(data = tripB[[1]]$trip,
label = "Trip B",
color = "orange") %>%
addPolylines(data = tripC[[1]]$trip,
label = "Trip C",
color = "dodgerblue") %>%
addPolylines(data = tripD[[1]]$trip,
label = "Trip D",
color = "red") %>%
addPolylines(data = tripE[[1]]$trip,
label = "Trip E",
color = "purple")Genetic Algorithm:
R will take a little time to run so you just be patient to wait.
The output show that the final distance is 20.0774 comparing to the above output is better.(20.0774 < 22.77717)
suggestion_route <- matrix(ncol = (max(route$node)-1), nrow = 100)
for (i in 1:100) {
try_route <- giant_path(demand = route$demand, distance = route_distance)
suggestion_route[i, ] <- try_route$route[2:max(route$node)] %>% as.numeric()
}
head(suggestion_route)## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14]
## [1,] 23 8 35 24 50 37 34 42 22 20 18 36 29 46
## [2,] 14 26 17 2 40 3 41 35 8 23 24 50 37 34
## [3,] 38 25 28 19 49 11 31 9 10 40 2 3 41 35
## [4,] 28 38 25 49 19 4 45 33 21 5 11 31 9 10
## [5,] 18 36 20 34 37 42 50 24 22 8 35 23 41 3
## [6,] 27 13 32 46 29 16 36 18 35 8 23 24 50 37
## [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24] [,25] [,26]
## [1,] 32 13 27 16 41 3 2 40 10 31 9 11
## [2,] 42 22 20 18 36 29 46 32 13 27 16 15
## [3,] 8 23 24 50 37 34 42 22 20 18 36 29
## [4,] 40 2 3 41 35 8 23 24 50 37 34 42
## [5,] 2 40 10 31 9 11 49 19 28 38 25 4
## [6,] 34 42 22 20 41 3 2 40 10 31 9 11
## [,27] [,28] [,29] [,30] [,31] [,32] [,33] [,34] [,35] [,36] [,37] [,38]
## [1,] 49 19 28 38 25 4 45 33 21 5 47 17
## [2,] 44 6 7 51 39 48 12 43 10 31 9 11
## [3,] 46 32 13 27 16 15 17 44 14 26 43 6
## [4,] 22 20 18 36 29 46 32 13 27 16 15 17
## [5,] 45 33 21 5 47 17 15 46 29 32 13 27
## [6,] 49 19 28 38 25 4 45 33 21 5 47 17
## [,39] [,40] [,41] [,42] [,43] [,44] [,45] [,46] [,47] [,48] [,49] [,50]
## [1,] 15 44 14 26 43 6 7 51 39 48 12 30
## [2,] 49 19 28 38 25 4 45 33 21 5 47 30
## [3,] 7 51 39 48 12 4 45 33 21 5 47 30
## [4,] 44 14 26 43 6 7 51 39 48 12 30 47
## [5,] 16 7 6 44 14 26 43 12 48 39 51 30
## [6,] 15 44 14 26 43 6 7 51 39 48 12 30library(GA)## Warning: package 'GA' was built under R version 4.2.3## Loading required package: foreach##
## Attaching package: 'foreach'## The following objects are masked from 'package:purrr':
##
## accumulate, when## Loading required package: iterators## Package 'GA' version 3.2.3
## Type 'citation("GA")' for citing this R package in publications.##
## Attaching package: 'GA'## The following object is masked from 'package:utils':
##
## dega_route <- ga(type = "permutation",
fitness = fitness,
capacity = 3000,
demand = route$demand,
distance = route_distance,
lower = 2,
upper = max(route$node),
mutation = gaperm_swMutation,
popSize = 100,
pmutation = 0.1,
maxiter = 5 * 10^4,
suggestions = suggestion_route,
monitor = F,
seed = 123)
summary(ga_route)## ── Genetic Algorithm ───────────────────
##
## GA settings:
## Type = permutation
## Population size = 100
## Number of generations = 50000
## Elitism = 5
## Crossover probability = 0.8
## Mutation probability = 0.1
## Suggestions =
## x1 x2 x3 x4 x5 x6 x7 x8 x9 x10 ... x49 x50
## 1 23 8 35 24 50 37 34 42 22 20 12 30
## 2 14 26 17 2 40 3 41 35 8 23 47 30
## 3 38 25 28 19 49 11 31 9 10 40 47 30
## 4 28 38 25 49 19 4 45 33 21 5 30 47
## 5 18 36 20 34 37 42 50 24 22 8 51 30
## 6 27 13 32 46 29 16 36 18 35 8 12 30
## 7 31 9 10 40 2 3 41 35 8 23 47 30
## 8 38 25 28 19 49 11 31 9 10 40 47 30
## 9 3 41 35 8 23 24 50 37 34 42 12 30
## 10 8 35 23 24 50 37 34 42 22 20 12 30
## ...
## 99 41 3 2 40 10 31 9 11 49 19 51 30
## 100 42 37 34 22 20 50 24 23 8 35 12 30
##
## GA results:
## Iterations = 50000
## Fitness function value = -20.55624
## Solution =
## x1 x2 x3 x4 x5 x6 x7 x8 x9 x10 ... x49 x50
## [1,] 8 23 24 50 42 37 34 20 22 3 36 18plot(ga_route)