目录
(3)model.cbLazy(lhs, sense, rhs):
(4)model.cbCut(lhs, sense, rhs):
(5)model.cbSetsolution(vars, solution):
在初始使用的时候,很多人觉得,最起码我是这么觉得,gurobi是个黑箱。即:只要输入决策变量,目标函数,和约束,直接调用model.optimize()即可,至于optimize是怎么求解,不了解不知道也没办法干预。
但其实不是的。
model.optimize()的完整形式是model.optimize(callback=none),既我们没有赋予括号中任何信息。但其实callback功能十分强大,在后期实现分支定价,分支切割算法的时候经常用。
callback可以实现监视,干预,也就是部分程度和gurobi求解器实现交流交互,管理gurobi的优化进程。更多callback的分析可以看gurobi自带的参考资料,即refman.pdf文件。
1.callback的使用
首先需要定义callback即回调函数,然后model.optimize(所定义的callback函数),即传入定义的callback函数即可。
而定义callback函数,有固定的模式,是半定制的形式(不像我们编程其它的函数,从头到尾都是自己定义;也不是全部定义好了,只需调用)。
callback里头的参数是where和what两个,where的取值直接决定了what的取值,即在什么地方(where)获取什么信息(what),因此二者要正确对应。where是callback函数的触发点,表明了我们在哪里进行callback的调用;而what是要进行的操作,即能够获取什么信息。其中where的取值有9种,what的取值也有很多种。比较常用的what取值是:
MIP_SOL:当前解的具体取值
MIP_OBJ:新解的目标值
MIP_OBJBST:当前最优的目标值
MIP_OBJBND:当前最优界
MIP_ODCNT:当前已搜索的节点数
MIP_SOLCNT:当前发现可行解的数量
MIP_CUTCNT:当前割平面使用次数
MIP_NODLFT:当前未搜索的节点数
MIP_ITRCNT:当前单纯形迭代步数
...
调用也很简单,如下述:
if where== grb.GRB.Callback.MULTIOBJ: # where
print(model.cbGet(grb.GRB.Callback.MULTIOBJ_OBJCNT)) # what
可以看到,需要用model.cbGet(what)获取信息。
2.callback其它常用函数
除了model.cbGet(what)之外,还有其它的一些常用函数来获取运行过程中的信息或修改模型运行状态。如:
(1)model.cbGetNodeRel(vars):
在callback函数中获取当前节点的线性松弛的解中决策变量的取值,即当前节点的松弛解。vars 为要查询的变量,需要注意的是,只有在where == GRB.Callback.MIPNODE并且 CRB.Callback.MIPNODE STATUS == GRB.OPTIMAL两个条件同时成立时才能使用。
# 查询变量在当前节点的松弛解
def mycallback(model, where):
if where == GRB.Callback.MIPNODE and model.cbGet(GRB.Callback.MIPNODE_STATUS) ==
GRBOPTIMAL:
print(model.cbGetNodeRel(model._vars))
model._vars= model.getVars()
model.optimize(mycallback)
(2)model.cbGetSolution(vars):
获取MIP的当前解(整数解),或者说查询变量在新可行解中的值。
(3)model.cbLazy(lhs, sense, rhs):
向MIP模型中添加一个lazy惰性约束,需要设置model.Params.lazyConstrains=1,即启用lazy约束。使用场景是:当模型中的约束集太大无法全部添加时(比如TSP破除子回路的约束,初始肯定不可能把所有的节点的真子集子回路避免约束全部添加进去),通常会使用lazy约束。通过只包含在分支割平面搜索过程中不满足条件的约束(即只有在被违反时才会起作用),有时也可以在只添加完整约束集的一部分时找到经验证的最优解。在添加lazy cut时应该先查询当前节点的解(通过cbGstSolution获取GRB.CB_MIPSOL或通过该cbGetNodeNodelRel获取GRB.CB_MIP_NODE)。
for m in model.getVars():
if (m.varName.startswith('x')):
a = (int)(m.varName.split('_')[1]) #m.varName为'x_10_9'.split('_')为['x','10','9']
b = (int)(m.varName.split('_')[2])
x_value[a][b] = model.cbGetSolution(m) # 获取model中变量的值
(4)model.cbCut(lhs, sense, rhs):
用于求解MIP问题时,在节点添加割平面。虽然割平面可以添加到树枝和切割树的任何节点,但是它们会增加在每个节点处求解松弛模型的大小,并且会显著降低节点处理的速度,因此应谨慎添加。除此以外,割平面通常用于切断当前松弛解,要在当前节点上检索松弛方案,因此应该先调用cbGetNodeRel。也必须将参数precrush设置为1,即model.Params.PreCrush = 1,意思是关闭gurobi预处理对模型约束的转化。
(5)model.cbSetsolution(vars, solution):
将一个已知解(完整的),或已知解的部分信息(部分解)传递给当前模型。如将启发式求解得到的解作为上界传给gurobi,gurobi可以在它的基础上求解,从而加快模型整体求解。如果要指定多组变量的值,可以多次调用该方法。
# Example usage:
def mycallback(model, where):
if where == GRB.Callback.MIPNODE:
model.cbSetSolution(vars, newsolution)
model.optimize(mycallback)
# 如:
if abs(subX[i + j -1].x) < 0.001: # 松弛解中该变量的取值在0附近
model.cbSetSolution(model._X[(i,j)], 0.0) # 将该变量的取值上界设为0
elif abs(subX[i + j-1].x - 1) < 0.001: # 在1附近
model.cbSetSolution(model._X[(i,j)], 1.0)
(6)model.cbUseSolution():
使用cbUseSolution导入一个解后,可以调用cbUseSolution来计算这个解对应的目标函数。
其实以上内容,可以参照gurobi的官方手册,以下是具体翻译内容:
callback:回调参数
Gurobi回调类。回调是一个用户函数,由gurobi优化器定期调用,以便允许用户查询或修改优化的状态。更准确地说,如果您为Model.optimize或Model.computeIIS传递一个(callback)函数,该函数接受两个参数(model和where)作为model的参数,您的函数将在优化期间被调用。然后你的回调函数可以调用Model.cbGet来查询优化器以获取优化状态的详细信息。
Gurobi回调既可以用来监视优化的进度,也可以用来修改Gurobi优化器的行为。一个简单的用户回调函数可能调用Model.cbGet生成自定义显示,或者提前终止优化(使用Model.terminate),或者继续进行计算的下一阶段(使用Model.cbProceed)。更复杂的MIP回调可能使用Model.cbGetNodeRel或Model.cbGetSolution从解决方案检索值到当前节点,然后使用Model.cbCut或Model.cbLazy很容易添加约束来切断该解决方案,或者适用Model.cbSetSolution导入从该解中构建的启发式解决方案。对于多目标问题,您可以使用Model.cbStopOneMultiObj在不停止分层优化过程的情况下,中断多目标MIP问题中某个优化步骤的优化过程。
Gurobi回调类提供了一组在用户回调函数中使用的常量。该类中的第一组常量列出了用户回调函数的where参数的选项。where参数指示从优化过程中的哪个位置调用用户回调。本文档的回调代码部分列出了选项。
该类中的另一组常量列出了Model.cbGet的what参数的选项。用户回调使用what参数来指示它想要检索的状态信息。完整的选项列表可以在回调代码部分找到。与where参数一样,您通过GRB.Callback引用一个what常量。例如,将使用GRB.Callback.SPX_OBJVAL请求单纯形目标值。
如果您希望将数据传递给回调函数,您可以通过Model对象这样做。例如,如果在优化之前您的程序包含语句model._value = 1,然后您的回调函数可以查询model._value的值。注意,用户数据字段的名称必须以下划线开头。
当使用多个线程求解模型时,用户回调只会从单个线程调用,因此您不需要担心回调的线程安全性。
注意,不支持在回调中更改参数,这样做可能会导致未定义的行为。
您可以查看示例目录中的callback.py,了解如何使用Gurobi回调的详细信息。
官方的callback.py文件:
# Copyright 2022, Gurobi Optimization, LLC
# This example reads a model from a file, sets up a callback that
# monitors optimization progress and implements a custom
# termination strategy, and outputs progress information to the
# screen and to a log file.
#
# The termination strategy implemented in this callback stops the
# optimization of a MIP model once at least one of the following two
# conditions have been satisfied:
# 1) The optimality gap is less than 10%
# 2) At least 10000 nodes have been explored, and an integer feasible
# solution has been found.
# Note that termination is normally handled through Gurobi parameters
# (MIPGap, NodeLimit, etc.). You should only use a callback for
# termination if the available parameters don't capture your desired
# termination criterion.
import sys
import gurobipy as gp
from gurobipy import GRB
# Define my callback function
def mycallback(model, where):
if where == GRB.Callback.POLLING:
# Ignore polling callback
pass
elif where == GRB.Callback.PRESOLVE:
# Presolve callback
cdels = model.cbGet(GRB.Callback.PRE_COLDEL)
rdels = model.cbGet(GRB.Callback.PRE_ROWDEL)
if cdels or rdels:
print('%d columns and %d rows are removed' % (cdels, rdels))
elif where == GRB.Callback.SIMPLEX:
# Simplex callback
itcnt = model.cbGet(GRB.Callback.SPX_ITRCNT)
if itcnt - model._lastiter >= 100:
model._lastiter = itcnt
obj = model.cbGet(GRB.Callback.SPX_OBJVAL)
ispert = model.cbGet(GRB.Callback.SPX_ISPERT)
pinf = model.cbGet(GRB.Callback.SPX_PRIMINF)
dinf = model.cbGet(GRB.Callback.SPX_DUALINF)
if ispert == 0:
ch = ' '
elif ispert == 1:
ch = 'S'
else:
ch = 'P'
print('%d %g%s %g %g' % (int(itcnt), obj, ch, pinf, dinf))
elif where == GRB.Callback.MIP:
# General MIP callback
nodecnt = model.cbGet(GRB.Callback.MIP_NODCNT)
objbst = model.cbGet(GRB.Callback.MIP_OBJBST)
objbnd = model.cbGet(GRB.Callback.MIP_OBJBND)
solcnt = model.cbGet(GRB.Callback.MIP_SOLCNT)
if nodecnt - model._lastnode >= 100:
model._lastnode = nodecnt
actnodes = model.cbGet(GRB.Callback.MIP_NODLFT)
itcnt = model.cbGet(GRB.Callback.MIP_ITRCNT)
cutcnt = model.cbGet(GRB.Callback.MIP_CUTCNT)
print('%d %d %d %g %g %d %d' % (nodecnt, actnodes,
itcnt, objbst, objbnd, solcnt, cutcnt))
if abs(objbst - objbnd) < 0.1 * (1.0 + abs(objbst)):
print('Stop early - 10% gap achieved')
model.terminate()
if nodecnt >= 10000 and solcnt:
print('Stop early - 10000 nodes explored')
model.terminate()
elif where == GRB.Callback.MIPSOL:
# MIP solution callback
nodecnt = model.cbGet(GRB.Callback.MIPSOL_NODCNT)
obj = model.cbGet(GRB.Callback.MIPSOL_OBJ)
solcnt = model.cbGet(GRB.Callback.MIPSOL_SOLCNT)
x = model.cbGetSolution(model._vars)
print('**** New solution at node %d, obj %g, sol %d, '
'x[0] = %g ****' % (nodecnt, obj, solcnt, x[0]))
elif where == GRB.Callback.MIPNODE:
# MIP node callback
print('**** New node ****')
if model.cbGet(GRB.Callback.MIPNODE_STATUS) == GRB.OPTIMAL:
x = model.cbGetNodeRel(model._vars)
model.cbSetSolution(model.getVars(), x)
elif where == GRB.Callback.BARRIER:
# Barrier callback
itcnt = model.cbGet(GRB.Callback.BARRIER_ITRCNT)
primobj = model.cbGet(GRB.Callback.BARRIER_PRIMOBJ)
dualobj = model.cbGet(GRB.Callback.BARRIER_DUALOBJ)
priminf = model.cbGet(GRB.Callback.BARRIER_PRIMINF)
dualinf = model.cbGet(GRB.Callback.BARRIER_DUALINF)
cmpl = model.cbGet(GRB.Callback.BARRIER_COMPL)
print('%d %g %g %g %g %g' % (itcnt, primobj, dualobj,
priminf, dualinf, cmpl))
elif where == GRB.Callback.MESSAGE:
# Message callback
msg = model.cbGet(GRB.Callback.MSG_STRING)
model._logfile.write(msg)
if len(sys.argv) < 2:
print('Usage: callback.py filename')
sys.exit(0)
# Turn off display and heuristics
gp.setParam('OutputFlag', 0)
gp.setParam('Heuristics', 0)
# Read model from file
model = gp.read(sys.argv[1])
# Open log file
logfile = open('cb.log', 'w')
# Pass data into my callback function
model._lastiter = -GRB.INFINITY
model._lastnode = -GRB.INFINITY
model._logfile = logfile
model._vars = model.getVars()
# Solve model and capture solution information
model.optimize(mycallback)
print('')
print('Optimization complete')
if model.SolCount == 0:
print('No solution found, optimization status = %d' % model.Status)
else:
print('Solution found, objective = %g' % model.ObjVal)
for v in model.getVars():
if v.X != 0.0:
print('%s %g' % (v.VarName, v.X))
# Close log file
logfile.close()