Python调用Gurobi求解简单线性规划问题

前言

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<=4x+y>=1x,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=0elif abs(yrelaxation-1)<=0.01:suby[i].lb=1submodel.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的学习，做一些复杂问题的建模求解。