Início / coding / Perguntas / 77937790
Accepted
Ryan
Asked: 2024-02-05 06:13:49 +0800 CST

Preso na formulação de uma restrição de adjacência com pyomo

  • 772

Estou tentando atribuir attributesuma matriz 3 x 3 com base em uma restrição de adjacência, mas estou preso na formulação da restrição de adjacência. Fico pensando que tenho a resposta, mas quando tento o resultado não é o esperado.

Portanto, elements são pré-organizados em uma element_mapmatriz 9 x 17. Cada um dos 9 elements tem 17 attributes, que são sempre zero ou um.

Os elements estão em posições fixas na matriz 3 x 3:

[[0, 3, 6
  1, 4, 7
  2, 5, 8]]

Por exemplo, o elementíndice 2at element_mapé fixado na posição (2,0)na matriz 3 x 3. E (2,0)é adjacente a (1,0), (1,1), (2,1).

As restrições são:

  1. cada um placepode ter no máximo 4 attributes atribuídos a ele ( feito )
  2. cada um placedeve ter um ou mais attributesatribuídos a ele ( feito )
  3. os s selecionados attributepara cada um placedevem ser iguais a zero ( concluído )
  4. Restrição de adjacência : explicada abaixo com um exemplo ( preciso de ajuda )
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)

Com um element_mapque se parece com isto:

[[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]]

... o resultado que obtenho é:

((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

O formato é: ((x,y), selected_attribute). (x,y)é a coordenada na matriz 3 x 3 e selected_attributeé um dos selecionados attributepara isso place.

Exemplo

Vamos considerar o elementat (1,1), pois ele é adjacente a todos os outros places. Os dados para isso elementsão:

[1 1 0 1 0 1 0 1 0 1 0 0 0 0 0 0 1]

Os attributes selecionados (1,1)no resultado são 4e 15. Dar uma olhada em element_map, attribute 4é uma boa escolha para this element, porque todos, exceto o penúltimo, elementtêm um valor 1para that attribute. O outro selecionado attribute, 15, não é bom, entretanto, porque o attributeíndice at 15é zero para o penúltimo elemento. Para elementat (1,1), attribute 4funcionaria em combinação com os seguintes attributeíndices: 3 or 9 or 16, porque para o penúltimo element, esses attributevalores são todos 1.

Definição de restrição de adjacência

Portanto, a restrição de adjacência deveria ser: atribuir attributes a a placetal que todos os places adjacentes tenham pelo menos um valor 1entre todos os attributes selecionados. E é um problema de minimização, portanto, no geral, quero o número mínimo de attributes atribuídos a cada um placeque atenda às restrições.

Espero que tenha ficado claro, comente se houver algo que eu possa esclarecer; obrigado por ler!

python

1 respostas

  1. Best Answer

    Você precisa de alguns componentes....

    Primeiro, você precisa de uma lista indexada de quais locais são adjacentes a cada local. Se você estiver usando seu sistema de indexação de um dígito, seria algo como:

    neighbors = { 0: [1, 3, 4], 1: [0, 3, 4, 5, 2], ...}

    estes índices serão a base da soma de uma restrição para garantir a cobertura. Algo como em pseudocódigo/não testado:

    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)

    Basicamente afirmando que para cada combo base-vizinho, você deve selecionar algum atributo placeque gere cobertura.

    Tente ver se você consegue fazer isso girar. Se travar, comente de volta!

    • 1

relate perguntas