# -*- coding: utf-8 -*- """ Created on Fri Jan 15 12:27:14 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 from multiprocessing import Pool from pathlib import Path import numpy as np from scipy.spatial import distance from itertools import repeat #%% 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" #LocationIdentifier Cap_Vehicles=50 TransferStations="TransferStationsLocations_"+Location+".csv" #Location of Transfer Stations NumSetsArray=np.array([[10, 20, 30,50,70,100,250,500,1000,3000], #This defines the number of CBS users and the number of iterations for each [50, 50, 50,50,50,50 ,25, 25, 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"] #%% 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 and process transfer station data (defined in pre-processing) TransferStationLocations = pd.read_csv (FilePath+"/InputFiles/"+Location+"_Inputs/OtherData/"+TransferStations) TransferStationLocations=TransferStationLocations[["POINT_X","POINT_Y"]] TS_Coords = [tuple(x) for x in TransferStationLocations.to_numpy()] TransferStationLocations = pd.read_csv (FilePath+"/InputFiles/"+Location+"_Inputs/OtherData/"+TransferStations) TransferStationPatterns=TransferStationLocations[["{Delta1}","{Delta2}","{Delta3}","{Delta4}","{Delta5}"]] TransferStationPatterns.columns = ['{Delta1Name}', '{Delta1Name}','{Delta1Name}','{Delta1Name}','{Delta1Name}'] #Names of different transfer stations density #%% Computation of Strategy C (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_","")) #Initialise results lists FinalRouteDistances_Distance100=[] FinalRouteDistances_Distance200=[] FinalRouteDistances_Distance300=[] FinalRouteDistances_Distance500=[] FinalRouteDistances_Distance1000=[] FinalRouteLoads_Distance100=[] FinalRouteLoads_Distance200=[] FinalRouteLoads_Distance300=[] FinalRouteLoads_Distance500=[] FinalRouteLoads_Distance1000=[] Runs_InfoList=[] for j in range(0,NumSets): #Import distance matrices from previous preprocessing 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") CBSUserLocations=pd.read_csv(FilePath+"/InputFiles/"+Location+"_Inputs/"+CBSUsersFileSets[i]+"/"+CBSUsersFileSets[i]+"_Set_"+str(j)+ "_CBSPoints.csv") CBSUserLocations=CBSUserLocations[["POINT_X","POINT_Y"]].values CBSUserLocations=[tuple(CBSUserLocations[l]) for l in range(0,len(CBSUserLocations))] DM_TS_Start=DistanceMatrix_Info.iloc[0,0]+DistanceMatrix_Info.iloc[1,0]+DistanceMatrix_Info.iloc[2,0]-1 DM_TS_End=DistanceMatrix_Info.iloc[0,0]+DistanceMatrix_Info.iloc[1,0]+DistanceMatrix_Info.iloc[2,0]+DistanceMatrix_Info.iloc[3,0]-1 TransferStationPatterns1=TransferStationPatterns.assign(Node=range(DM_TS_Start,DM_TS_End)) for t in range(0,5): #For each of the different transfer station densities my_file = Path(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyC/"+"RouteOptimisation_"+str(NumUsers)+"_t"+str(t)+".csv") if my_file.is_file(): # file exists pass else: #Assign transfer station locations ListTransferStations=[] [ListTransferStations.append(TransferStationPatterns1.iloc[x,5]) for x in range(0,len(TransferStationPatterns1)) if TransferStationPatterns1.iloc[x,t]==1] ListTransferCoords=[] [ListTransferCoords.append(TS_Coords[x]) for x in range(0,len(TransferStationPatterns1)) if TransferStationPatterns1.iloc[x,t]==1] ArrayTransferCoords=np.array(ListTransferCoords) #Assigns transfer station loads ListTransferWeights= [0 for i in range(0,len(ListTransferStations))] ListTransferDistances= [] for x in range(0, len(CBSUserLocations)): ListTransferWeights[closest_node(CBSUserLocations[x],ArrayTransferCoords)]+=1 ListTransferDistances.append(closest_node_distance(CBSUserLocations[x],ArrayTransferCoords)) Final_ListTransferWeights=[] [Final_ListTransferWeights.append(ListTransferWeights[x]) for x in range(0,len(ListTransferWeights)) if ListTransferWeights[x]!=0] Final_ListTransferStations=[] [Final_ListTransferStations.append(ListTransferStations[x]) for x in range(0,len(ListTransferWeights)) if ListTransferWeights[x]!=0] #Adds the entrance/exit line Final_ListTransferWeights=[0]+Final_ListTransferWeights Final_ListTransferStations=[0]+Final_ListTransferStations DistanceMatrix_Tailored=DistanceMatrix.iloc[Final_ListTransferStations,Final_ListTransferStations] DistanceMatrix_Altered=DistanceMatrix_Tailored #This ensures that if the node load exceeds the vehicle capacity that multiple vehicles are dispatched to the same point if max(Final_ListTransferWeights)>Cap_Vehicles: for p in range(0, len(Final_ListTransferStations)): if Final_ListTransferWeights[p] > Cap_Vehicles: NVisits=int(np.ceil(Final_ListTransferWeights[p]/Cap_Vehicles)) StationRow=DistanceMatrix_Tailored.iloc[p,:] StationRows=pd.DataFrame([StationRow,]*(NVisits-1)) DistanceMatrix_Altered=DistanceMatrix_Altered.append(StationRows) Final_ListTransferWeights[p] = Final_ListTransferWeights[p] % Cap_Vehicles Final_ListTransferWeights.extend(repeat(Cap_Vehicles, (NVisits-1))) ## ## 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) Num_Vehicle=int(np.ceil((NumUsers*1.1)/Cap_Vehicles)) Cap_Vehicles_Input=[Cap_Vehicles] * 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 # yapf: disable data['num_vehicles'] = Num_Vehicle data['depot'] = 0 #This has been scripted from the 'Exit Point' defined in Script 1 data['demands'] = Final_ListTransferWeights data['vehicle_capacities'] = Cap_Vehicles_Input return data #OutputRouteDistances=[] #OutputRouteLoads=[] def print_solution(data, manager, routing, solution): """Prints solution on console.""" global OutputRouteDistances global OutputRouteLoads 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 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) 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 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): 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 #search_parameters.solution_limit = 100 # 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 #Postprocessing depending on transfer station density Runs_InfoList.append([j]+[t]+[Cap_Vehicles]+[Num_Vehicle]+Final_ListTransferStations+Final_ListTransferWeights) if t==0: FinalRouteDistances_Distance100.append(OutputRouteDistances) FinalRouteLoads_Distance100.append(OutputRouteLoads) elif t==1: FinalRouteDistances_Distance200.append(OutputRouteDistances) FinalRouteLoads_Distance200.append(OutputRouteLoads) elif t==2: FinalRouteDistances_Distance300.append(OutputRouteDistances) FinalRouteLoads_Distance300.append(OutputRouteLoads) elif t==3: FinalRouteDistances_Distance500.append(OutputRouteDistances) FinalRouteLoads_Distance500.append(OutputRouteLoads) elif t==4: FinalRouteDistances_Distance1000.append(OutputRouteDistances) FinalRouteLoads_Distance1000.append(OutputRouteLoads) ## #Postprocessing for each transfer station density FinalRouteDistances_100_df=pd.DataFrame(FinalRouteDistances_Distance100) FinalRouteDistances_100_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_100_df] FinalRouteLoads_100_df=pd.DataFrame(FinalRouteLoads_Distance100) FinalRouteLoads_100_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_100_df] Export_df_100 = pd.concat([FinalRouteDistances_100_df, FinalRouteLoads_100_df], axis=1) Export_df_100.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyC/"+"RouteOptimisation_"+str(NumUsers)+"_t0.csv") FinalRouteDistances_200_df=pd.DataFrame(FinalRouteDistances_Distance200) FinalRouteDistances_200_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_200_df] FinalRouteLoads_200_df=pd.DataFrame(FinalRouteLoads_Distance200) FinalRouteLoads_200_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_200_df] Export_df_200 = pd.concat([FinalRouteDistances_200_df, FinalRouteLoads_200_df], axis=1) Export_df_200.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyC/"+"RouteOptimisation_"+str(NumUsers)+"_t1.csv") FinalRouteDistances_300_df=pd.DataFrame(FinalRouteDistances_Distance300) FinalRouteDistances_300_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_300_df] FinalRouteLoads_300_df=pd.DataFrame(FinalRouteLoads_Distance300) FinalRouteLoads_300_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_300_df] Export_df_300 = pd.concat([FinalRouteDistances_300_df, FinalRouteLoads_300_df], axis=1) Export_df_300.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyC/"+"RouteOptimisation_"+str(NumUsers)+"_t2.csv") FinalRouteDistances_500_df=pd.DataFrame(FinalRouteDistances_Distance500) FinalRouteDistances_500_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_500_df] FinalRouteLoads_500_df=pd.DataFrame(FinalRouteLoads_Distance500) FinalRouteLoads_500_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_500_df] Export_df_500 = pd.concat([FinalRouteDistances_500_df, FinalRouteLoads_500_df], axis=1) Export_df_500.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyC/"+"RouteOptimisation_"+str(NumUsers)+"_t3.csv") FinalRouteDistances_1000_df=pd.DataFrame(FinalRouteDistances_Distance1000) FinalRouteDistances_1000_df.columns = ['FinalDistanceRoute_'+str(c) for c in FinalRouteDistances_1000_df] FinalRouteLoads_1000_df=pd.DataFrame(FinalRouteLoads_Distance1000) FinalRouteLoads_1000_df.columns = ['FinalLoadRoute_'+str(c) for c in FinalRouteLoads_1000_df] Export_df_1000 = pd.concat([FinalRouteDistances_1000_df, FinalRouteLoads_1000_df], axis=1) Export_df_1000.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyC/"+"RouteOptimisation_"+str(NumUsers)+"_t4.csv") Runs_Info_df=pd.DataFrame(Runs_InfoList) Runs_Info_df.to_csv(FilePath+"/OutputFiles/"+Location+"_Outputs/StrategyC/"+"RouteOptimisation_"+str(NumUsers)+"Info.csv") if __name__ == '__main__': pool = Pool() #os.cpu_count() - 1 pool.map(PrimaryFunction, range(0,len(CBSUsersFileSets)))