### Prepared for paper submission

parent 7bca90b5
This diff is collapsed.
 ... ... @@ -9,19 +9,19 @@ from time import time from lambdasStrategyFunctions import * def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas=None): def Solver(filename, s = None, p = None, S_max=10, nb_Gen=5, UB=None, GivenR_max=None, GivenLambdas=None): # Récupérer le jeu de données (n, T, t0, Rmax, A, R) = DataReader(filename) T = min(T, UB) Rmax = GivenR_max if GivenR_max else Rmax if s == None: s = (1, 1) if p == None: p = (n-2, n-2) A = np.array(A) R = np.array(R) class State: def __init__(self, ids, t, r, parent, valueFunction, **args): def __init__(self, ids, t, r, parent, valueFunction): self.valueFunction = valueFunction self.ids = ids self.parent = parent ... ... @@ -38,7 +38,7 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas return self.t <= other.t and self.cumulr <= other.cumulr def value(self, **args): return self.valueFunction(self.arc, self.t, self.cumulr, **args) return self.valueFunction(self.t, self.cumulr, **args) def dist(self, v): return np.sqrt(pow(self.ids-v, 2)+pow(self.ids-v, 2)) ... ... @@ -75,7 +75,7 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas nextMove.append((x, y-1)) return nextMove A = np.array(A) pq = PriorityQueue() shortest_paths = np.full((n, n), sys.maxsize) pq.put((0, p)) ... ... @@ -88,7 +88,7 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas shortest_paths[v] = shortest_paths[u] + A[e] pq.put((shortest_paths[v], v)) def L(u, path, A, Rmax, T, S_max, nextNodeF, valueF, filtrageF, speedF, lbdas, Deltas): def L(u, path, A, Rmax, T, S_max, nextNodeF, valueF, filtrageF, speedF, lbdas, nb_Gen, Deltas): n = u.arc # Nombre d'arcs déjà parcourus t1 = u.t # Entrer dans l'arc suivant à la date t1 = date d'arrivée en u expansion = dict() ... ... @@ -99,21 +99,26 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas # Moyennes et déviation if vIds not in Deltas: Deltas[vIds] = 0 nb = 0 Delta = 0 Deltas[vIds] = [0, 0, 0] nbDeltas = [0, 0, 0] deltas = [0, 0, 0] # Vérifie si au moins une date à été insérée atLeastOne = False # Pour toutes les dates d'arrivées possibles for t2 in speedF(u, vIds, T, A, R, RDeriv, lbdas): for (fromLbda, t2) in speedF(u, vIds, T, A, R, RDeriv, lbdas, nb_Gen): # Etat d'arrivée potentiel : v v = State(vIds, t2, RtC(T, t1, t2, A[e], np.array(R[e])), u, valueF, lbdas=lbdas) v = State(vIds, t2, RtC(T, t1, t2, A[e], R[e]), u, valueF) # Mise à jour des moyennes et déviation Delta = (nb*Delta + v.cumulr/Rmax-v.dist(s)*Rmax/(v.dist(s)+v.dist(p))) / (nb+1) nb += 1 iD = 1 if fromLbda == 0: iD = 0 elif fromLbda == nb_Gen-1: iD = 2 deltas[iD] = (nbDeltas[iD]*deltas[iD] + v.cumulr/Rmax-v.dist(s)*Rmax/(v.dist(s)+v.dist(p))) / (nbDeltas[iD]+1) nbDeltas[iD] += 1 # Insertion si état v valide if v.cumulr < Rmax: ... ... @@ -137,15 +142,15 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas if atLeastOne == False: # Ajouter la première vitesse valide t2 = t1 idx = bisect_right(np.array(R[e])[:,0], t1)-1 idx = bisect_left(R[e], (t1,)) K = R[e][idx] while t2 < t1+A[e] or (t2 < T and u.cumulr + r > Rmax): t2 += 1 idx += (idx+1 < len(R[e]) and R[e][idx+1] == t2) K += R[e][idx] r = pow(A[e] / (t2 - t1), 2) * K expansion[vIds].append((State(vIds, t2, r, u, valueF, lbdas=lbdas), path+[u])) Deltas[vIds] = (Deltas[vIds]+Delta) / 2 expansion[vIds].append((State(vIds, t2, r, u, valueF), path+[u])) Deltas[vIds] = [(Deltas[vIds][i]+deltas[i]) / 2 for i in range(3)] keys = list(expansion.keys()) for ids in keys: ... ... @@ -153,14 +158,14 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas del expansion[ids] return expansion lbdas = GivenLambdas if GivenLambdas else getLambdas(n, T, Rmax) lbdas = GivenLambdas if GivenLambdas else getLambdas(UB, Rmax) Deltas = dict() valueF = lbdaValue filtrageF = lbdaFiltrage speedF = lbdaSpeed updateF = lbdaUpdate sNode = (0, (State(s, t0, 0, None, valueF, lbdas=lbdas), [])) sNode = (0, (State(s, t0, 0, None, valueF), [])) # dict with key = nodeId, value = [(nodeValue, (node, path w/o node)), ...] Lstar = dict() Lstar[s] = [] ... ... @@ -182,12 +187,12 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas (e, (u, path)) = heappop(Lstar[nextNodeId]) expansion = L(u, path, A, Rmax, T, S_max, nextNodeF, valueF, filtrageF, speedF, lbdas, Deltas) valueF, filtrageF, speedF, lbdas, nb_Gen, Deltas) for vIds in expansion.keys(): if vIds not in Lstar: Lstar[vIds] = [] S_maxs[vIds] = S_max S_maxs[vIds] = 2*S_max for (v, path) in expansion[vIds]: if S_maxs[vIds] != 0: heappush(Lstar[vIds], (v.estimatedValue+(T if v.cumulr>Rmax*v.t/(v.t + shortest_paths[v.ids]) else 0), (v, path))) ... ... @@ -197,9 +202,10 @@ def Solver(filename, s = None, p = None, S_max=10, GivenR_max=None, GivenLambdas nextNodeId = None for vIds in Lstar.keys(): if S_maxs[vIds]: if len(Lstar[vIds]) and\ (not nextNodeId or\ Lstar[vIds].estimatedValue < Lstar[nextNodeId].estimatedValue): if len(Lstar[vIds]) and \ (not nextNodeId or \ Lstar[vIds] < Lstar[nextNodeId] or \ (nextNodeId == p and Lstar[vIds] <= Lstar[nextNodeId])): nextNodeId = vIds return Lstar[p] if len(Lstar[p]) else (State(s, t0, 0, None, valueF, lbdas=lbdas), []) \ No newline at end of file return Lstar[p] if len(Lstar[p]) else (State(s, t0, 0, None, valueF), []) \ No newline at end of file
