Блог

Что я сделал

Во-первых, это может быть не самый лучший форум, так что приношу свои извинения, если это так. Я создаю модель Pyomo, в которой я создал двоичную матрицу следующим образом:

 model.binMat = Var(range(6),range(6),domain=Binary)
  

Моя модель решает для этой матрицы с типичным выводом, подобным этому:

 binaryMatrix =  [[0 1 0 1 0 0]
                 [1 0 1 0 0 0]
                 [0 1 0 0 0 1]
                 [1 0 0 0 1 0]
                 [0 0 0 1 0 1]
                 [0 0 1 0 1 0]]
  

Результаты интерпретируются как координаты единиц, т.е. (1,2),(1,4),(2,1),(2,3),(3,2),(3,6),(4,1),(4,5),(5,4),(5,6),(6,3),(6,5) в этот пример.

Затем это рассматривается в терминах групп подключенных элементов. В этом случае будет только 1 уникальная группа: (1,2,3,4,5,6).

Что мне нужно

Я хотел бы получить помощь в создании нового ограничения, позволяющего использовать только 2 уникальные группы одинакового размера, ссылаясь на значения в model.binMat.

Примером того, как могли бы выглядеть эти конечные группы, являются: (1,5,6) и (2,3,4). Соответствующими координатами для этого могут быть: (1,5),(1,6),(2,3),(2,4),(3,2),(3,4),(4,2),(4,3),(5,1),(5,6),(6,1),(6,5)

В настоящее время я пытаюсь решить эту проблему с помощью наборов Pyomo, но поскольку они для меня новы, мне не повезло.

Редактировать

Для тех, кто интересуется альтернативными подходами к той же проблеме, я также разместил это здесь

Может быть более простой способ, но лучший способ, который я мог придумать, — это добавить двоичные ограничения для проверки каждого возможного такого набора и принудительно выбрать один из этих наборов уникальных компонентов одинакового размера. Обратите внимание, что этот подход приводит к экспоненциальному числу ограничений, поэтому он не является хорошим решением для более крупных проблем.

 import pyomo.environ as pyo
import itertools

nodes = set(range(6))
# the possible sets of components of length 3
full_comp_list = [(set(i),nodes-set(i)) for i in itertools.combinations(nodes,3)]
# only take the first half because it's symmetric with six nodes and equal size
comp_list = full_comp_list[:int(len(full_comp_list)/2)]

num_poss_component_sets = len(comp_list)

#%% Build model
model = pyo.ConcreteModel()
model.binMat = pyo.Var(nodes,nodes,domain=pyo.Binary)

#%% Additional Variables
# binaries to track if each component connected
model.comp1_connected= pyo.Var(range(num_poss_component_sets),within=pyo.Binary)
model.comp2_connected= pyo.Var(range(num_poss_component_sets),within=pyo.Binary)
# tracks if the two components are disjoint
model.comps_disjoint = pyo.Var(range(num_poss_component_sets),within=pyo.Binary)
# tracks if the criteria met for any set of components
model.meet_criteria = pyo.Var(range(num_poss_component_sets),within=pyo.Binary)

#%% Additional constraints
def is_comp1_connected_rule(model,comp_num):
    ''' The component is complete iff the number of (directed) edges between ==6 (all three undirected edges selected)'''
    return(sum(model.binMat[i,j] for i,j in itertools.combinations(comp_list[comp_num][0],2))
    >=3*model.comp1_connected[comp_num])
   
model.comp1_connected_constraint = pyo.Constraint(range(num_poss_component_sets),
                                                  rule=is_comp1_connected_rule)

# Check if each component set is a complete graph
def is_comp2_connected_rule(model,comp_num):
    ''' The component is complete iff the number of (directed) edges between == 6 (all three undirected edges selected)'''
    return(sum(model.binMat[i,j] for i,j in itertools.combinations(comp_list[comp_num][1],2))
    >= 3*model.comp2_connected[comp_num])
   
model.comp2_connected_constraint = pyo.Constraint(range(num_poss_component_sets),
                                                  rule=is_comp2_connected_rule)

# Check if components are separate from each other (no edges between)
def are_both_disjoint_rule(model,comp_num):
    '''Disjoint if no edges between any nodes in different component
    If there are ANY edges between, then not disjoint (model.both_comps_connected must be 0)
    '''
    return(sum([model.binMat[i,j] for i in comp_list[comp_num][0] for j in comp_list[comp_num][1]])
    <= 9 * (1-model.comps_disjoint[comp_num]))
   
model.comps_disjoint_constraint = pyo.Constraint(range(num_poss_component_sets),
                                                      rule=are_both_disjoint_rule)

# Determines if a given set of components meet the rule
def meet_criteria_rule(model,comp_num):
    '''Rule met if both components are connected and separate from each other'''
    return(model.comp1_connected[comp_num]   model.comp2_connected[comp_num]
      model.comps_disjoint[comp_num] >= 3 * model.meet_criteria[comp_num])
   
model.comp_meets_criteria_constraint = pyo.Constraint(range(num_poss_component_sets),
                                                rule=meet_criteria_rule)

# at least one component must meet rule that theyre separate connected components
model.must_meet_criteria_constraint = pyo.Constraint(expr = sum(model.meet_criteria[comp_num]
for comp_num in range(num_poss_component_sets)) >= 1)

### New constraint to make adjacency matrix symmetric (binMat_{i,j} == binMat_{j,i})
def edges_symmetric_rule(model,node1,node2):
    '''Rule requiring both directions for edges to be the same'''
    return(model.binMat[node1,node2] == model.binMat[node2,node1])
model.edges_symmetric_constraint = pyo.Constraint(nodes,nodes,rule=edges_symmetric_rule)

#%% Add objective and solve
des_edges = [(4,0),(1,2),(1,3),(2,1),(2,3),(3,1),(3,2)]
pos_c_dict = {e:1 for e in des_edges}
c = [[pos_c_dict.get((i,j),-1) for i in nodes] for j in nodes]
model.obj = pyo.Objective(expr = sum([c[i][j]*model.binMat[i,j] for i in nodes for j in nodes]),
                          sense=pyo.maximize)

solver = pyo.SolverFactory('glpk')
res = solver.solve(model)

# get the components and the index for what's chosen
[comp_list[i] for i in range(len(comp_list)) if pyo.value(model.meet_criteria[i])]
# [({0, 4, 5}, {1, 2, 3})]
[i for i in range(len(comp_list)) if pyo.value(model.meet_criteria[i])]
# 9

# View the final binMat
final_binMat = pd.DataFrame({'source':list(nodes)*len(nodes),
                             'target':[j for i in nodes for j in [i]*len(nodes)]})
final_binMat['value'] = [pyo.value(model.binMat[i,j]) for i,j in final_binMat.values]
final_binMat['cost'] = [c[i][j] for i,j in final_binMat[['source','target']].values]
final_binMat_wide = pd.pivot(data=final_binMat,index='source',columns='target',values='value')

# target    0    1    2    3    4    5
# source                              
# 0       0.0  0.0  0.0  0.0  1.0  1.0
# 1       0.0  0.0  1.0  1.0  0.0  0.0
# 2       0.0  1.0  0.0  1.0  0.0  0.0
# 3       0.0  1.0  1.0  0.0  0.0  0.0
# 4       1.0  0.0  0.0  0.0  0.0  1.0
# 5       1.0  0.0  0.0  0.0  1.0  0.0