前言

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的学习，做一些复杂问题的建模求解。

# python # 开发语言 # 后端
Logo
华为开发者空间

华为开发者空间，是为全球开发者打造的专属开发空间，汇聚了华为优质开发资源及工具，致力于让每一位开发者拥有一台云主机，基于华为根生态开发、创新。

更多推荐

cover

开源for Huawei，Beam适配GaussDB实践案例分享

avatar 华为开发者空间
cover

华为云数据库亮相创原会：拥抱AI，共创未来

avatar 华为开发者空间
cover

GaussDB典型SQL调优点之自诊断和语句下推调优

avatar 华为开发者空间