This diff is collapsed.
This source diff could not be displayed because it is too large. You can view the blob instead.
 ... ... @@ -23,3 +23,4 @@ RUN jupyter nbextension enable code_prettify/autopep8 RUN pip install scanf RUN pip install tqdm RUN pip install joblib
This diff is collapsed.
 ... ... @@ -122,8 +122,8 @@ "filename = \"04.dat\"\n", "(n, T, t0, Rmax, A, R) = DataReader(filename)\n", "(OptPath, OptTs, OptR) = SolReader(filename)\n", "A = np.array(A)\n", "R = np.array(R)\n", "\n", "", "s = (1, 1)\n", "p = (n-2, n-2)\n", "# Nombre de décisions max dans la programmation dynamique\n", ... ...
This diff is collapsed.
 ... ... @@ -8,15 +8,16 @@ from bisect import bisect_left import matplotlib.pyplot as plt from dataReader import DataReader, SolReader from time import time from ppp.lambdasStrategyFunctions import * from ppp.idealStrategyFunctions import * from lambdasStrategyFunctions import * class State: def __init__(self, u, t, r, valueFunction = None): def __init__(self, ids, t, r, parent, valueFunction = None): self.valueFunction = valueFunction self.arc = 1 + (u.arc if u else -1) self.cumulr = r + (u.cumulr if u else 0) self.ids = ids self.parent = parent self.arc = 1 + (parent.arc if parent else -1) self.cumulr = r + (parent.cumulr if parent else 0) self.t = t self.r = r ... ... @@ -26,73 +27,78 @@ class State: def __ge__(self, other): return self.t <= other.t and self.cumulr <= other.cumulr def value(self, **args): return self.valueFunction(self.arc, self.t, self.cumulr, **args) def value(self, lbdas): return self.valueFunction(self.t, self.cumulr, lbdas) def __str__(self): return "(%u, %f)" % (self.t, self.cumulr) def Solver(T, t0, Rmax, D, R, S_max, strategy="lambdas", GivenLambdas=None): def Solver(path, t0, T, Rmax, A, R, S_max, nb_Gen, UB, GivenLambdas=None): # Récupérer le jeu de données n = len(D) n = len(path) totalD = sum([A[e] for e in path]) partialD = 0 V = np.empty((n+1, 2), dtype=int) V = 0 V[n] = T for a in range(1, n+1): V[a] = V[a-1] + D[a-1] V[n-a] = V[n+1-a] - D[n-a] lbdas = GivenLambdas if GivenLambdas else (getLambdas(n, T, Rmax, V) if strategy == "lambdas" else None) valueF = lbdaValue if strategy == "lambdas" else idealValue filtrageF = lbdaFiltrage if strategy == "lambdas" else idealFiltrage speedF = lbdaSpeed if strategy == "lambdas" else idealSpeed updateF = lbdaUpdate if strategy == "lambdas" else idealUpdate lbdas = GivenLambdas if GivenLambdas else getLambdas(UB, Rmax) valueF = lbdaValue filtrageF = lbdaFiltrage speedF = lbdaSpeed updateF = lbdaUpdate from_state = [] # [(value, (node, path w/o node)), ...] heappush(from_state, (0, (State(None, t0, 0, valueF), []))) heappush(from_state, (0, (State(path, t0, 0, None, valueF), []))) to_state = [] # [(node, path), ...] # Jusqu'à la fin de chemin for a in range(n): for a in range(1, n): to_state = [] # Destination vIds = path[a] partialD += A[path[a]] # Pour tous les états à considérer for (e, (u, path)) in from_state: # Etat courant : u for (e, (u, incomplete_path)) in from_state: # Etat courant : u # Arc suivant edge = (u.ids+vIds, u.ids+vIds) # Entrer dans l'arc suivant à la date t1 = date d'arrivée en u t1 = u.t # déviation nb = 0 Delta = 0 nbDeltas = [0, 0, 0] deltas = [0, 0, 0] # Vérifie si au moins une date à été insérée atLeastOne = False # Pour toutes les dates d'arrivées possibles for t2 in speedF(u, D[a], R[a], T, RDerivWrapper(a, D, R), lbdas=lbdas, n=n, Rmax=Rmax, V=V): for (fromLbda, t2) in speedF(u, vIds, T, A, R, RDeriv, lbdas, nb_Gen): # Etat d'arrivée potentiel : v v = State(u, t2, RtC(t1, t2, D[a], R[a]), valueF) v = State(vIds, t2, RtC(T, t1, t2, A[edge], R[edge]), u, valueF) # Mise à jour des moyennes et déviation Delta = (nb*Delta + v.cumulr/Rmax-sum(D[:a+1])/sum(D)) / (nb+1) nb += 1 iD = 1 if fromLbda == 0: iD = 0 elif fromLbda == nb_Gen-1: iD = 2 deltas[iD] = (nbDeltas[iD]*deltas[iD] + v.cumulr/Rmax-partialD/totalD) / (nbDeltas[iD]+1) nbDeltas[iD] += 1 # Insertion si état v valide if v.cumulr < Rmax: atLeastOne = True # Si la liste est vide => on insert sans vérification if len(to_state) == 0: to_state.append((v, path+[u])) to_state.append((v, incomplete_path+[u])) else: # Sinon, on insert dans l'ordre d'arrivée et on vérifie les dominations # Trouver le plus grands état arrivant avant v i = bisect_left(np.array(to_state)[:, 0], v) # Insertion si état v non dominé if not to_state[i-1] >= v: to_state.insert(i, (v, path+[u])) to_state.insert(i, (v, incomplete_path+[u])) # Supprimer les états dominés par v while i < len(to_state)-1 and v >= to_state[i+1]: ... ... @@ -101,27 +107,28 @@ def Solver(T, t0, Rmax, D, R, S_max, strategy="lambdas", GivenLambdas=None): # Si aucune vitesse sélectionnée n'est valide if atLeastOne == False: # Ajouter la première vitesse valide t2 = t1 + D[a] K = sum(R[a][t1:t1+D[a]+1]) r = pow(D[a]/(t2-t1),2)*K while t2 < len(R[a])-1 and u.cumulr + r > Rmax: t2 = t1 idx = bisect_left(R[edge], (t1,)) K = R[edge][idx] while t2 < t1+A[edge] or (t2 < T and u.cumulr + r > Rmax): t2 += 1 K += R[a][t2] r = pow(D[a]/(t2-t1),2)*K to_state.append((State(u, t2, r, valueF), path+[u])) idx += (idx+1 < len(R[edge]) and R[edge][idx+1] == t2) K += R[edge][idx] r = pow(A[edge] / (t2 - t1), 2) * K to_state.append((State(vIds, t2, r, u, valueF), incomplete_path+[u])) # Si le chemin n'est pas fini : évaluer et trier, sinon prendre le chemin le plus rapide if a < n - 1: # Evaluer tous les nouveaux états et les ranger par ordre croissants from_state = [] for (v, path) in to_state: heappush(from_state, (v.value(lbdas=lbdas, n=n, Rmax=Rmax, V=V), (v, path))) for (v, incomplete_path) in to_state: heappush(from_state, (v.value(lbdas=lbdas), (v, incomplete_path))) #Filtrage des nouveaux états from_state = filtrageF(from_state, S_max, lbdas, Delta) from_state = filtrageF(from_state, S_max, lbdas, deltas) # Mettre à jour les hyperparamètres si besoin updateF(lbdas, Delta) updateF(lbdas, deltas) return to_state \ No newline at end of file
