前言
Gurobi是一款功能强大的商用求解器,支持Python、C、C++、Java等多种语言调用,相比于Cplex,Gurobi封装更高,更加方便,但是对于初学者而言会更难一些。Gurobi与Cplex存在兼容,Gurobi生成的mps或者lp文件可以在Cplex中运行。另外,在校学生和教师可以免费使用Gurobi的学术版,求解规模和求解速度不受限制。本文将简单记录python调用Gurobi配置及其求解的过程。
安装
用户可以通过gurobi官网下载并安装软件:gurobi官网,下载软件后学生可通过学生认证获取免费证书,获取证书后在电脑上验证即可。获取license后,在终端运行:grbgetkey xxxxxxxx,如下图所示。
运行完成后会得到gurobi.lic文件,可以放在指定目录下。之后在anaconda对应的环境中安装gurobi包,安装语句为:
conda install -c gurobi gurobi
安装完成后,每次运行程序需要激活许可证,即上面提到的gurobi.lic文件,输入以下语句:
export GRB_LICENSE_FILE=/Users/liuxinwei/gurobi.lic
许可证验证完成后即可调用gurobi包进行数学规划等操作。
简单示例(一)
'''
线性规划模型
max x+y+2z
s.t.
x+2y+3z<=4
x+y>=1
x,y,z=0 or 1
'''
from gurobipy import *
try:
#模型
model=Model('mip')
#变量
x=model.addVar(vtype=GRB.BINARY,name='x')
y=model.addVar(vtype=GRB.BINARY,name='y')
z=model.addVar(vtype=GRB.BINARY,name='z')
#目标函数
model.setObjective(x+y+2*z,GRB.MAXIMIZE)
#约束
model.addConstr(x+2*y+3*z<=4,name='c1')
model.addConstr(x+y>=1,name='c2')
#求解
model.setParam('outPutFlag',0)#不输出求解日志
model.optimize()
#输出
print('obj=',model.objVal)
for v in model.getVars():
print(v.varName,':',v.x)
except GurobiError as e:
print('Error code '+str(e.errno)+':'+str(e))
except AttributeError:
print('Encountered an attribute error')
简单示例(二)
'''
营养配方模型
目标:以最低的价格满足人体所需营养
约束:人体所需营养必须有上限和下限,每种食物有营养成分和价格的属性
'''
from gurobipy import *
#数据
categories,minNutrition,maxNutrition=multidict({
'calories':[1800,2200],
'protein':[91,GRB.INFINITY],
'fat':[0,65],
'sodium':[0,1779]
})
foods,cost=multidict({
'hamburger':2.49,
'chicken':2.89,
'hot dog':1.50,
'fries':1.89,
'macaroni':2.09,
'pizza':1.99,
'salad':2.49,
'milk':0.89,
'ice cream':1.59
})
nutritionValues = {
('hamburger', 'calories'): 410,
('hamburger', 'protein'): 24,
('hamburger', 'fat'): 26,
('hamburger', 'sodium'): 730,
('chicken', 'calories'): 420,
('chicken', 'protein'): 32,
('chicken', 'fat'): 10,
('chicken', 'sodium'): 1190,
('hot dog', 'calories'): 560,
('hot dog', 'protein'): 20,
('hot dog', 'fat'): 32,
('hot dog', 'sodium'): 1800,
('fries', 'calories'): 380,
('fries', 'protein'): 4,
('fries', 'fat'): 19,
('fries', 'sodium'): 270,
('macaroni', 'calories'): 320,
('macaroni', 'protein'): 12,
('macaroni', 'fat'): 10,
('macaroni', 'sodium'): 930,
('pizza', 'calories'): 320,
('pizza', 'protein'): 15,
('pizza', 'fat'): 12,
('pizza', 'sodium'): 820,
('salad', 'calories'): 320,
('salad', 'protein'): 31,
('salad', 'fat'): 12,
('salad', 'sodium'): 1230,
('milk', 'calories'): 100,
('milk', 'protein'): 8,
('milk', 'fat'): 2.5,
('milk', 'sodium'): 125,
('ice cream', 'calories'): 330,
('ice cream', 'protein'): 8,
('ice cream', 'fat'): 10,
('ice cream', 'sodium'): 180}
#模型
model=Model('diet')
#变量
buy=model.addVars(foods,name='buy')
#目标函数
model.setObjective(buy.prod(cost),GRB.MINIMIZE)
#model.setObjective(sum(buy[f]*cost[f] for f in foods),GRB.MINIMIZE)
#约束(expr==[lb,ub]-->lb<=expr<=ub)
model.addConstrs((quicksum(nutritionValues[f,c]*buy[f] for f in foods)\
==[minNutrition[c],maxNutrition[c]] for c in categories),'con')
#求解
model.setParam('OutputFlag', 0)
model.optimize()
#输出
if model.status==GRB.Status.OPTIMAL:
print('obj=',model.objVal)
buyx=model.getAttr('x',buy)#取buy的x属性
for f in foods:
if buyx[f]>0.0001:
print(f,':',buyx[f])
简单示例(三)
'''
指派模型
目标:成本最小化,时间最小化
'''
from gurobipy import *
import random
#数据
N=3
random.seed(1)
T={(i,j):random.randint(0,100) for i in range(N) for j in range(N)}
C={(i,j):random.randint(0,100) for i in range(N) for j in range(N)}
#模型
model=Model('assignment')
#添加变量
x=model.addVars([(i,j) for i in range(N) for j in range(N)],vtype='B',name='x')
#目标函数
model.setObjectiveN(x.prod(T),index=0,weight=0.1,name='Time Objective')
model.setObjectiveN(x.prod(C),index=1,weight=0.5,name='Cost Objective')
#约束
model.addConstrs(((quicksum(x[i,j] for i in range(N))==1) for j in range(N)),name='Job Constraingt')
model.addConstrs(((quicksum(x[i,j] for j in range(N))==1) for i in range(N)),name='Worker Constraingt')
#求解
model.Params.OutputFlag=0
model.optimize()
#输出
for i in range(model.NumObj):
model.setParam('ObjNumber',i)
print('Obj{}:'.format(i),model.objVal)
for i in range(N):
for j in range(N):
if x[i,j].x>0:
print('job {}-->worker {}'.format(j+1,i+1))
model.write('assignment.lp')
Callback函数调用
'''
callback函数调用
'''
from gurobipy import *
import random
#callback函数
def RINScallback(model,where):
if where==GRB.Callback.MIPNODE:
#节点数量%100==0时,调用callback函数
if model.cbGET(GRB.Callback.MINPNODE_NODCNT)%100==0 and\
model.cbGET(GRB.Callback.MIPNODE.STATUS)==GRB>OPTIMAL:
submodel=model.copy()
suby=submodel.getVars()
#获得节点松弛解
yrelaxation=model.cbGetNodeRel(model._y)
#固定变量取值
for i in range(model._N):
if abs(yrelaxation[i])<0.01:
suby[i].ub=0
elif abs(yrelaxation-1)<=0.01:
suby[i].lb=1
submodel.setParam(GRB.Param.TimeLimit,30)
submodel.optimize()
if submodel.objval>model.cbGet(GRB.Callback.MIPNODE_OBJBST):
#将解传递给原模型
for i in range(model._N):
if abs(suby[i].x)<=0.001:
model.cbSetSolution(model._y[i],0.0)
elif abs(suby[i].x-1)<=0.001:
model.cbSetSolution(model._y[i],1.0)
#数据
random.seed(1)
N=10
C={(i,j):random.randint(0,100) for i in range(N) for j in range(N)}
#模型
model=Model('MaxCut')
#变量
y=model.addVars(N,vtype='B',name='y')
#目标函数
obj=QuadExpr()
for i in range(N):
for j in range(N):
obj+=C[i,j]*(2*y[i]-1)*(2*y[j]-1)
model.setObjective(obj,GRB.MAXIMIZE)
#设置求解时间
model.Params.TimeLimit=600
#外部变量
model._y=y
model._N=N
#求解
model.optimize(RINScallback)
广义约束及其线性化
#模型
model=Model('gurobi')
#变量
x=model.addVar(vtype='C',name='x')
y=model.addVar(vtype='C',name='y')
z=model.addVar(vtype='C',name='z')
#目标
model.setObjective(x+y+z,GRB.MAXIMIZE)
#约束
model.addConstr(x+y<=10)
model.addConstr(z==min_(x,y,6))
#求解
model.Params.OutputFlag=0
model.optimize()
#输出
print('obj=',model.objVal)
for v in model.getVars():
print(v.varName+'='+str(v.x))
另外,对于取绝对值的情况,例如y=|x|,可写成model.addConstr(y==abs_(x))。
如上约束也可写成:
min z
x<=My
y={0,1}
z>=cx+k+(y-1)M
对于分段式目标函数,我们也可通过变量的转化将其进行线性化处理:
对于逻辑或的约束,即两个约束至少有一个成立的问题,可按照如下方式进行线性化:
总结
后续将继续跟进Gurobi的学习,做一些复杂问题的建模求解。