Commit 9eb9696d authored by Michael Keller's avatar Michael Keller
Browse files

Savepoint

parent b4af023e
from cmath import sqrt
import random
import numpy as np
import gurobipy as gp
from gurobipy import GRB
np.set_printoptions(threshold=np.inf)
# Toy Example
# X = 1
# Y = 4
# Y = 6
# R = lambda a, b: 1 if a == b else 0
# dist = [2,2]
# dist = [1,1,1,1,1,1]
# X = 3
# Y = 3
# R = lambda a, b: -1 if abs(a - b) == 0 else 0
# dist = [4,5]
X = 1
X = 5
Y = 5
R = lambda a, b: 1 if abs(a - b) <= 1 else 0
dist = [1,1,1,1,1]
# X = 5
# Y = 5
# R = lambda a, b: [[-1, 0, 0.5, -0.5, 1, 0, 0.5, 1, 0.5, 0],
# [0, -1, 0, 1, 0, 0, 0.5, 0.5, 0.5, 0],
# [0.5, 0, -1, 0, 0, 0, 0.5, 0.5, 0, -0.5],
# [-0.5, 1, 0, -1, 0, 0, -1, 0.5, 0.5, 0.5],
# [0, 0, 0, 0, -1, 0, 0, 0, 0, 0],
# [0, 0, 0, 0, 0, -1, 1, 0.5, 0, 0],
# [0.5, 0.5, 0.5, -0.5, 0, 0.5, -1, 0.5, 0.5, 0.5],
# [0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
# [0.5, 0.5, 0, 1, 0, 0, 0.5, 0.5, -1, 0],
# [0, 0, -0.5, 0.5, 0, 0, 1, 0.5, 0.5, -1]][a][b]
# dist = [3, 3, 3, 3, 3, 2, 2, 2, 2, 2]
R = lambda a, b: [[-1, 0, 0.5, -0.5, 1, 0, 0.5, 1, 0.5, 0],
[0, -1, 0, 1, 0, 0, 0.5, 0.5, 0.5, 0],
[0.5, 0, -1, 0, 0, 0, 0.5, 0.5, 0, -0.5],
[-0.5, 1, 0, -1, 0, 0, -1, 0.5, 0.5, 0.5],
[0, 0, 0, 0, -1, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, -1, 1, 0.5, 0, 0],
[0.5, 0.5, 0.5, -0.5, 0, 0.5, -1, 0.5, 0.5, 0.5],
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
[0.5, 0.5, 0, 1, 0, 0, 0.5, 0.5, -1, 0],
[0, 0, -0.5, 0.5, 0, 0, 1, 0.5, 0.5, -1]][a][b]
dist = [2,2,2,2,2,3,3,3,3,3]
# auto calculate some vars
numb_crops = len(dist)
......@@ -72,6 +76,9 @@ def main():
if get_ind(pA, pB, cA, cB) is None:
obj.append(-(R(cA, cB) + R(cB, cA)))
if ind < 150:
print(ind, (pA, pB, cA, cB))
# memoization
dp[str(pA) + "|" + str(pB) + "|" + str(cA) + "|" + str(cB)] = ind
dp[str(pB) + "|" + str(pA) + "|" + str(cB) + "|" + str(cA)] = ind
......@@ -117,6 +124,21 @@ def main():
totals_one = np.array(totals_one)
totals_one_b = np.array(totals_one_b)
m.addConstr(totals_one @ x == totals_one_b, "totals one")
# all pixel amounts sum to 1
pixels_totals_one = []
pixels_totals_one_b = []
for pA in pixels:
pB = get_neighbors(pA).pop()
constr = [0 for _ in obj]
for cA in range(numb_crops):
for cB in range(numb_crops):
constr[get_ind(pA, pB, cA, cB)] = 1
pixels_totals_one.append(constr)
pixels_totals_one_b.append(1)
pixels_totals_one = np.array(pixels_totals_one)
pixels_totals_one_b = np.array(pixels_totals_one_b)
m.addConstr(pixels_totals_one @ x == pixels_totals_one_b, "pixels totals one")
# Enforce the global crop distribution
global_crops = []
......@@ -133,40 +155,99 @@ def main():
global_crops_b = np.array(global_crops_b)
m.addConstr(global_crops @ x == global_crops_b, "global crops")
# Quadratic constraint, the linear version (C = 2)
quad_lin = []
quad_lin_b = []
for pA in pixels:
neighbors = get_neighbors(pA)
for pB in neighbors:
for cA in range(numb_crops):
for cB in range(numb_crops):
constr = [0 for _ in obj]
if cA == cB:
continue
constr[get_ind(pA, pB, cA, cA)] = -1
constr[get_ind(pA, pB, cB, cB)] = -1
constr[get_ind(pA, pB, cA, cB)] = 1
constr[get_ind(pA, pB, cB, cA)] = 1
quad_lin.append(constr)
quad_lin_b.append(0)
quad_lin = np.array(quad_lin)
quad_lin_b = np.array(quad_lin_b)
m.addConstr(quad_lin @ x == quad_lin_b, "quad lin")
# quadratic constraint (relaxed)
# m.params.NonConvex = 2
# for pA in pixels:
# neighbors = get_neighbors(pA)
# for pB in neighbors:
# for cA in range(numb_crops):
# for cB in range(numb_crops):
# O = np.zeros((n, n))
# val = 0
# for cAn in range(numb_crops):
# for cBn in range(numb_crops):
# O[get_ind(pA, pB, cA, cAn)][get_ind(pB, pA, cB, cBn)] += 1
# # O[get_ind(pB, pA, cB, cBn)][get_ind(pA, pB, cA, cAn)] += 1
# val += x[get_ind(pA, pB, cA, cB)] + x[get_ind(pB, pA, cB, cBn)]
# print(O)
# print(np.all(np.linalg.eigvals(O) >= -0.00000000001))
# O = np.zeros((n,n))
# extra = 0
# for cAn in range(numb_crops):
# for cBn in range(numb_crops):
# O[get_ind(pA, pB, cA, cAn)][get_ind(pB, pA, cB, cBn)] += 1
# O[get_ind(pB, pA, cB, cBn)][get_ind(pA, pB, cA, cAn)] += 1
# m.addConstr(x @ O @ x == val, "quad constraint")
# # for i in range(n):
# # if O[i][i] == 0:
# # for j in range(n):
# # if O[i][j] > 0:
# # O[i][i] = 1
# # extra += x[i]
# # break
# print(O)
# O[4][4] = 1
# O[4][8] = 1
# O[8][4] = 1
# O[4][12] = 1
# O[12][4] = 1
# O[8][8] = 1
# O[12][12] = 1
# quadratic constraints (strict)
m.params.NonConvex = 2
for pA in pixels:
neighbors = get_neighbors(pA)
for pB in neighbors:
for cA in range(numb_crops):
for cB in range(numb_crops):
O = np.zeros((n,n))
for cAn in range(numb_crops):
for cBn in range(numb_crops):
O[get_ind(pA, pB, cA, cAn)][get_ind(pB, pA, cB, cBn)] = 1
# O[1][1] = 1
# O[1][2] = 1
# O[1][3] = 1
# O[2][1] = 1
# O[2][2] = 1
# O[2][3] = 1
# O[3][1] = 1
# O[3][2] = 1
# O[3][3] = 1
# print(O)
m.addConstr(x @ O @ x == x[get_ind(pA, pB, cA, cB)], "quad constraint")
# w, v = np.linalg.eigh(O)
# print(np.all(w >= -0.00000001))
# # print(v)
# # exit()
# m.addConstr(x @ O @ x <= 2 * x[get_ind(pA, pB, cA, cB)] + extra, "quad constraint")
# m.addConstr(extra <= 0.5, "quad constraint")
# quadratic constraints (strict)
# m.params.NonConvex = 2
# for pA in pixels:
# neighbors = get_neighbors(pA)
# for pB in neighbors:
# for cA in range(numb_crops):
# for cB in range(numb_crops):
# O = np.zeros((n,n))
# for cAn in range(numb_crops):
# for cBn in range(numb_crops):
# O[get_ind(pA, pB, cA, cAn)][get_ind(pB, pA, cB, cBn)] = 1
# # print(O)
# m.addConstr(x @ O @ x == x[get_ind(pA, pB, cA, cB)], "quad constraint")
# opt_sol = np.array([1,0,0,0,0,1,0,0,0,0,0,1])
# print("Optimal Score:", obj @ opt_sol)
# print(quad_lin @ opt_sol)
# print(opt_sol)
# print(quad_lin)
# ask the oracle
print("Ask Oracle")
......@@ -184,9 +265,15 @@ def main():
pB = get_neighbors(pA).pop()
for cB in range(numb_crops):
cA_amount += sol[get_ind(pA, pB, cA, cB)]
print(round(cA_amount, 2), end=", ")
print(cA_amount, end=", ")
print("")
print(sol)
# print(sol)
print(obj @ sol)
# for i in range(len(obj)):
# if obj[i] * sol[i] < 0:
# print(i, obj[i], sol[i])
# print(quad_lin @ sol)
......
from numpy import linalg as LA
import numpy as np
X = 5
Y = 5
R = lambda a, b: [[-1, 0, 0.5, -0.5, 1, 0, 0.5, 1, 0.5, 0],
[0, -1, 0, 1, 0, 0, 0.5, 0.5, 0.5, 0],
[0.5, 0, -1, 0, 0, 0, 0.5, 0.5, 0, -0.5],
[-0.5, 1, 0, -1, 0, 0, -1, 0.5, 0.5, 0.5],
[0, 0, 0, 0, -1, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, -1, 1, 0.5, 0, 0],
[0.5, 0.5, 0.5, -0.5, 0, 0.5, -1, 0.5, 0.5, 0.5],
[0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
[0.5, 0.5, 0, 1, 0, 0, 0.5, 0.5, -1, 0],
[0, 0, -0.5, 0.5, 0, 0, 1, 0.5, 0.5, -1]][a][b]
dist = [2,2,2,2,2,3,3,3,3,3]
# matrix = [
# [2,1,1,1,0,0,1,0,0],
# [1,1,1,1,0,0,1,0,0],
# [1,1,1,1,0,0,1,0,0],
# [1,1,1,1,0,0,1,0,0],
# [0,0,0,0,0,0,0,0,0],
# [0,0,0,0,0,0,0,0,0],
# [1,1,1,1,0,0,1,0,0],
# [0,0,0,0,0,0,0,0,0],
# [0,0,0,0,0,0,0,0,0]
# ]
matrix = [[0 for _ in range(len(dist))] for _ in range(len(dist))]
for x in range(len(dist)):
for y in range(len(dist)):
matrix[x][x] += (R(x, y) + R(y, x)) / 2 + 1
matrix[x][y] += (R(x, y) + R(y, x)) / 2 + 1
matrix[y][x] += (R(x, y) + R(y, x)) / 2 + 1
matrix[y][y] += (R(x, y) + R(y, x)) / 2 + 1
w, v = LA.eig(matrix)
print(w)
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment