import numpy as np
import shed
import matplotlib.pyplot as plt
from scipy.spatial import Delaunay
import networkx as nx
from scipy.spatial import distance
from scipy.optimize import linprog
class Path:
def __init__ (self,C,f):
self.C = C
self.flow = f
def getK (self,p):
d = 0
for i in range(len(self.C)-1):
d += distance.euclidean(p[self.C[i]],p[self.C[i+1]])
return d/distance.euclidean(p[self.C[0]],p[self.C[-1]])
def __repr__(self):
return " %s<-%s " % (self.C, self.flow)
def __str__(self):
return " %s<-%s " % (self.C, self.flow)
def __eq__(self, other):
if len(self.C) != len(other.C):
return False
for i in range(len(self.C)):
if self.C[i]!=other.C[i]:
return False
return True
def getMinFlow(self):
if len(self.C)==1:
return 0
min_flow = np.inf
for i in range(len(self.C)-1):
if G[self.C[i]][self.C[i+1]]['capacity']<min_flow:
min_flow = G[self.C[i]][self.C[i+1]]['capacity']
return min_flow
def complete(self,v,G):
for i in range(len(self.C)-1):
if G[self.C[i]][self.C[i+1]]['problem']:
return False
if self.getMinFlow()>v:
return True
else:
return False
def builtRandomTransit (n):
f = []
for i in range(n):
f.append([0] * n)
for i in range(n):
for j in range(n):
f[i][j]=np.random.randint(10)+1
if i==j:
f[i][j]=0
return f
def check(lst, sub):
for i in range(0, len(lst)):
if lst[i:i+len(sub)] == sub:
return True
return False
def simplexNewPaths (G, loupe_paths, loupe_edges, V,p):
print ("\n Work simplex mathod")
k = [i.getK(p) for i in loupe_paths]
loupe_paths = [i.C for i in loupe_paths]
size = len(loupe_paths)
A_e = []
A_u = []
b_e = []
b_u = [d['capacity'] for u,v,d in loupe_edges]
d = {}
for i in loupe_paths:
if d.get((i[0],i[-1])) is None:
d[(i[0],i[-1])] = []
if i not in d[(i[0],i[-1])]:
d[(i[0],i[-1])].append(i)
for key, value in d.items():
if len(value)>1:
print("{0}: {1}".format(key,value))
A_e.append([0]*size)
for j in value:
xxx = list(d.keys()).index(key)
yyy = loupe_paths.index(j)
A_e[xxx][yyy] = 1
b_e.append(V[key[0]][key[1]])
print (loupe_paths)
shed.ppprint(A_e, len(A_e), size)
print(b_e)
for i in loupe_edges:
A_u.append([0]*size)
for j in loupe_paths:
if check(j,[i[0], i[1]]):
A_u[-1][loupe_paths.index(j)]=1
shed.ppprint(A_u, len(A_u), size)
print(b_u)
print (k)
res = linprog(k, A_ub=A_u, b_ub=b_u, A_eq = A_e, b_eq=b_e,bounds=(0, None))
print(res.x)
return res.x
# return old_paths, new_paths
#Создание и отрисовка случайного графа при помощи триангуляции Делоне
N = 6
p = shed.builtPoints(N,10,10)
print(p)
#p = [[6, 6], [8, 4], [2, 1], [1, 6], [1, 9], [4, 2]]
tri = Delaunay(p)
e,f = shed.getEdgesDelaunay(tri)
G = nx.DiGraph()
for edge in e:
G.add_edge(edge[0],edge[1],weight = round(distance.euclidean(p[edge[0]],p[edge[1]]),3), capacity = np.random.randint(50)+1,flow=[], problem = True)
G.add_edge(edge[1],edge[0],weight = round(distance.euclidean(p[edge[0]],p[edge[1]]),3), capacity = np.random.randint(50)+1,flow=[], problem = True)
#G.add_edge(2,5,weight = round(distance.euclidean(p[2],p[5]),3), capacity = 9,flow=[], problem = True)
#G.add_edge(2,3,weight = round(distance.euclidean(p[2],p[3]),3), capacity = 6,flow=[], problem = True)
#G.add_edge(5,2,weight = round(distance.euclidean(p[5],p[2]),3), capacity = 46,flow=[], problem = True)
#G.add_edge(5,3,weight = round(distance.euclidean(p[5],p[3]),3), capacity = 2,flow=[], problem = True)
#G.add_edge(5,0,weight = round(distance.euclidean(p[5],p[0]),3), capacity = 50,flow=[], problem = True)
#G.add_edge(5,1,weight = round(distance.euclidean(p[5],p[1]),3), capacity = 42,flow=[], problem = True)
#G.add_edge(3,2,weight = round(distance.euclidean(p[3],p[2]),3), capacity = 19,flow=[], problem = True)
#G.add_edge(3,5,weight = round(distance.euclidean(p[3],p[5]),3), capacity = 24,flow=[], problem = True)
#G.add_edge(3,0,weight = round(distance.euclidean(p[3],p[0]),3), capacity = 3,flow=[], problem = True)
#G.add_edge(3,4,weight = round(distance.euclidean(p[3],p[4]),3), capacity = 34,flow=[], problem = True)
#G.add_edge(0,5,weight = round(distance.euclidean(p[0],p[5]),3), capacity = 16,flow=[], problem = True)
#G.add_edge(0,3,weight = round(distance.euclidean(p[0],p[3]),3), capacity = 9,flow=[], problem = True)
#G.add_edge(0,1,weight = round(distance.euclidean(p[0],p[1]),3), capacity = 28,flow=[], problem = True)
#G.add_edge(0,4,weight = round(distance.euclidean(p[0],p[4]),3), capacity = 7,flow=[], problem = True)
#G.add_edge(1,0,weight = round(distance.euclidean(p[1],p[0]),3), capacity = 3,flow=[], problem = True)
#G.add_edge(1,5,weight = round(distance.euclidean(p[1],p[5]),3), capacity = 23,flow=[], problem = True)
#G.add_edge(4,0,weight = round(distance.euclidean(p[4],p[0]),3), capacity = 21,flow=[], problem = True)
#G.add_edge(4,3,weight = round(distance.euclidean(p[4],p[3]),3), capacity = 44,flow=[], problem = True)
positions_vertexes = [(p[i][0], p[i][1]) for i in range(N)]
#Создание случайной матрциы перевозок
V = builtRandomTransit(N)
print (V)
#V = [[0, 1, 3, 3, 5, 3], [5, 0, 3, 3, 4, 5], [2, 4, 0, 4, 1, 3], [3, 3, 5, 0, 3, 2], [3, 3, 5, 4, 0, 1], [3, 3, 4, 1, 5, 0]]
# shed.pprint (V)
result_flows = [[[] for i in range(N)] for j in range(N)]
#Поиск всех кратчайших путей
path = nx.all_pairs_shortest_path(G)
result_paths = []
for i in range(N):
for j in range(N):
if path[i][j] is not None:
for k in range(len(path[i][j])-1):
result_paths.append(Path(path[i][j],V[i][j]))
G[path[i][j][k]][path[i][j][k+1]]['flow'].append(result_paths[-1])
else:
if V[i][j]>0:
print ("Задача неразрешима")
exit(1)
# edge_labels=dict([((u,v,),(d['capacity'])) for u,v,d in G.edges(data=True)])
# nx.draw_networkx_edge_labels(G,positions_vertexes,edge_labels=edge_labels,label_pos=0.75)
# nx.draw_networkx(G, positions_vertexes, with_labels=True, arrows=True, node_color='Red')
# plt.show()
print ("\nПути\n",path)
problem_edges = []
for u,v,d in G.edges(data=True):
sum_flow = sum([paths_edge.flow for paths_edge in d['flow'] ])
if sum_flow<=d['capacity']:
d['problem']=False
else:
problem_edges.append([u,v,d])
if len(problem_edges) == 0:
print ("Решение получено")
exit(1)
print ("\nВсе ребра")
for u,v,d in G.edges(data=True):
print (d['problem'], u,v,d['flow'],d['capacity'])
print ("\nРезультурующие пути\n",result_paths)
simplex_result_paths = []
print("\nПроблемные ребра\n", problem_edges)
capacities = []
for u,v,d in problem_edges:
G.remove_edge(u,v)
for pa in d['flow']:
if pa not in simplex_result_paths:
simplex_result_paths.append(pa)
tmp_pp = Path(nx.shortest_path(G,pa.C[0],pa.C[-1]),None)
if tmp_pp not in simplex_result_paths:
simplex_result_paths.append(tmp_pp)
G.add_edge(u,v,d)
print()
resx =simplexNewPaths(G,simplex_result_paths,problem_edges, V, p)
for i in range(len(resx)):
simplex_result_paths[i].flow = resx[i]
for i in range(len(resx)):
if i%2 == 0:
if simplex_result_paths[i].flow > 0.0:
for j in range(len(simplex_result_paths[i].C) - 1):
if simplex_result_paths[i] in G[simplex_result_paths[i].C[j]][simplex_result_paths[i].C[j + 1]]['flow']:
index = G[simplex_result_paths[i].C[j]][simplex_result_paths[i].C[j + 1]]['flow'].index(simplex_result_paths[i])
G[simplex_result_paths[i].C[j]][simplex_result_paths[i].C[j + 1]]['flow'][index].flow = simplex_result_paths[i].flow
else:
for j in range(len(simplex_result_paths[i].C) - 1):
if simplex_result_paths[i] in G[simplex_result_paths[i].C[j]][simplex_result_paths[i].C[j + 1]]['flow']:
G[simplex_result_paths[i].C[j]][simplex_result_paths[i].C[j + 1]]['flow'].remove(simplex_result_paths[i])
else:
if simplex_result_paths[i].flow != 0.0:
for j in range(len(simplex_result_paths[i].C)-1):
if simplex_result_paths[i] not in G[simplex_result_paths[i].C[j]][simplex_result_paths[i].C[j + 1]]['flow']:
G[simplex_result_paths[i].C[j]][simplex_result_paths[i].C[j+1]]['flow'].append(simplex_result_paths[i])
for u,v,d in G.edges(data=True):
sum_flow = sum([paths_edge.flow for paths_edge in d['flow'] ])
if sum_flow<=d['capacity']:
d['problem']=False
print ("\nВсе ребра")
for u,v,d in G.edges(data=True):
print (d['problem'], u,v,d['flow'],d['capacity'])