I'm trying to assign attributes into a 3 x 3 matrix based on an adjacency constraint, but I'm stuck formulating the adjacency constraint. I keep thinking I've got the answer, then when I try it the result isn't as expected.
So, elements are pre-arranged into an element_map which is a 9 x 17 matrix. Each of the 9 elements has 17 attributes, which are always zero or one.
The elements are in fixed positions in the 3 x 3 matrix:
[[0, 3, 6
1, 4, 7
2, 5, 8]]
For example, the element at index 2 in the element_map is fixed at position (2,0) in the 3 x 3 matrix. And (2,0) is adjacent to (1,0), (1,1), (2,1).
The constraints are:
place can have at most 4 attributes assigned to it (done)place must have one or more attributes assigned to it (done)attributes for each place must be equal to zero (done)import pyomo.environ as pyo
import numpy as np
"""
fit elements into matrix based on adjacency rules
"""
class Element:
"""a convenience to hold the rows of attribute values"""
def __init__(self, row):
self.attributes = tuple(row)
def attribute(self, idx):
return self.attributes[idx]
def __repr__(self):
return str(self.attributes)
class Position:
"""a convenience for (x, y) positions that must have equality & hash defined for consistency"""
def __init__(self, x, y):
self.x = x
self.y = y
def __repr__(self):
return f'({self.x}, {self.y})'
def __hash__(self):
return hash((self.x, self.y))
def __eq__(self, other):
if isinstance(other, Position):
return (self.x, self.y) == (other.x, other.y)
return False
# each 'row' corresponds to an element
# each 'column' corresponds to an attribute of the various elements
# here, each element has attributes which are always 0 or 1
element_map = np.array([[0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0],
[1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1],
[0, 0, 0, 1, 1, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 1],
[0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1],
[1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1],
[1, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1],
[0, 0, 0, 0, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1],
[0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1],
[1, 0, 0, 1, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 1]])
print(element_map, 'map')
print (element_map.shape)
matrix_a_rows = 3
matrix_a_cols = 3
matrix_a = np.zeros((matrix_a_rows, matrix_a_cols))
mask = np.arange(1,(matrix_a_rows*matrix_a_cols)+1).reshape(matrix_a_cols, matrix_a_rows).T
def get_element(position):
#get the element vector at a position
x,y =position.x, position.y
idx = mask[x,y]-1
return idx #return 0-indexed element idx
def adj_xy(mat, p: Position):
x, y = p.x, p.y
res = []
rows = len(mat) - 1
cols = len(mat[0]) - 1
for i in range(x - 1, x + 2):
for j in range(y - 1, y + 2):
if all((0 <= i <= rows, 0 <= j <= cols, (i, j) != (x, y))):
res.append(Position(i, j))
return res
# SET UP ILP
m = pyo.ConcreteModel('matrix_fitter')
# SETS
elements = np.array([np.array(row) for row in element_map])
m.P = pyo.Set(initialize=[Position(x, y) for x in range(len(matrix_a)) for y in range(len(matrix_a[0]))],
doc='positions')
m.A = pyo.Set(initialize=list(range(len(element_map[0]))), doc='attribute')
# VARS
# place element e in position p based on attribute a being 0...
m.place = pyo.Var(m.P, m.A, domain=pyo.Binary, doc='place')
# OBJ: minimize attributes assigned to each position
m.obj = pyo.Objective(expr=pyo.sum_product(m.place), sense=pyo.minimize)
#each place must have 4 or fewer attributes assigned to it
m.attr_constraint = pyo.ConstraintList()
for p in m.P:
s = 0
for a in m.A:
s += m.place[p,a]
m.attr_constraint.add(s <= 4)
#each place must have one or more attributes assigned to it
m.enzyme_constraint = pyo.ConstraintList()
for p in m.P:
s = 0
for a in m.A:
s += m.place[p,a]
m.attr_constraint.add(s >= 1)
#the selected attribute for each position must be equal to zero
m.cut_constraint = pyo.ConstraintList()
for p in m.P:
e_idx = get_element(p)
element = elements[e_idx]
for i,a in enumerate(m.A):
if element[i] == 1:
m.cut_constraint.add(m.place[p,a] == 0)
#adjacency constraint
#doesn't work as expected
#This is where I need help...
m.adjacency_constraint = pyo.ConstraintList()
for p in m.P:
plate_element = elements[get_element(p)]
neighbor_positions = adj_xy(matrix_a, p)
print (p, 'place')
for i,a in enumerate(m.A):
s = 0
for j,aa in enumerate(m.A):
for neighbor in neighbor_positions:
neighbor_element = elements[get_element(neighbor)]
neighbor_element_attribute = neighbor_element[i]
s += m.place[neighbor,aa]*neighbor_element_attribute
m.adjacency_constraint.add(s >= len(neighbor_positions)*m.place[p,a])
solver = pyo.SolverFactory('cbc')
results = solver.solve(m, tee=True)
print(results, 'results')
if results.solver.termination_condition == pyo.TerminationCondition.optimal:
for idx in m.place.index_set():
if m.place[idx].value == 1:
s = 0
att = idx[1]
neighs = adj_xy(matrix_a, idx[0])
for i in neighs:
place_element = elements[get_element(i)]
s += place_element[att]
print(idx, 'idx', s)
if pyo.value(m.obj) == matrix_a_rows * matrix_a_cols:
# all positions were filled
print('success!')
else:
print(f'the number of elements that can be placed is {pyo.value(m.obj)} / {matrix_a_rows * matrix_a_cols}')
else:
print('Problem with model...see results')
print (element_map)
With an element_map that looks like this:
[[0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 1 0]
[1 0 0 0 1 0 0 1 0 1 0 0 0 0 0 1 1]
[0 0 0 1 1 0 0 1 0 1 0 1 0 0 0 1 1]
[0 1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 1]
[1 1 0 1 0 1 0 1 0 1 0 0 0 0 0 0 1]
[1 0 0 0 1 0 0 1 0 0 0 0 0 0 0 0 1]
[0 0 0 0 1 0 1 1 0 1 0 0 0 0 0 0 1]
[0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 0 1]
[1 0 0 1 1 1 0 0 0 1 0 0 0 0 0 1 1]]
...the result I get is:
((0, 0), 0) idx 2
((0, 0), 1) idx 2
((0, 0), 9) idx 2
((0, 0), 16) idx 3
((0, 1), 5) idx 2
((0, 1), 6) idx 2
((0, 1), 9) idx 4
((0, 2), 3) idx 2
((1, 0), 1) idx 2
((1, 0), 5) idx 2
((1, 1), 4) idx 7
((1, 1), 15) idx 4
((1, 2), 1) idx 2
((2, 0), 0) idx 3
((2, 1), 3) idx 4
((2, 2), 7) idx 2
The format is: ((x,y), selected_attribute). (x,y) is the coordinate in the 3 x 3 matrix and selected_attribute is one of the selected attributes for that place.
Example
Let's consider the element at (1,1), because it is adjacent to all other places. The data for this element is:
[1 1 0 1 0 1 0 1 0 1 0 0 0 0 0 0 1]
The selected attributes for (1,1) in the result are 4 and 15. Taking a look at the element_map, attribute 4 is a good choice for this element, because all except for the second-to-last element have a value of 1 for that attribute. The other selected attribute, 15, is not good, however, because the attribute at index 15 is a zero for the second-to-last element. For the element at (1,1), attribute 4 would work in combination with the following attribute indices: 3 or 9 or 16, because for the second-to-last element, those attribute values are all 1.
Adjacency Constraint Definition
So the adjacency constraint should be: assign attributes to a place such that all adjacent places have at least one value of 1 among all the selected attributes. And it is a minimization problem, so overall I want the minimum number of attributes assigned to each place that meet the constraints.
I hope this was clear, please comment if there's something I can clarify; thanks for reading!
You need a couple of components....
First, you need an indexed list of which places are adjacent to each place. If you are using your single-digit indexing system, that would be something like:
neighbors = { 0: [1, 3, 4], 1: [0, 3, 4, 5, 2], ...}
these indices will be the basis of summation for a constraint to ensure coverage. Something like in pseudocode/untested:
def cover(model, base, neighbor):
return sum(model.place[base, a] * element_map[neighbor, a] for a in attributes) >= 1
# make some flat list of all the possible base-neighbor combos...
bn = [(b, n) for b in neighbors for n in neighbors[b] ]
model.neighbor_coverage = pyo.Constraint(bn, rule=cover)
Basically stating that for each base-neighbor combo, you must select some attribute with place that generates coverage.
Try to see if you can get that spinning. If stuck, comment back!