from vrptw_base import VrptwGraph from multiple_ant_colony_system import MultipleAntColonySystem import numpy as np import csv import math import random import networkx as nx import matplotlib.pyplot as plt class Node: """ Attributes ---------- index : int index ubigeo : str 6 digits lat : float Latitud (angulo de los Paralelos) en grados sexagesimales lon : float Longitud (angulo de los Meridianos) en grados sexagesimales is_depot : bool demand : int ready_time : float due_time : float service_time : float x : float Aproximacion local, en Km respecto y : float Aproximacion local, en Km respecto ("Lima") Notes ----- Web Mercator projection (BAD): x = floor(256 / (2 * math.pi) * 2**(zoom_level) * (lon + math.pi)) y = floor(265 / (2 * math.pi) * 2**(zoom_level) * (math.pi - math.ln(math.tan( math.pi / 4 + lat / 2 )))) x = R * lon y = R * ln(tan(pi/4 + lat/2) Both `lon` and `lat` in radians. "Lima": -12.04591952,-77.03049615 (lat, long) """ def __init__(self, index: int, ubigeo, lat, lon, is_depot, demand, ready_time, due_time, service_time): super() self.index = index self.ubigeo = ubigeo if is_depot: self.is_depot = True else: self.is_depot = False earth_radius_km = 6371 # Avg. radius lima_lat = -12.04591952 lima_lon = -77.03049615 self.lat = lat self.lon = lon self.x = (lon - lima_lon) * (math.pi / 180) * earth_radius_km self.y = (lat - lima_lat) * (math.pi / 180) * earth_radius_km self.x = round(self.x, 3) self.y = round(self.y, 3) self.demand = demand self.ready_time = 0 # ready_time self.due_time = due_time self.service_time = service_time def _deg2rad(degrees): return degrees * math.pi / 180 def distance_between_coordinates(lat1, lon1, lat2, lon2): """Returns the distance in km""" earth_radius_km = 6371 # Avg. radius lat1 = _deg2rad(lat1) lat2 = _deg2rad(lat2) d_lat = lat2 - lat1 d_lon = lon2 - lon1 a = math.sin(d_lat/2) * math.sin(d_lat/2) \ + math.sin(d_lon/2) * math.sin(d_lon/2) * math.cos(lat1) * math.cos(lat2) c = 2 * math.atan2(math.sqrt(a), math.sqrt(1-a)) return earth_radius_km * c def _read_nodes(depot_num): """ Read nodes from input file. Attributes ---------- depot_num : int Number of depots Returns ------- list List of Nodes. The first `depot_num` nodes are the depots. """ depots = [] no_depots = [] with open('odiparpack/inf226.oficinas_mod.csv', newline='') as csvfile: orders = csv.reader(csvfile) count = 1 index_customer = depot_num index_depot = 0 for row in orders: if count == 1: count += 1 continue ubigeo, dept, prov, lat, lon, region_natural, is_depot = row[:7] demand, ready_time, due_time, service_time = row[7:] n = Node( -1, ubigeo, float(lat), float(lon), int(is_depot), demand=float(demand), ready_time=0, due_time=float(due_time), service_time=60 ) if n.is_depot: n.index = index_depot depots.append(n) index_depot += 1 else: n.index = index_customer no_depots.append(n) index_customer += 1 count += 1 return depots + no_depots def _read_tramos(): """ Lee archivo de tramos Returns ------- list Lista de tuplas con los tramos (ubigeo1, ubigeo2) """ tramos = [] with open('odiparpack/inf226.tramos.v.2.0.csv', newline='') as f: rows = csv.reader(f) count = 1 for row in rows: if count >= 2: tramos.append([row[0], row[1]]) count += 1 return tramos def make_complete_customer_node_graph(nodes, tramos): """ Nodes + Tramos make a network. An undirected graph. But we only care about the nodes that are either have demand (aka. customers), or are depots. This function creates a sub-graph of only depot-or-customer nodes, and finds the nearest path + path_lenght between them (on the original network). Parameters ---------- nodes : list List of Node objects tramos : list List of tuples of "ubigeo" Returns ------- list List of depot-or-customer nodes, Shortest distance matrix (np.narray), Dict with shortest distance and path object for all pairs of nodes """ elist = [] # Load sparse graph nodes_d = {n.ubigeo: n for n in nodes} for ubigeo1, ubigeo2 in tramos: n1 = nodes_d[ubigeo1] n2 = nodes_d[ubigeo2] dist = round(distance_between_coordinates(n1.lat, n1.lon, n2.lat, n2.lon), 0) elist.append((ubigeo1, ubigeo2, dist)) G = nx.Graph() G.add_weighted_edges_from(elist) # Calculate missing edges using the shortest path # (Thus making the graph complete) len_path = dict(nx.all_pairs_dijkstra(G, weight='weight')) # for node, (dist, path) in len_path: # print(f"{node}: {dist} [{path}]") # Prune all nodes that have no demand and "reset index" of nodes depot_customer_nodes = [] i = 0 for n in nodes: if n.demand > 0 or n.is_depot: n.index = i depot_customer_nodes.append(n) i += 1 node_num = i # Create distance matrix node_dist_mat = np.zeros((node_num, node_num)) for n1 in depot_customer_nodes: for n2 in depot_customer_nodes: node_dist_mat[n1.index, n2.index] = \ len_path[n1.ubigeo][0][n2.ubigeo] return depot_customer_nodes, node_dist_mat, len_path def lab(): depot_num = 3 nodes = _read_nodes(depot_num) tramos = _read_tramos() nodes, node_dist_mat, len_path = \ make_complete_customer_node_graph(nodes, tramos) # nodes with demand for n in nodes[:8]: print(n.__dict__) # shortest distance between nodes in network print(node_dist_mat[:8, :8]) # shortest path for n1 in nodes[:1]: for n2 in nodes[:8]: shortest_path = len_path[n1.ubigeo][1][n2.ubigeo] print(f"{n1.ubigeo} => {n2.ubigeo}: [{shortest_path}]") def main(): #file_path = './solomon-100/c101.txt' file_path = './odiparpack/pedidosperu195.txt' ants_num = 10 beta = 2 # 5 q0 = 0.1 # 0.5 ? show_figure = True graph = VrptwGraph(file_path) macs = MultipleAntColonySystem(graph, ants_num=ants_num, beta=beta, q0=q0, whether_or_not_to_show_figure=show_figure) macs.run_multiple_ant_colony_system() if __name__ == '__main__': lab()