使用Gurobi创建具有最高边缘连通性的networkx.Graph

Question

我有以下图表

G

G

使用以下代码创建

import networkx as nx
import matplotlib.pyplot as plt

G = nx.hoffman_singleton_graph()
pos = nx.spring_layout(G)
nx.draw(G, pos=pos)
nx.draw_networkx_labels(G=G, pos=pos)
plt.show()

G

由 50 个节点组成。我只想包含 25 个节点。此外，我只想包含最大化节点

A

(=5) 和节点

B

(=20) 之间连接的节点（和边）。

我写了以下代码：

import numpy as np
import gurobipy as grb
from networkx.algorithms.connectivity import local_edge_connectivity


A = 5
B = 20
nodes = list(G.nodes)
n_nodes = len(nodes)
edges = nx.to_numpy_array(G, nodelist=nodes)
thresh_nodes = 25


model = grb.Model()
F = model.addMVar(n_nodes, vtype=grb.GRB.BINARY)

model.addConstr(F.sum() == thresh_nodes)
model.addConstr(F[nodes.index(A)] == 1)
model.addConstr(F[nodes.index(B)] == 1)


E = F * edges

model.setObjective(local_edge_connectivity(nx.from_numpy_array(A=E), A, B), grb.GRB.MAXIMIZE)
model.optimize()

这会导致错误，因为 nx.from_numpy_array() 无法处理

E

的数据类型。如何创建

np.ndarray

的 .X 值 ({0, 1}) 的临时

E

并使用它来确定解？

Answer 1

这可以通过 MIP（混合整数程序）来解决，就像您尝试做的那样。 下面是

pyomo

中的实现，我相信如果需要的话，它可以轻松地直接翻译为

gurobipy

。 解决得很快。

备注：

• 局部边缘连通性与最大单位流同义，因此我们可以将其构造为具有单位容量的最大流
• 我们需要增强无向图以使其成为有向图（本质上是将边加倍）
• 设置流量平衡和流量限制...
• 本地连通性为7，也可以通过将阈值降低到15来实现

代码：

import pyomo.environ as pyo
import networkx as nx
import matplotlib.pyplot as plt

G = nx.hoffman_singleton_graph()

edges = list(G.edges())
# augment the graph with reverse edges
edges.extend([(n1, n0) for (n0, n1) in edges])
node_edges = {n:[e for e in edges if e[0]==n or e[1]==n] for n in G.nodes()}

m = pyo.ConcreteModel('graph model')

# SETS
m.N = pyo.Set(initialize=list(G.nodes()))
m.E = pyo.Set(initialize=edges)

# VARS
m.active_node = pyo.Var(m.N, domain=pyo.Binary, doc='Node Active')
m.max_flow = pyo.Var(doc='Max Flow')
m.flow = pyo.Var(m.E, domain=pyo.NonNegativeReals, doc='flow')

s=5
t=20
threshold=25  # <- this could be cleaned up as we know s, t MUST be active

# OBJ: maximize flow
m.obj = pyo.Objective(expr=m.max_flow, sense=pyo.maximize)

# CONSTRAINTS
# 1.  Flow balance
@m.Constraint(m.N)
def flow_balance(m, n):
    inbounds = [e for e in m.E if e[1]==n]
    outbounds = [e for e in m.E if e[0]==n]
    if n == s:
        return m.max_flow <= sum(m.flow[e] for e in outbounds) - sum(m.flow[e] for e in inbounds)
    elif n == t:
        return m.max_flow <= sum(m.flow[e] for e in inbounds) - sum(m.flow[e] for e in outbounds)
    else:
        return sum(m.flow[e] for e in outbounds) == sum(m.flow[e] for e in inbounds)

# 2.  Limit Flow to selection
@m.Constraint(m.E)
def flow_limit_1(m, *e):
    return m.flow[e] <= m.active_node[e[1]]

@m.Constraint(m.E)
def flow_limit_2(m, *e):
    return m.flow[e] <= m.active_node[e[0]]

# 3.  Limit selections
m.node_limit = pyo.Constraint(expr=sum(m.active_node[n] for n in m.N) <= threshold)

# Check it...
m.pprint()

# solve it...
opt = pyo.SolverFactory('appsi_highs')
res = opt.solve(m, tee=True)

m.active_node.display()

print(f'Conectivity: {pyo.value(m.max_flow)}')

# prune it for display
active_nodes = [n for n in m.N if pyo.value(m.active_node[n])]
active_nodes.extend([s, t])
for node in m.N:
    if node not in active_nodes:
        G.remove_node(node)


pos = nx.spring_layout(G)
pos[s] = [-10,0]
pos[t] = [10,0]
nx.draw(G, pos=pos)
nx.draw_networkx_labels(G=G, pos=pos)
plt.show()

输出（部分）