# Delivery Vehicle Routing¶

Consider the following problem:

A delivery person is given a list of dropoff locations for his/her daily route. Our goal then is to reach all of our target destinations in as short of a distance as possible. Alternatively, we could also consider other objectives to optimize for, either separately or at the same time:

• shortest time

• number of left turns

• avoiding highways

• eco-friendly routes

• avoid residential areas

```import pandas as pd
import math
import numpy as np
import random
import time
import matplotlib.pyplot as plt
import osmnx as ox
import networkx as nx
from tqdm.notebook import tqdm
from smart_mobility_utilities.poi import poi, drawPOIS, drawRouteOrder
```
```# You can use different excel data
file_path = '../../data/rideshare_GA/rideshare_dropoff.csv'
numberOfStops = len(df)
```
```POIS = []
for _ in df.itertuples():
POIS
```
```[Name: Bay Street ID: 1601697722,
Name: 653 ID: 76112568,
Name: Castle Frank ID: 43804186,
Name: 719-721 ID: 363349504,
Name: Coxwell Avenue ID: 420539103,
Name: Don Mills Road ID: 31653252,
Name: 3945 ID: 76130372,
Name: Fit for Life ID: 1561890202,
Name: Dupont ID: 5320349317,
Name: Finch ID: 561044357,
Name: 1380 ID: 66246309,
Name: Scarborough Southwest ID: 685674284,
Name: Donlands Avenue ID: 420539141,
Name: Kipling Station ID: 32052148,
Name: 1870 ID: 67739333,
```
```drawPOIS(POIS, zoom=10)
```
```G = nx.DiGraph()
```
```# this will take 3-4 minutes because you need to find
# two routes between each pair of POIs
durations = []
pbar = tqdm(total=len(G.nodes())**2)
for source in G.nodes():
duration = []
for destination in G.nodes():
pbar.update()
if source == destination:
duration.append(0)
continue
route = source.route_to(destination)
duration.append(route['duration'])
durations.append(duration)
```
```trip = np.array(durations)
trip
```
```array([[   0. , 1244. ,  211.3,  270.4,  649.3, 1163.4, 1116.8,  192.9,
224.3, 1383.1,  739.7, 1234.2,  511.4, 1254.3,  736.8,  525.9],
[1230.4,    0. , 1019.3, 1268.5, 1140.8,  230.1,  610.8, 1349.8,
1221.9,  434.5,  701.1,  894.6, 1075.3, 1596.3, 1343.4,  872.3],
[ 211.1, 1128.5,    0. ,  481.5,  438. ,  952.1, 1283.8,  337.4,
399. , 1448.5,  907.5, 1022.9,  300.1, 1325.2,  947.9,  641.9],
[ 270.4, 1205.9,  481.7,    0. ,  919.7, 1406.4,  972.6,  463.3,
223.3, 1394.7,  553.8, 1504.6,  781.8, 1156.2,  466.4,  693.1],
[ 649.1, 1113.3,  438. ,  919.5,    0. ,  936.9, 1510.8,  764.9,
837. , 1433.3, 1241.2,  848.4,  137.9, 1607.9, 1385.2, 1031.8],
[1089.6,  229.9,  878.5, 1360. , 1000. ,    0. ,  761.3, 1209. ,
1277.5,  635.2,  834.5,  827.8,  934.5, 1732.6, 1479.7, 1008.6],
[1116.2,  628.2, 1276.3,  971.6, 1528.3,  841.2,    0. , 1309.1,
927.1,  620.3,  404.2, 1282.1, 1462.8, 1368.3, 1087.9,  811.6],
[ 192.9, 1318.7,  340.6,  463.3,  765. , 1142.3, 1292.9,    0. ,
400.4, 1551.9,  915.8, 1213.1,  640.7, 1125.2,  796.9,  694.7],
[ 219.9, 1139.1,  398.8,  226.6,  836.8, 1316.6,  928.1,  406.5,
0. , 1327.9,  531.7, 1387.4,  698.9, 1176.8,  620.5,  499.7],
[1334.4,  450.3, 1347.9, 1370.8, 1469.4,  621.5,  562.8, 1524.2,
1324.2,    0. ,  803.4, 1244.8, 1403.9, 1698.6, 1445.7,  952. ],
[ 720.3,  652.1,  886.9,  567.4, 1238.5,  852.6,  418.8,  913.2,
522.9,  840.9,    0. , 1293.4, 1107.5, 1121.3,  683.7,  422.2],
[1279. ,  994.4, 1067.9, 1549.4,  918.3,  912.7, 1391.7, 1397.7,
1466.9, 1314.4, 1464.9,    0. , 1056.2, 2145.9, 1923.2, 1263.3],
[ 511.2, 1050.3,  300.1,  781.6,  137.9,  873.9, 1447.8,  637.5,
699.1, 1370.3, 1108.6,  944.7,    0. , 1572. , 1248. ,  899.2],
[1300.3, 1439.4, 1370.6, 1147.4, 1662.6, 1652.4, 1295.4, 1174. ,
1178.6, 1628.2, 1118.9, 2036.9, 1621.1,    0. ,  681. , 1432.7],
[ 736.8, 1322.6,  948.1,  466.4, 1386.1, 1535.6, 1106.6,  793.7,
619.2, 1511.4,  687.8, 1782.1, 1248.2,  689.8,    0. ,  977.9],
[ 509.3,  854.4,  651.1,  740.5, 1022.6, 1067.4,  801.6,  699.1,
530.2,  990.5,  425.3, 1158.4,  891.6, 1485. , 1032.7,    0. ]])
```
```# define a function that can compute the fitness value of each solution in the population
def cal_pop_objective_1(m, pop):
M = np.zeros([8,numberOfStops-1])
fitness = np.zeros(8)
for i in range(8):
for j in range(numberOfStops-1):
M[i,j] = m[pop[i,j]-1,pop[i,j+1]-1]
for k in range(8):
fitness[k] = np.sum(M[k])
return fitness
```
```# define a function that can select the best individuals in the current generation as the parents to produce the offsprings
def select_mating(pop, fitness, num_parents):
parents = np.zeros([num_parents, pop.shape[1]])
for parent_num in range(num_parents):
min_fitness_index = np.where(fitness == np.min(fitness))
min_fitness_index = min_fitness_index[0][0]
parents[parent_num, :] = pop[min_fitness_index, :]
fitness[min_fitness_index] = 999999999
return parents
```
```# implementation of order 1 crossover
size = len(mum)

# select random start/end position for crossover
alice, bob = [-1] * size, [-1] * size
start, end = sorted([random.randrange(size) for _ in range(2)])

# replicate mum's sequence for alice
alice_inherited = []
for i in range(start, end + 1):
alice[i] = mum[i]
alice_inherited.append(mum[i])

fixed_pos = list(range(start, end + 1))
i = 0
while i < size:
if i in fixed_pos:
i += 1
continue
test_alice = alice[i]
if test_alice==-1:
i += 1
return alice
```
```# implementation of swap mutation
def mutation(offspring_crossover):
a = random.randint(0,numberOfStops-1)
b = random.randint(0,numberOfStops-1)
for i in range(offspring_crossover.shape[0]):
mut_1 = offspring_crossover[i,a]
mut_2 = offspring_crossover[i,b]
offspring_crossover[i,a] = mut_2
offspring_crossover[i,b] = mut_1
return offspring_crossover
```
```"implementation of genetic algorithm on ride sharing problem"

# ga parameters
sol_per_pop = 8
num_parents_mating = 4
num_offsprings = 4

# define population size and offspring size
pop_size = (sol_per_pop, trip.shape[1])
offspring_size = (num_offsprings, trip.shape[1])

# generate the initial population
new_population = np.zeros(pop_size)
for i in range(8):
new_population[i] = np.arange(1,numberOfStops+1)
random.shuffle(new_population[i])
new_population = new_population.astype(int)

best_outputs = []
num_of_generation = 30000

running_time = np.zeros(num_of_generation)

for generation in tqdm(range(num_of_generation), total=num_of_generation):

# record the starting time
start_time = time.time()

# compute the fitness of each individual in the population
fitness = cal_pop_objective_1(trip, new_population)

# record the best fitness in the current generation
best_outputs.append(np.min(fitness))

# select the best 4 individuals in the population as parents for mating
parents = select_mating(new_population, fitness, num_parents_mating)

# generate offsprings using crossover
offspring_crossover = np.zeros(offspring_size)
r_1 = np.random.random()
if r_1 < 0.7: # crossover probability = 0.7
offspring_crossover[0] = crossover(parents[0], parents[1])
offspring_crossover[1] = crossover(parents[1], parents[0])
offspring_crossover[2] = crossover(parents[2], parents[3])
offspring_crossover[3] = crossover(parents[3], parents[2])
else:
offspring_crossover[0] = parents[0]
offspring_crossover[1] = parents[1]
offspring_crossover[2] = parents[2]
offspring_crossover[3] = parents[3]

# offsprings mutation
r_2 = np.random.random()
if r_2 < 0.05: # mutation probability = 0.05
offspring_mutation = mutation(offspring_crossover)
else:
offspring_mutation = offspring_crossover

# generate the new population based on parents and offsprings
new_population[0:parents.shape[0], :] = parents
new_population[parents.shape[0]:, :] = offspring_mutation

# compute and store the running time
running_time[generation] = time.time() - start_time
```
```# get the best solution after iterating all the generations
# compute the fitness of each individual in the final generation
fitness = cal_pop_objective_1(trip, new_population)
# return the index of the solution corresponding to the best fitness
best_index = np.where(fitness == np.min(fitness))

running_time = running_time.cumsum()

# to get the index of iteration when the fitness reaches its minimum value (convergence point of the algorithm)
convergence_index = np.where(best_outputs == np.min(best_outputs))
convergence_time = running_time[convergence_index]

print("Best solution:", new_population[best_index, :][0][0])
print("Best fitness:", fitness[best_index][0])
print("The algorithm converges in", np.min(convergence_index), "iterations.")
print("The algorithm converges in", np.min(convergence_time), "seconds.")
```
```Best solution: [10  2 12  5 13  3  6  7 11 16 14 15  4  9  1  8]
Best fitness: 8509.5
The algorithm converges in 5246 iterations.
The algorithm converges in 0.838292121887207 seconds.
```
```fig, ax = plt.subplots(2, figsize = (8,numberOfStops-2))
ax[0].plot(best_outputs, 'green')
ax[0].set_xlabel("Iteration", fontsize = 10)
ax[0].set_ylabel("Fitness", fontsize = 10)
ax[0].set_title("GA (Fitness vs # of iterations)   1 station", fontsize = 12)

ax[1].plot(running_time, best_outputs, 'red')
ax[1].set_xlabel("Running time/s", fontsize = 10)
ax[1].set_ylabel("Fitness", fontsize = 10)
ax[1].set_title("GA (Fitness vs Running time)   1 station", fontsize = 12)
plt.show()
```
```#Visualization

optimalRoute = new_population[best_index, :]
optimalRoute = optimalRoute[0][0]
currentRoute = G[POIS[optimalRoute[0]]][POIS[optimalRoute[1]]]['route']['coords']
for i in range(1,len(optimalRoute)-1):
currentRoute += G[POIS[optimalRoute[i]-1]][POIS[optimalRoute[i+1]-1]]['route']['coords']
drawRouteOrder(route = currentRoute, POIS = list(G.nodes), order= list(new_population[best_index][0]))
```
Make this Notebook Trusted to load map: File -> Trust Notebook