# -*- coding: utf-8 -*- """ Created on Fri Jan 22 14:08:13 2021 @author: Charlie Ferguson """ ### Packages for ORTools from __future__ import print_function from ortools.constraint_solver import routing_enums_pb2 from ortools.constraint_solver import pywrapcp #General Packages import pandas as pd #import os import numpy as np from multiprocessing import Pool from scipy.spatial import distance from pathlib import Path from itertools import repeat from collections import Counter #%% Define Type of Solver Solver="local" #or "server" if Solver=="local": FilePath="{Path to local server}" elif Solver=="server": FilePath="{Path to remote server}" del Solver #%% Define Variables Location="Nairobi" #Location Identifier Cap_Vehicle=50 #The number of filled containers able to be carried by the collection vehicle NumSetsArray=np.array([[10, 20, 30,50,70,100,250,500,1000,3000, 5000], #This defines the number of CBS users and the number of iterations for each [50, 50, 50,50,50,50 ,25, 25, 10, 10,10]]) CBSUsersFileSets=["CBSUsers_Num_10","CBSUsers_Num_20","CBSUsers_Num_30","CBSUsers_Num_50","CBSUsers_Num_70", #This defines the names of the folders (following from previous Python Script 1) "CBSUsers_Num_100","CBSUsers_Num_250","CBSUsers_Num_500","CBSUsers_Num_1000", "CBSUsers_Num_3000","CBSUsers_Num_5000"] #%% Define Functions def closest_node(node, nodes): closest_index = distance.cdist([node], nodes).argmin() return closest_index def closest_node_distance(node, nodes): closest_distance = min(min(distance.cdist([node], nodes,"euclidean"))) return closest_distance #%% Call preprocessed data #This is the transfer centre locations TransferCentreLocations=pd.read_csv(FilePath+"/InputFiles/"+Location+"_Inputs/OtherData/InternalTransferCentres_"+Location+".csv") TransferCentreLocations=TransferCentreLocations[["POINT_X","POINT_Y"]] TC_Coords = [tuple(x) for x in TransferCentreLocations.to_numpy()] #%% Computation of Strategy B (this is written to run in a HPC, can be converted in for loop for local machine use) def PrimaryFunction(i): #Determines how many folders to consider inputs from NumSets=NumSetsArray[1,i] NumUsers=int(CBSUsersFileSets[i].replace("CBSUsers_Num_","")) Num_Vehicle=int(np.ceil(NumUsers/Cap_Vehicle)) Runs_InfoList=[] FinalRouteDistances_Centre0=[] FinalRouteLoads_Centre0=[] FinalRouteDistances_Centre1=[] FinalRouteLoads_Centre1=[] FinalRouteDistances_Centre2=[] FinalRouteLoads_Centre2=[] FinalRouteDistances_Centre3=[] FinalRouteLoads_Centre3=[] FinalRouteNodes_Centre0=[] FinalRouteNodes_Centre1=[] FinalRouteNodes_Centre2=[] FinalRouteNodes_Centre3=[] for j in range(0,NumSets): #Bring in Distance Matrix created in earlier python scripts DistanceMatrix=pd.read_csv(FilePath+"/InputFiles/"+Location+"_Inputs/"+CBSUsersFileSets[i]+"/"+CBSUsersFileSets[i]+"_Set_"+str(j)+ "_DistanceMatrix.csv",header=None) DistanceMatrix_Info=pd.read_csv(FilePath+"/InputFiles/"+Location+"_Inputs/"+CBSUsersFileSets[i]+"/"+CBSUsersFileSets[i]+"_Set_"+str(j)+ "_DistanceMatrixComponents.csv") DM_Households=DistanceMatrix_Info.iloc[1,0] #Find locations of CBS user 'nearest points' CBSUserRoadPoints=pd.read_csv(FilePath+"/InputFiles/"+Location+"_Inputs/"+CBSUsersFileSets[i]+"/"+CBSUsersFileSets[i]+"_Set_"+str(j)+ "_PointsCSV.csv") CBSUserRoadPoints=CBSUserRoadPoints[["POINT_X","POINT_Y"]] RoadPoints_Cords = [tuple(x) for x in CBSUserRoadPoints.to_numpy()] #Find locations of depots TransferCentreIndexes=[RoadPoints_Cords.index(TC_Coords[x]) for x in range(0,len(TC_Coords))] CBSUserPoints=pd.read_csv(FilePath+"/InputFiles/"+Location+"_Inputs/"+CBSUsersFileSets[i]+"/"+CBSUsersFileSets[i]+"_Set_"+str(j)+ "_CBSPoints.csv") CBSUserPoints=CBSUserPoints[["POINT_X","POINT_Y"]] CBSUserPoints = [tuple(x) for x in CBSUserPoints.to_numpy()] CBSUserRoadPoints_Sample=CBSUserRoadPoints.iloc[1:DistanceMatrix_Info.iloc[1,0]+1,:] for t in range (0, 4): #For each of the different depots #t=0 my_file = Path(FilePath+"/InputFiles/"+Location+"_Inputs/"+CBSUsersFileSets[i]+"/"+CBSUsersFileSets[i]+"_Set_"+str(j)+ "_RouteOptimisationResultStrategyB_Centre"+str(t)+".csv") if my_file.is_file(): # The file exists pass else: DM_Indices=[TransferCentreIndexes[t]] +list(range(1,DM_Households+1)) # Import depot location as row 0 DistanceMatrix_Tailored=DistanceMatrix.iloc[DM_Indices, DM_Indices] ClosestDemands= [0 for i in range(0,DistanceMatrix_Info.iloc[1,0])] for x in range(0, len(CBSUserPoints)): ClosestDemands[closest_node(CBSUserPoints[x],CBSUserRoadPoints_Sample)-1]+=1 #This ensures that if the node load exceeds the vehicle capacity that multiple vehicles are dispatched to the same point ClosestDemands1=ClosestDemands ExtraVehicleControl=[] ExtraVehicleControl.append([0,len(ClosestDemands1)]) DistanceMatrix_Altered=DistanceMatrix_Tailored if max(ClosestDemands1)>Cap_Vehicle: for p in range(0, len(ClosestDemands1)): if ClosestDemands1[p] > Cap_Vehicle: NVisits=int(np.ceil(ClosestDemands1[p]/Cap_Vehicle)) StationRow=DistanceMatrix_Tailored.iloc[p,:] StationRows=pd.DataFrame([StationRow,]*(NVisits-1)) DistanceMatrix_Altered=DistanceMatrix_Altered.append(StationRows) ClosestDemands1[p] = ClosestDemands1[p] % Cap_Vehicle ClosestDemands1.extend(repeat(Cap_Vehicle, (NVisits-1))) ExtraVehicleControl.append([p,len(ClosestDemands1)]) ## ## DistanceMatrix_Altered=DistanceMatrix_Altered.reset_index(drop=True) DistanceMatrix_Altered.columns = range(DistanceMatrix_Altered.shape[1]) DistanceMatrix_AlteredT=DistanceMatrix_Altered.transpose() DistanceMatrix_Altered=DistanceMatrix_Altered.combine_first(DistanceMatrix_AlteredT) DistanceMatrix_Altered = DistanceMatrix_Altered.fillna(999999) ## DistanceMatrix_Tailored=DistanceMatrix_Altered.values np.fill_diagonal(DistanceMatrix_Tailored,0) ClosestDemands1=[0]+ClosestDemands1 Vehicle_Capacities=[Cap_Vehicle for i in range(Num_Vehicle)] ###### #Routing Exercise - THIS IS BASED ON GOOGLE OR-TOOLS CAPACITY-DRIVEN VEHICLE ROUTING PROBLEM (SEE FOR DETAILS/DOCUMENTATION) https://developers.google.com/optimization/routing/cvrp def create_data_model(): """Stores the data for the problem.""" data = {} data['distance_matrix'] = DistanceMatrix_Tailored data['num_vehicles'] = Num_Vehicle data['depot'] = 0 #This is the line with the depot point data['demands'] = ClosestDemands1 data['vehicle_capacities'] = Vehicle_Capacities return data global OutputRouteDistances global OutputRouteLoads global NodeRoutes NodeRoutes=[] OutputRouteDistances=[] OutputRouteLoads=[] def print_solution(data, manager, routing, solution): """Prints solution on console.""" global OutputRouteDistances global OutputRouteLoads global NodeRoutes NodeRoutes=[] OutputRouteDistances=[] OutputRouteLoads=[] total_distance = 0 total_load = 0 for vehicle_id in range(data['num_vehicles']): index = routing.Start(vehicle_id) plan_output = 'Route for vehicle {}:\n'.format(vehicle_id) route_distance = 0 route_load = 0 previous_index=0 while not routing.IsEnd(index): node_index = manager.IndexToNode(index) route_load += data['demands'][node_index] plan_output += ' {0} Load({1}) -> '.format(node_index, route_load) NodeRoutes.append([previous_index,index]) previous_index = index index = solution.Value(routing.NextVar(index)) route_distance += routing.GetArcCostForVehicle( previous_index, index, vehicle_id) plan_output += ' {0} Load({1})\n'.format(manager.IndexToNode(index), route_load) plan_output += 'Distance of the route: {}m\n'.format(route_distance) OutputRouteDistances.append(route_distance) plan_output += 'Load of the route: {}\n'.format(route_load) OutputRouteLoads.append(route_load) total_distance += route_distance total_load += route_load return OutputRouteDistances, OutputRouteLoads,NodeRoutes #These outputs (alterations from base Google OR Tools script) are stored for later evaluation def main(): data = create_data_model() manager = pywrapcp.RoutingIndexManager(len(data['distance_matrix']), data['num_vehicles'], data['depot']) routing = pywrapcp.RoutingModel(manager) def distance_callback(from_index, to_index): """Returns the distance between the two nodes.""" # Convert from routing variable Index to distance matrix NodeIndex. from_node = manager.IndexToNode(from_index) to_node = manager.IndexToNode(to_index) return data['distance_matrix'][from_node][to_node] transit_callback_index = routing.RegisterTransitCallback(distance_callback) routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index) def demand_callback(from_index): """Returns the demand of the node.""" # Convert from routing variable Index to demands NodeIndex. from_node = manager.IndexToNode(from_index) return data['demands'][from_node] demand_callback_index = routing.RegisterUnaryTransitCallback( demand_callback) routing.AddDimensionWithVehicleCapacity( demand_callback_index, 0, # null capacity slack data['vehicle_capacities'], # vehicle maximum capacities True, # start cumul to zero 'Capacity') # Setting first solution heuristic. search_parameters = pywrapcp.DefaultRoutingSearchParameters() search_parameters.first_solution_strategy = ( routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC) search_parameters.local_search_metaheuristic = ( routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH) search_parameters.time_limit.FromSeconds(1) search_parameters.time_limit.seconds = 300 # Solve the problem. solution = routing.SolveWithParameters(search_parameters) global routes def get_routes(solution, routing, manager): routes = [] for route_nbr in range(routing.vehicles()): index = routing.Start(route_nbr) route = [manager.IndexToNode(index)] while not routing.IsEnd(index): index = solution.Value(routing.NextVar(index)) route.append(manager.IndexToNode(index)) routes.append(route) return routes routes = get_routes(solution, routing, manager) if solution: print_solution(data, manager, routing, solution) if __name__ == '__main__': main() #END OF GOOGLE OR-TOOLS SCRIPT NodeRoutes1 = NodeRoutes for x in range(0, len(NodeRoutes1)): if NodeRoutes1[x][0]==0: NodeRoutes1[x][1]=0 for x in range(0, len(NodeRoutes1)): if NodeRoutes1[x]==[0,0]: NodeRoutes1[x+1][0]=0 for x in range(0, len(NodeRoutes1)): if NodeRoutes1[x]==[0,0]: NodeRoutes1[x][0]=NodeRoutes1[x-1][1] NodeRoutes1=NodeRoutes1[1:] NodeRoutes2=pd.DataFrame(NodeRoutes1) ExtraVehicleControlReversed=list(reversed(ExtraVehicleControl)) ExtraVehicleControlReversed=ExtraVehicleControlReversed[:-1] if NodeRoutes2.values.max()>len(CBSUserRoadPoints): for t in range(0,len(ExtraVehicleControl)): NodeRoutes2[NodeRoutes2>ExtraVehicleControl[t][1]] = ExtraVehicleControl[t][0] CBSUserRoadPoints_Corrected=CBSUserRoadPoints CBSUserRoadPoints_Corrected.loc[-1] = [0,0] # adding a row CBSUserRoadPoints_Corrected.index = CBSUserRoadPoints_Corrected.index + 1 # shifting index CBSUserRoadPoints_Corrected = CBSUserRoadPoints_Corrected.sort_index() NodeRoutes_Locations = pd.DataFrame(index=NodeRoutes2.index,columns=NodeRoutes2.columns) NodeRoutes_Locations=pd.concat([NodeRoutes_Locations,NodeRoutes_Locations],axis=1) NodeRoutes_Locations.columns=['Ax','Ay','Bx','By'] for x in range(0, len(NodeRoutes_Locations)): NodeRoutes_Locations.iloc[x,0]=CBSUserRoadPoints_Corrected .iloc[NodeRoutes2.iloc[x,0],0] NodeRoutes_Locations.iloc[x,1]=CBSUserRoadPoints_Corrected .iloc[NodeRoutes2.iloc[x,0],1] NodeRoutes_Locations.iloc[x,2]=CBSUserRoadPoints_Corrected .iloc[NodeRoutes2.iloc[x,1],0] NodeRoutes_Locations.iloc[x,3]=CBSUserRoadPoints_Corrected .iloc[NodeRoutes2.iloc[x,1],1] ### ## NodeRoutes_Locations_List=NodeRoutes_Locations.values.tolist() ################# #Compute the results depending on depot Runs_InfoList.append([j]+[t]+[Cap_Vehicle]+[Num_Vehicle]+routes) if t==0: FinalRouteDistances_Centre0.append(OutputRouteDistances) FinalRouteLoads_Centre0.append(OutputRouteLoads) FinalRouteNodes_Centre0.extend(NodeRoutes_Locations_List) elif t==1: FinalRouteDistances_Centre1.append(OutputRouteDistances) FinalRouteLoads_Centre1.append(OutputRouteLoads) FinalRouteNodes_Centre1.extend(NodeRoutes_Locations_List) elif t==2: FinalRouteDistances_Centre2.append(OutputRouteDistances) FinalRouteLoads_Centre2.append(OutputRouteLoads) FinalRouteNodes_Centre2.extend(NodeRoutes_Locations_List) elif t==3: FinalRouteDistances_Centre3.append(OutputRouteDistances) FinalRouteLoads_Centre3.append(OutputRouteLoads) FinalRouteNodes_Centre3.extend(NodeRoutes_Locations_List) Runs_Info_df=pd.DataFrame(Runs_InfoList) #Rename columns of route tables FinalRouteDistances_Centre0_df=pd.DataFrame(FinalRouteDistances_Centre0) FinalRouteDistances_Centre0_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_Centre0_df] FinalRouteDistances_Centre1_df=pd.DataFrame(FinalRouteDistances_Centre1) FinalRouteDistances_Centre1_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_Centre1_df] FinalRouteDistances_Centre2_df=pd.DataFrame(FinalRouteDistances_Centre2) FinalRouteDistances_Centre2_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_Centre2_df] FinalRouteDistances_Centre3_df=pd.DataFrame(FinalRouteDistances_Centre3) FinalRouteDistances_Centre3_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_Centre3_df] FinalRouteLoads_Centre0_df=pd.DataFrame(FinalRouteLoads_Centre0) FinalRouteLoads_Centre0_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_Centre0_df] FinalRouteLoads_Centre1_df=pd.DataFrame(FinalRouteLoads_Centre1) FinalRouteLoads_Centre1_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_Centre1_df] FinalRouteLoads_Centre2_df=pd.DataFrame(FinalRouteLoads_Centre2) FinalRouteLoads_Centre2_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_Centre2_df] FinalRouteLoads_Centre3_df=pd.DataFrame(FinalRouteLoads_Centre3) FinalRouteLoads_Centre3_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_Centre3_df] #Combine load/route outputs Export_df_Centre0 = pd.concat([FinalRouteDistances_Centre0_df, FinalRouteLoads_Centre0_df], axis=1) Export_df_Centre1 = pd.concat([FinalRouteDistances_Centre1_df, FinalRouteLoads_Centre1_df], axis=1) Export_df_Centre2= pd.concat([FinalRouteDistances_Centre2_df, FinalRouteLoads_Centre2_df], axis=1) Export_df_Centre3 = pd.concat([FinalRouteDistances_Centre3_df, FinalRouteLoads_Centre3_df], axis=1) #Export to CSV Runs_Info_df.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteOptimisation_"+str(NumUsers)+"Info.csv") Export_df_Centre0.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteOptimisation_"+str(NumUsers)+"_Centre0.csv") Export_df_Centre1.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteOptimisation_"+str(NumUsers)+"_Centre1.csv") Export_df_Centre2.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteOptimisation_"+str(NumUsers)+"_Centre2.csv") Export_df_Centre3.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteOptimisation_"+str(NumUsers)+"_Centre3.csv") FinalRouteNodes_Centre0_df=pd.DataFrame(FinalRouteNodes_Centre0) FinalRouteNodes_Centre0_df_Count=FinalRouteNodes_Centre0_df.groupby(FinalRouteNodes_Centre0_df.columns.tolist(),as_index=False).size() FinalRouteNodes_Centre0_df_Count.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteNodes_"+str(NumUsers)+"Capacity_"+str(Cap_Vehicle)+"_Centre0.csv") FinalRouteNodes_Centre1_df=pd.DataFrame(FinalRouteNodes_Centre1) FinalRouteNodes_Centre1_df_Count=FinalRouteNodes_Centre1_df.groupby(FinalRouteNodes_Centre1_df.columns.tolist(),as_index=False).size() FinalRouteNodes_Centre1_df_Count.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteNodes_"+str(NumUsers)+"Capacity_"+str(Cap_Vehicle)+"_Centre1.csv") FinalRouteNodes_Centre2_df=pd.DataFrame(FinalRouteNodes_Centre2) FinalRouteNodes_Centre2_df_Count=FinalRouteNodes_Centre2_df.groupby(FinalRouteNodes_Centre2_df.columns.tolist(),as_index=False).size() FinalRouteNodes_Centre2_df_Count.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteNodes_"+str(NumUsers)+"Capacity_"+str(Cap_Vehicle)+"_Centre2.csv") FinalRouteNodes_Centre3_df=pd.DataFrame(FinalRouteNodes_Centre3) FinalRouteNodes_Centre3_df_Count=FinalRouteNodes_Centre3_df.groupby(FinalRouteNodes_Centre3_df.columns.tolist(),as_index=False).size() FinalRouteNodes_Centre3_df_Count.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyB/"+"RouteNodes_"+str(NumUsers)+"Capacity_"+str(Cap_Vehicle)+"_Centre3.csv") if __name__ == '__main__': pool = Pool() #os.cpu_count() - 1 pool.map(PrimaryFunction, range(0,len(CBSUsersFileSets))) #0,len(CBSUsersFileSets)