Pyscipopt是scip求解器的Python库，可以通过Python调用Pyscipopt中的函数、建模、求解。
scip求解器介绍
SCIP是用于混合整数规划和混合整数非线性规划的最快的非商业求解器之一。SCIP是一个混合整数规划求解器和一个用于分支定界以及分支定价的框架。由德国柏林Zuse 研究所（ZIB）团队开发维护。
python中使用scip
python中可以直接通过安装pyscipopt库来使用scip。
安装pyscipopt
与安装其他库一样，使用pip直接安装。
pip install pyscipopt
基础使用
pyscipopt可以直接进行建模，非常方便，先给一个基本的案例。
from pyscipopt import Model#导入包
model = Model("Example") # 命名模型
x = model.addVar("x")
y = model.addVar("y")
model.setObjective(x+y)#目标值
model.addCons(x>= 1)
model.addCons(y>= 1)
model.addCons(2*x - x*y >= 2)
model.addCons(3*x + y*y <= 10)
model.optimize()#默认最小值
sol = model.getBestSol()
print("x: {}".format(sol[x]))
print("y: {}".format(sol[y]))
tsp问题
tsp问题是一个np-hard问题，也是运筹学里很经典的一个问题，这里使用scip也可以通过建模求解。
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx
import time
import pyscipopt as opt
def print_fun_run_time(func):
def wrapper(*args, **kwargs):
local_time = time.time()
func(*args, **kwargs)
print('current function [%s] run time is %.4f'%(func.__name__, time.time()-local_time))
return wrapper
def prepare_data():
# 测试数据: 邻接矩阵标识图结构，其中极大值1000表示不连通, 测试数据包含11个节点
matrix = np.array([[0, 2, 8, 1, 1000, 1000, 1000, 1000, 1000, 1000, 1000],
[2, 0, 6, 1000, 1, 1000, 1000, 1000, 1000, 1000, 1000],
[8, 6, 0, 7, 5, 1, 2, 1000, 1000, 1000, 1000],
[1, 1000, 7, 0, 1000, 1000, 9, 1000, 1000, 1000, 1000],
[1000, 1, 5, 1000, 0, 3, 1000, 2, 1000, 1000, 1000],
[1000, 1000, 1, 1000, 3, 0, 4, 1000, 6, 1000, 1000],
[1000, 1000, 2, 9, 1000, 4, 0, 1000, 3, 1, 1000],
[1000, 1000, 1000, 1000, 2, 1000, 1000, 0, 7, 1000, 9],
[1000, 1000, 1000, 1000, 1000, 6, 3, 7, 0, 1, 2],
[1000, 1000, 1000, 1000, 1000, 1000, 1, 1000, 1, 0, 4],
[1000, 1000, 1000, 1000, 1000, 1000, 1000, 9, 2, 4, 0]])
# 转换数据（只是为了方便理解，不为效率）, 获取不包含本节点且连通的点及弧
row, col = np.where((matrix != 1000) & (matrix != 0))
Nodes = list(np.unique(np.unique(np.append(row, col)) + 1).astype('str'))
Arcs = {}
for i in range(len(row)):
Arcs[str(row[i] + 1), str(col[i] + 1)] = matrix[row[i], col[i]]
print('Nodes:\n', Nodes)
print('Arcs:\n', Arcs)
return Nodes, Arcs
@print_fun_run_time
def SCIP_solver(Nodes, Arcs, source, target):
model = opt.Model('SPP')
# ==========定义变量==========
arc = {}
for idx in Arcs.keys():
arc[idx] = model.addVar(vtype='B', name=str(idx))
# ==========定义约束==========
# 其实类似VRP的流量约束，源节点一定是流出，结束节点一定是流入，其他节点有流入一定有流出
# 约束1: 源点一定流出
model.addCons(opt.quicksum(arc[idx] for idx in Arcs.keys() if idx[0] == source) == 1, name='source_flow')
# 约束2: 终点一定流入
model.addCons(opt.quicksum(arc[idx] for idx in Arcs.keys() if idx[1] == target) == 1, name='target_flow')
# 约束3: 其他点流平衡
for node in Nodes:
if node != source and node != target:
model.addCons(opt.quicksum(arc[idx] for idx in Arcs.keys() if idx[0] == node and idx[1] != source) -
opt.quicksum(arc[idx] for idx in Arcs.keys() if idx[1] == node and idx[0] != target) == 0)
# ==========定义目标==========
model.setObjective(opt.quicksum(arc[arc_idx] * arc_weight for arc_idx, arc_weight in Arcs.items()))
model.setMinimize()
model.optimize()
# ==========输出结果==========
distance = model.getObjVal()
# 转换结果
visit_arcs = []
for idx in Arcs.keys():
if model.getVal(arc[idx]) > 0:
visit_arcs.append(idx)
visit_arcs_dict = dict([list(i) for i in visit_arcs])
idx = source
path = []
for i in range(len(Nodes)):
next_node = visit_arcs_dict.get(idx, None)
if next_node is not None:
path.append(next_node)
idx = next_node
print(f'节点 {source} 到 {target} 的路径={path}, 最短距离={distance}')
return True
if __name__ == '__main__':
source = '1' # starting node
target = '11' # ending node
# 测试数据
Nodes, Arcs = prepare_data()
# 方式3: python调用开源求解器SCIP求解最短路径问题
SCIP_solver(Nodes, Arcs, source, target)
基础设施选址问题
官网提供了一个基础设施选址问题的代码案例。
from pyscipopt import Model, quicksum, multidict
def flp(I,J,d,M,f,c):
"""flp -- model for the capacitated facility location problem
Parameters:
- I: set of customers
- J: set of facilities
- d[i]: demand for customer i
- M[j]: capacity of facility j
- f[j]: fixed cost for using a facility in point j
- c[i,j]: unit cost of servicing demand point i from facility j
Returns a model, ready to be solved.
"""
model = Model("flp")
x,y = {},{}
for j in J:
y[j] = model.addVar(vtype="B", name="y(%s)"%j)
for i in I:
x[i,j] = model.addVar(vtype="C", name="x(%s,%s)"%(i,j))
for i in I:
model.addCons(quicksum(x[i,j] for j in J) == d[i], "Demand(%s)"%i)
for j in M:
model.addCons(quicksum(x[i,j] for i in I) <= M[j]*y[j], "Capacity(%s)"%i)
for (i,j) in x:
model.addCons(x[i,j] <= d[i]*y[j], "Strong(%s,%s)"%(i,j))
model.setObjective(
quicksum(f[j]*y[j] for j in J) +
quicksum(c[i,j]*x[i,j] for i in I for j in J),
"minimize")
model.data = x,y
return model
def make_data():
"""creates example data set"""
I,d = multidict({1:80, 2:270, 3:250, 4:160, 5:180}) # demand
J,M,f = multidict({1:[500,1000], 2:[500,1000], 3:[500,1000]}) # capacity, fixed costs
c = {(1,1):4, (1,2):6, (1,3):9, # transportation costs
(2,1):5, (2,2):4, (2,3):7,
(3,1):6, (3,2):3, (3,3):4,
(4,1):8, (4,2):5, (4,3):3,
(5,1):10, (5,2):8, (5,3):4,
}
return I,J,d,M,f,c
if __name__ == "__main__":
I,J,d,M,f,c = make_data()
model = flp(I,J,d,M,f,c)
model.optimize()
EPS = 1.e-6
x,y = model.data
edges = [(i,j) for (i,j) in x if model.getVal(x[i,j]) > EPS]
facilities = [j for j in y if model.getVal(y[j]) > EPS]
print("Optimal value:", model.getObjVal())
print("Facilities at nodes:", facilities)
print("Edges:", edges)
try: # plot the result using networkx and matplotlib
import networkx as NX
import matplotlib.pyplot as P
P.clf()
G = NX.Graph()
other = [j for j in y if j not in facilities]
customers = ["c%s"%i for i in d]
G.add_nodes_from(facilities)
G.add_nodes_from(other)
G.add_nodes_from(customers)
for (i,j) in edges:
G.add_edge("c%s"%i,j)
position = NX.drawing.layout.spring_layout(G)
NX.draw(G,position,node_color="y",nodelist=facilities)
NX.draw(G,position,node_color="g",nodelist=other)
NX.draw(G,position,node_color="b",nodelist=customers)
P.show()
except ImportError:
print("install 'networkx' and 'matplotlib' for plotting")
pyscipopt在建模和求解方面表现还是不错的，虽然比gurobi、cplex稍逊一筹，但scip开源免费，也是一个非常良心的求解器了。