问题描述
我做了什么
首先,这可能不是最好的论坛,如果是这种情况,我们深表歉意。我正在创建一个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的坐标,即(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)。
我需要什么
我希望通过参考model.binMat中的值来创建一个新约束,以仅允许 2个大小相同的唯一组。
这些最终组的外观示例如下:(1,6)和(2,4)。为此的相应坐标可以是:(1,5),(1,6),(2,3),(2,4),(3,4),(4,2), (4,3),(5,1),(5,1),(6,5)
我目前正在尝试使用Pyomo套件解决此问题,但是由于这些对我来说是新的,所以我没有任何运气。
编辑
对于那些对解决同一问题的替代方法感兴趣的人,我还发布了此here
解决方法
也许有一种更简单的方法,但是我能想到的最好方法是添加二进制约束来检查每个可能的这样的集合,并强制选择那些大小相等的唯一组件集。请注意,这种方法会导致成倍的约束,因此对于较大的问题并不是一个好的解决方案。
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 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,rule=edges_symmetric_rule)
#%% Add objective and solve
des_edges = [(4,0),(1,2),3),(2,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