added the optimizer
This commit is contained in:
		
							
								
								
									
										206
									
								
								backend/src/optimizer.py
									
									
									
									
									
										Normal file
									
								
							
							
						
						
									
										206
									
								
								backend/src/optimizer.py
									
									
									
									
									
										Normal file
									
								
							| @@ -0,0 +1,206 @@ | ||||
| from scipy.optimize import linprog | ||||
| import numpy as np | ||||
| from scipy.linalg import block_diag | ||||
|  | ||||
|  | ||||
| # Defines the landmark class (aka some place there is to visit) | ||||
| class landmark : | ||||
|     def __init__(self, name: str, attractiveness: int, loc: tuple): | ||||
|         self.name = name | ||||
|         self.attractiveness = attractiveness | ||||
|         self.loc = loc | ||||
|  | ||||
| # Convert the result (edges from j to k like d_25 = edge between vertex 2 and vertex 5) into the list of indices corresponding to the landmarks | ||||
| def untangle(res: list) : | ||||
|     N = len(res)                # length of res | ||||
|     L = int(np.sqrt(N))         # number of landmarks. CAST INTO INT but should not be a problem because N = L**2 by def. | ||||
|     n_landmarks = res.sum()     # number of visited landmarks | ||||
|     visit_order = [] | ||||
|     cnt = 0 | ||||
|  | ||||
|     if n_landmarks % 2 == 1 :                                     # if odd number of visited checkpoints | ||||
|         for i in range(L) : | ||||
|             for j in range(L) : | ||||
|                 if res[i*L + j] == 1 :              # if index is 1 | ||||
|                     cnt += 1                        # increment counter | ||||
|                     if cnt % 2 == 1 :               # if counter odd | ||||
|                         visit_order.append(i) | ||||
|                         visit_order.append(j) | ||||
|     else :                                   # if even number of ones | ||||
|         for i in range(L) : | ||||
|             for j in range(L) : | ||||
|                 if res[i*L + j] == 1 :              # if index is one | ||||
|                     cnt += 1                        # increment counter | ||||
|                     if j % (L-1) == 0 :             # if last node | ||||
|                         visit_order.append(j)       # append only the last index | ||||
|                         return visit_order          # return | ||||
|                     if cnt % 2 == 1 :  | ||||
|                         visit_order.append(i) | ||||
|                         visit_order.append(j) | ||||
|     return visit_order | ||||
|  | ||||
| # Just to print the result | ||||
| def print_res(res: list, P) : | ||||
|     X = abs(res.x) | ||||
|     order = untangle(X) | ||||
|  | ||||
|     # print("Optimal value:", -res.fun)  # Minimization, so we negate to get the maximum | ||||
|     # print("Optimal point:", res.x) | ||||
|     # N = int(np.sqrt(len(X))) | ||||
|     # for i in range(N): | ||||
|     #     print(X[i*N:i*N+N]) | ||||
|     # print(order) | ||||
|  | ||||
|     if (X.sum()+1)**2 == len(X) :  | ||||
|         print('\nAll landmarks can be visited within max_steps, the following order is most likely not the fastest') | ||||
|     else : | ||||
|         print('Could not visit all the landmarks, the following order could be the fastest but not sure') | ||||
|     print("Order of visit :") | ||||
|     for i, elem in enumerate(landmarks) :  | ||||
|         if i in order : print('- ' + elem.name) | ||||
|  | ||||
|     steps = path_length(P, abs(res.x)) | ||||
|     print("\nSteps walked : " + str(steps)) | ||||
|  | ||||
| # Constraint to use only the upper triangular indices for travel | ||||
| def break_sym(landmarks, A_eq, b_eq): | ||||
|     L = len(landmarks) | ||||
|     l = [0]*L*L | ||||
|     for i in range(L) : | ||||
|         for j in range(L) : | ||||
|             if i >= j : | ||||
|                 l[j+i*L] = 1 | ||||
|  | ||||
|     A_eq = np.vstack((A_eq,l)) | ||||
|     b_eq.append(0) | ||||
|  | ||||
|     return A_eq, b_eq | ||||
|  | ||||
| # Constraint to respect max number of travels | ||||
| def respect_number(landmarks, A_ub, b_ub): | ||||
|     h = [] | ||||
|     for i in range(len(landmarks)) : h.append([1]*len(landmarks)) | ||||
|     T = block_diag(*h) | ||||
|     return np.vstack((A_ub, T)), b_ub + [1]*len(landmarks) | ||||
|  | ||||
| # Constraint to tie the problem together and have a connected path | ||||
| def respect_order(landmarks: list, A_eq, b_eq):  | ||||
|     N = len(landmarks) | ||||
|     for i in range(N-1) :     # Prevent stacked ones | ||||
|         if i == 0 : | ||||
|             continue | ||||
|         else :  | ||||
|             l = [0]*N | ||||
|             l[i] = -1 | ||||
|             l = l*N | ||||
|             for j in range(N) : | ||||
|                 l[i*N + j] = 1 | ||||
|  | ||||
|             A_eq = np.vstack((A_eq,l)) | ||||
|             b_eq.append(0) | ||||
|  | ||||
|     return A_eq, b_eq | ||||
|  | ||||
| # Compute manhattan distance between 2 locations | ||||
| def manhattan_distance(loc1: tuple, loc2: tuple): | ||||
|     x1, y1 = loc1 | ||||
|     x2, y2 = loc2 | ||||
|     return abs(x1 - x2) + abs(y1 - y2) | ||||
|  | ||||
| # Constraint to not stay in position | ||||
| def init_eq_not_stay(landmarks):  | ||||
|     L = len(landmarks) | ||||
|     l = [0]*L*L | ||||
|     for i in range(L) : | ||||
|         for j in range(L) : | ||||
|             if j == i : | ||||
|                 l[j + i*L] = 1 | ||||
|     #A_eq = np.array([np.array(xi) for xi in A_eq])                  # Must convert A_eq into an np array | ||||
|     l = np.array(np.array(l)) | ||||
|     return [l], [0] | ||||
|  | ||||
| # Initialize A and c. Compute the distances from all landmarks to each other and store attractiveness | ||||
| # We want to maximize the sightseeing :  max(c) st. A*x < b   and   A_eq*x = b_eq | ||||
| def init_ub_dist(landmarks: list, max_steps: int): | ||||
|     # Objective function coefficients. a*x1 + b*x2 + c*x3 + ... | ||||
|     c = [] | ||||
|     # Coefficients of inequality constraints (left-hand side) | ||||
|     A = [] | ||||
|     for i, spot1 in enumerate(landmarks) : | ||||
|         dist_table = [0]*len(landmarks) | ||||
|         c.append(-spot1.attractiveness) | ||||
|         for j, spot2 in enumerate(landmarks) : | ||||
|             dist_table[j] = manhattan_distance(spot1.loc, spot2.loc) | ||||
|         A.append(dist_table) | ||||
|     c = c*len(landmarks) | ||||
|     A_ub = [] | ||||
|     for line in A : | ||||
|         A_ub += line | ||||
|     return c, A_ub, [max_steps] | ||||
|  | ||||
| # Go through the landmarks and force the optimizer to use landmarks where attractiveness is set to -1 | ||||
| def respect_user_mustsee(landmarks: list, A_eq: list, b_eq: list) : | ||||
|     L = len(landmarks) | ||||
|     for i, elem in enumerate(landmarks) : | ||||
|         if elem.attractiveness == -1 : | ||||
|             l = [0]*L*L | ||||
|             if elem.name != "arrivée" : | ||||
|                 for j in range(L) : | ||||
|                     l[j +i*L] = 1 | ||||
|             else :                          # This ensures we go to goal | ||||
|                 for k in range(L-1) : | ||||
|                         l[k*L+L-1] = 1 | ||||
|             A_eq = np.vstack((A_eq,l)) | ||||
|             b_eq.append(1) | ||||
|     return A_eq, b_eq | ||||
|  | ||||
| # Computes the path length given path matrix (dist_table) and a result | ||||
| def path_length(P: list, resx: list) : | ||||
|     return np.dot(P, resx) | ||||
|  | ||||
| # Initialize all landmarks (+ start and goal). Order matters here | ||||
| landmarks = [] | ||||
| landmarks.append(landmark("départ", -1, (0, 0))) | ||||
| landmarks.append(landmark("concorde", -1, (5,5))) | ||||
| landmarks.append(landmark("tour eiffel", 99, (1,1)))                           # PUT IN JSON | ||||
| landmarks.append(landmark("arc de triomphe", 99, (2,3))) | ||||
| landmarks.append(landmark("louvre", 70, (4,2))) | ||||
| landmarks.append(landmark("montmartre", 20, (0,2))) | ||||
| landmarks.append(landmark("arrivée", -1, (0, 0))) | ||||
|  | ||||
|  | ||||
|  | ||||
| # CONSTRAINT TO RESPECT MAX NUMBER OF STEPS | ||||
| max_steps = 25 | ||||
|  | ||||
| # SET CONSTRAINTS FOR INEQUALITY | ||||
| c, A_ub, b_ub = init_ub_dist(landmarks, max_steps)              # Add the distances from each landmark to the other | ||||
| P = A_ub                                                        # store the paths for later. Needed to compute path length | ||||
| A_ub, b_ub = respect_number(landmarks, A_ub, b_ub)              # Respect max number of visits.  | ||||
|  | ||||
| # SET CONSTRAINTS FOR EQUALITY | ||||
| A_eq, b_eq = init_eq_not_stay(landmarks)                       # Force solution not to stay in same place | ||||
| A_eq, b_eq = respect_user_mustsee(landmarks, A_eq, b_eq)       # Check if there are user_defined must_see. Also takes care of start/goal | ||||
| A_eq, b_eq = break_sym(landmarks, A_eq, b_eq)                  # break the symmetry. Only use the upper diagonal values | ||||
| A_eq, b_eq = respect_order(landmarks, A_eq, b_eq)              # Respect order of visit (only works when max_steps is limiting factor) | ||||
|  | ||||
| # Bounds for variables (x can only be 0 or 1) | ||||
| x_bounds = [(0, 1)] * len(c) | ||||
|  | ||||
| # Solve linear programming problem | ||||
| res = linprog(c, A_ub=A_ub, b_ub=b_ub, A_eq=A_eq, b_eq = b_eq, bounds=x_bounds, method='highs', integrality=3) | ||||
|  | ||||
| # Raise error if no solution is found | ||||
| if not res.success : | ||||
|     raise ValueError("No solution has been found, please adapt your max steps") | ||||
|  | ||||
| # Print result | ||||
| print_res(res, P) | ||||
|  | ||||
|  | ||||
|  | ||||
|  | ||||
|  | ||||
|  | ||||
|  | ||||
|  | ||||
		Reference in New Issue
	
	Block a user