This diff is collapsed.
This diff is collapsed.
 import numpy as np from heapq import heappush from scanf import scanf def StringToVector(line): return [int(nbr) for nbr in line.split(',')] return np.array([int(nbr) for nbr in line.split(',')]) def StringToTupleVector(line): splitted = line.split(',') v = [] for i in range(0, len(splitted), 2): v.append((int(splitted[i][1:]), int(splitted[i+1][:-1]))) heappush(v, (int(splitted[i][1:]), int(splitted[i+1][:-1]))) return v def DataReader(fname): ... ... @@ -25,7 +27,7 @@ def DataReader(fname): R.append([None]*(2*n-1)) for j in range(n-(i%2==0)): R[-1][2*j+(i%2==0)] = StringToTupleVector(scanf("%s", f.readline())) return n, T, t0, Rmax, A, R return n, T, t0, Rmax, np.array(A), np.array(R) def SolReader(fname): ... ...
 ... ... @@ -2,49 +2,25 @@ import sys import numpy as np from dijkstar import Graph, find_path from PathPhasingProblem import Solver as pppSolver from lambdasStrategyFunctions import RtC, getLambdas from dataReader import DataReader, SolReader from time import time class PPP: def __init__(self, spath, T, A, R): self.spath = spath self.D = [] self.risk_planning = [] for i in range(len(spath)-1): edge = (spath[i]+spath[i+1], spath[i]+spath[i+1]) self.D.append(A[edge][edge]) self.risk_planning.append([]) t = 0 h = 0 for (l, nexth) in R[edge][edge] : while t < l: self.risk_planning[-1].append(h) t += 1 h = nexth while t <= T: self.risk_planning[-1].append(h) t += 1 def solve(self, t0, T, Rmax, S_max, strategy="lambdas", GivenLambdas=None): return pppSolver(T, t0, Rmax, self.D, self.risk_planning, S_max, strategy, GivenLambdas) def RtC(self, a, t1, t2): return pow(self.D[a] / (t2 - t1), 2) * sum(self.risk_planning[a][t1:t2]) def Solver(filename, s=None, p=None, S_max=10, strategy="lambdas", GivenR_max=None, GivenLambdas=None): def Solver(filename, s=None, p=None, S_max=10, nb_Gen=5, UB=None, GivenR_max=None, GivenLambdas=None): (n, T, t0, Rmax, A, R) = DataReader(filename) T = min(T, UB) Rmax = GivenR_max if GivenR_max else Rmax lbdas = GivenLambdas if GivenLambdas else getLambdas(UB, Rmax) if s == None: s = (1, 1) if p == None: p = (n-2, n-2) if Rmax == 0: return [], [], 0 R = np.array(R) G = Graph(undirected=True) for i in range(n): for j in range(n): ... ... @@ -55,19 +31,21 @@ def Solver(filename, s=None, p=None, S_max=10, strategy="lambdas", GivenR_max=No # Shortest path spath = find_path(G, s, p).nodes # Solve Path Phasing Problem ppp = PPP(spath, T, A, R) (end_node, path_node) = ppp.solve(t0, T, Rmax, S_max, strategy=strategy, GivenLambdas=GivenLambdas) (end_node, path_node) = pppSolver(spath, t0, T, Rmax, A, R, S_max, nb_Gen, UB, GivenLambdas=lbdas) path_node = path_node+[end_node] phase = [node.t for node in path_node] cumulrs = [node.cumulr for node in path_node] risks = [cumulrs[i]-cumulrs[i-1] for i in range(1,len(cumulrs))] def neighbour(ppp, phase): risks = [ppp.RtC(a, phase[a], phase[a+1]) for a in range(len(phase)-1)] def neighbour(path, phase): risks = [] for a in range(len(phase)-1): e = (path[a]+path[a+1], path[a]+path[a+1]) risks.append(RtC(T, phase[a], phase[a+1], A[e], R[e])) edge = np.argsort(risks)[::-1] for a in range(len(edge)): (ui, uj) = ppp.spath[a] (vi, vj) = ppp.spath[a+1] (ui, uj) = path[a] (vi, vj) = path[a+1] (ei, ej) = (ui+vi, uj+vj) G.remove_edge((ui, uj), (vi, vj)) try: ... ... @@ -84,8 +62,7 @@ def Solver(filename, s=None, p=None, S_max=10, strategy="lambdas", GivenR_max=No mPath = spath mPhase = phase visited = {hash(tuple(spath))} ppp = PPP(spath, T, A, R) res = ppp.solve(t0, min(mPhase[-1], T), Rmax, S_max, strategy=strategy, GivenLambdas=GivenLambdas) res = pppSolver(spath, t0, min(mPhase[-1], T), Rmax, A, R, S_max, nb_Gen, UB, GivenLambdas=lbdas) if res != None: (end_node, path_node) = res path_node = path_node+[end_node] ... ... @@ -96,14 +73,14 @@ def Solver(filename, s=None, p=None, S_max=10, strategy="lambdas", GivenR_max=No mPath = spath mPhase = phase path = spath while not STOP: neighbourhood = neighbour(ppp, phase) neighbourhood = neighbour(path, phase) for path in neighbourhood: h = hash(tuple(path)) if h not in visited: visited.add(h) ppp = PPP(path, T, A, R) res = ppp.solve(t0, min(mPhase[-1], T), Rmax, S_max, strategy=strategy, GivenLambdas=GivenLambdas) res = pppSolver(path, t0, mPhase[-1], Rmax, A, R, S_max, nb_Gen, UB, GivenLambdas=lbdas) if res != None: (end_node, path_node) = res path_node = path_node+[end_node] ... ...
 ... ... @@ -16,8 +16,7 @@ def Solver(filename, s = None, p = None, GivenR_max=None): s = (1, 1) if p == None: p = (n-2, n-2) R = np.array(R) G = Graph(undirected=True) for i in range(n): for j in range(n): ... ... @@ -30,8 +29,7 @@ def Solver(filename, s = None, p = None, GivenR_max=None): def Solver(T, t0, Rmax, D, R): arcs = len(D) R = np.array(R) V = np.empty((arcs+1, 2), dtype=int) V = 0 V[arcs] = T ... ...
 from bisect import bisect_right from bisect import bisect_left import numpy as np def RtC(T, t1, t2, d, r): idl = bisect_right(r[:,0], t1)-1 idr = bisect_right(r[:,0], t2, lo=idl) K = sum([0 if t1 > T else (T-max(r[i], t1))*r[i] if t2 > T else (min(r[i+1], t2)-max(r[i], t1))*r[i] for i in range(idl, min(len(r)-1,idr))]) idl = bisect_left(r, (t1,)) idl -= t1 != r[idl] idr = bisect_left(r, (t2,), lo=idl) K = 0 for i in range(idl, min(len(r)-1,idr)): if t1 < T: K += (T-max(r[i], t1))*r[i] if t2 > T else (min(r[i+1], t2)-max(r[i], t1))*r[i] return RtCK(t1, t2, d, K) RtCK = lambda t1, t2, d, K : pow(d / (t2 - t1), 2) * K def getLambdas(n, T, Rmax, V = None): gamma = [T / 4 / Rmax, T / 2 / Rmax, T / Rmax] return [.2, .5, .8] return [gamma / (1 + gamma), gamma / (1 + gamma), gamma / (1 + gamma)] def getLambdas(UB, Rmax): lbda = Rmax / UB return [lbda / 4, lbda, lbda * 4] # Dérivé de la fonction de risque def RDeriv(t1, t2, d, r): # Rechercher le début du plateau de risque lu = bisect_right(r[:,0], t1)-1 lk = bisect_right(r[:,0], t2, lo=lu)-1 K = sum([(r[i+1]-max(r[i], t1))*r[i] for i in range(lu, lk-1)]) return pow(d/(t2-t1), 2)*(-r[lk] if lk == lu else (-2/(t2-t1)*K+r[lk]*(-2*(t2-r[lk])/(t2-t1)+1))) idl = bisect_left(r, (t1,)) idl -= t1 != r[idl] idr = bisect_left(r, (t2,), lo=idl) K = 0 for i in range(idl, min(len(r)-1,idr+1)): K += (min(r[i+1], t2)-max(r[i], t1))*r[i] return pow(d/(t2-t1), 2)*(-2/(t2-t1)*K+r[idr]) # Recherche dichotomique def dicho(t1, a, b, f, e, *args): ... ... @@ -37,34 +43,39 @@ def dicho(t1, a, b, f, e, *args): ecart = fin - debut return m def lbdaValue(a, t, r, lbdas, **args): return lbdas * r + (1 - lbdas) * t def lbdaValue(t, r, lbdas, **args): return r + lbdas * t def lbdaUpdate(lbdas, Delta): return def lbdaUpdate(lbdas, Deltas): #print("Delta=",Delta) for i in range(len(lbdas)): lbdas[i] = max(0, min(1, lbdas[i]+(lbdas[i] if Delta < 0 else 1-lbdas[i]) * Delta)) lbdas[i] = lbdas[i]*(1-.2*Deltas[i]) def lbdaDeriv(t1, t2, d, r, RDeriv, lbda): return RDeriv(t1, t2, d, r) * lbda + (1 - lbda) return RDeriv(t1, t2, d, r) + lbda def lbdaWrapper(RDeriv, lbda): return lambda t1, t2, d, r: lbdaDeriv(t1, t2, d, r, RDeriv, lbda)