Gurobi中的Callback函数用法：如何干预优化求解过程

Word count: 1.8k   |   Reading time≈ 7 min

本文主要分为三个部分，第一部分介绍callback函数的用途以及调用callback函数方法；第二部分参考gurobi官方的用户手册，介绍callback函数具体能干预的函数方法；第三部分解释gurobi官方的tsp.py中使用callback函数添加subtour约束的方法。

什么是Callback & 什么时候可以调用Callback

在optimize的函数介绍下有callback的描述：

通常使用Gurobi建完模型后重要一步是optimize，其可选参数就是callback函数，用于在指定条件触发情况下对优化求解的过程进行一定干预。

Callback函数主要由两个输入，分别是model和where，在求解过程中会周期性或者由条件触发该函数。其中model是需要被优化的模型，where就是触发callback函数的阶段，在每一个where下由很多的属性可以调用，称为what，如下图所示

比如表中的presolve阶段有PRE_COLDEL和PRE_ROWDEL属性，可以对其进行判断并进行对应的求解干预，比如终止优化、添加行添加列等等。下面的代码展示了在presolve阶段调用callback

1
2
3
4
5
6
7
8
9
def mycallback(model, where):
    if 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))

model.optimize(mycallback)

这段代码只执行了print操作，并没有实际干预求解过程，所以callback都能干预求解的哪些参数呢。

Callback怎样干预求解过程

参考gurobi10.0的用户手册，发现有如下方法可以调用

Model.cbGet()：查询信息

1
2
3
4
def mycallback(model, where):
    if where == GRB.Callback.SIMPLEX:
        print(model.cbGet(GRB.Callback.SPX_OBJVAL))
model.optimize(mycallback)

查询单纯形法里面的目标函数值（Current simplex objective value）。

Model.cbGetNodeRel()：查询当前节点的松弛解

本方法只应用于where是MIPNODE并且求得最优解的情况下。

1
2
3
4
5
6
7
def mycallback(model, where):
    if where == GRB.Callback.MIPNODE:
        status = model.cbGet(GRB.Callback.MIPNODE_STATUS)
        if status == GRB.OPTIMAL:
            print(model.cbGetNodeRel(model._vars))
model._vars = model.getVars()
model.optimize(mycallback)

该callback含义：在当前节点取得最优值的情况下，print所有变量的松弛解

Model.cbCut()：添加割平面

1
2
3
4
5
6
7
8
9
10
def mycallback(model, where):
    if where == GRB.Callback.MIPNODE:
        status = model.cbGet(GRB.Callback.MIPNODE_STATUS)
        if status == GRB.OPTIMAL:
        rel = model.cbGetNodeRel([model._vars[0], model._vars[1]])
        if rel[0] + rel[1] > 1.1:
            model.cbCut(model._vars[0] + model._vars[1] <= 1)

model._vars = model.getVars()
model.optimize(mycallback)

如果当前节点取到最优解，并且如果第前两个变量的松弛解之和大于1.1，则添加前两个变量之和小于等于1的cut。

Model.cbGetSolution()：查询可行解的值

1
2
3
4
5
def mycallback(model, where):
    if where == GRB.Callback.MIPSOL:
        print (model.cbGetSolution(model._vars))
model._vars = model.getVars()
model.optimize(mycallback)

print所有可行解的取值

Model.cbLazy()：添加lazycut

Lazy constraint是在求解完毕时候进行检验，如果违反了再把约束加到模型中，因为Lazy constraint可能影响求解速度。

如果通过回调函数创建Lazy constraint，则必须设置参数LazyConstraint = 1。否则，Gurobi可能会应用对Lazy constraint无效的dual presolve reductions。

如果在不添加任何Lazy constraint的情况下设置LazyConstraint = 1参数，则Gurobi的算法行为将发生变化（因为禁用了dual presolve reductions），但最优解的值不会变。

1
2
3
4
5
6
7
8
def mycallback(model, where):
    if where == GRB.Callback.MIPSOL:
        sol = model.cbGetSolution([model._vars[0], model._vars[1]]) #model.cbGetSolution()获得新MIP solution函数值；()内的var是指定变量，表示想知道这些变量在solution中的值；返回括号内指定变量在solution中的值
            if sol[0] + sol[1] > 1.1:
                model.cbLazy(model._vars[0] + model._vars[1] <= 1)
model._vars = model.getVars()
model.Params.lazyConstraints = 1
model.optimize(mycallback)

如果求解的前两个变量的和大于1.1，则增加lazycut让前两个变量求和小于等于1

Model.cbProceed():

1
2
3
4
5
6
7
def mycallback(model, where):
    if where == GRB.Callback.MIPSOL:
        phase = model.cbGet(GRB.Callback.MIPSOL_PHASE)
        obj = model.cbGet(GRB.Callback.MIPSOL_OBJ)
        if phase == GRB.PHASE_MIP_NOREL and obj < target_obj:
            model.cbProceed()
model.optimize(mycallback)

用于跳转到下一个求解阶段。如果目前是处于MIP搜索阶段，并且当前目标小于target_obj，则跳过当前阶段（理解为continue？）

Model.cbSetSolution()，Model.cbUseSolution()：向当前节点导入解，并快速求解

1
2
3
4
5
def mycallback(model, where):
    if where == GRB.Callback.MIPNODE:
        model.cbSetSolution(vars, newsolution)
        objval = model.cbUseSolution()
model.optimize(mycallback)

把vars变量赋值为newsolution，通过model.cbUseSolution()快速求解一个obj

Model.cbStopOneMultiObj()：

用于中断多目标求解。
Interrupt the optimization process of one of the optimization steps in a multi-objective MIP problem without
stopping the hierarchical optimization process. Only available for multi-objective MIP models and when the where
member variable is not equal to GRB.Callback.MULTIOBJ

1
2
3
4
5
6
7
8
9
10
11
12
13
def mycallback(model, where):
    if where == GRB.Callback.MULTIOBJ:
        # get current objective number
        model._objcnt = model.cbGet(GRB.Callback.MULTIOBJ_OBJCNT)
        # reset start time to current time
        model._starttime = time.time()
    # See if we want to stop current multiobjective step
    else if time.time() - model._starttime > 1000: # or solution is good enough
        # stop only this optimization step
        model.cbStopOneMultiObj(model._objcnt)
model._objcnt = 0
model._starttime = time.time()
model.optimize(mycallback)

如果求解多目标的时间大于1000，则不在进行多目标求解（把求解中止的目标数目设置为当前已经求解的数目）

使用callback为TSP添加subtour约束

参考gurobi官方tsp案例，其中使用callback添加subtour约束的函数如下（n是节点数目）：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
# Callback - use lazy constraints to eliminate sub-tours
def subtourelim(model, where):
    if where == GRB.Callback.MIPSOL:
        vals = model.cbGetSolution(model._vars)
        # find the shortest cycle in the selected edge list
        tour = subtour(vals)
        if len(tour) < n:
            # add subtour elimination constr. for every pair of cities in tour
            model.cbLazy(gp.quicksum(model._vars[i, j]
                                     for i, j in combinations(tour, 2))
                         <= len(tour)-1)


# Given a tuplelist of edges, find the shortest subtour

def subtour(vals):
    # make a list of edges selected in the solution
    edges = gp.tuplelist((i, j) for i, j in vals.keys()
                         if vals[i, j] > 0.5)
    unvisited = list(range(n))
    cycle = range(n+1)  # initial length has 1 more city
    while unvisited:  # true if list is non-empty
        thiscycle = []
        neighbors = unvisited
        while neighbors:
            current = neighbors[0]
            thiscycle.append(current)
            unvisited.remove(current)
            neighbors = [j for i, j in edges.select(current, '*')
                         if j in unvisited]
        if len(cycle) > len(thiscycle):
            cycle = thiscycle
    return cycle

如果求解生成的回路没有包含全部节点，判定他为一个subtour，并加约束让当前子回路中的几个边不同时出现。

总结

gurobi参考手册里面的第16小节就是 Callback Codes，专门描述了callback可以适用的阶段（where）和可以提取的参数（what）。使用的阶段从presolve到simplex、MIP再到barrier，可以调用的参数也是数不胜数，根据实际的问题需求合理的加速求解，或者在合适时机和合适条件下加入cut来干预优化求解的过程，可以起到不错的效果。不过这个cut是在求解过程中去加的，之前的decomposition method大多数是分为了主问题和一个或多个子问题，并不适用这种callback添加cut的方法。之后如果做到需要加速单次求解，或者在求解过程中添加cut，那callback则必不可少。

（本来很早就计划把gurobi的reference manual学习一遍了，一直拖延症没看，没有ddl的前提下提高执行力确实很困难。

人长期往上走的首要条件是，稳住下限，让低质量的决策尽量减少发生，然后再慢慢变好。不水幻想一夜暴富或力挽狂澜。长久布局和持之以恒，才是无法被打败的铜墙铁壁。

Donate
  • Copyright： Copyright is owned by the author. For commercial reprints, please contact the author for authorization. For non-commercial reprints, please indicate the source.
  • Copyrights © 2015-2024 galaxy
  • Visitors: | Views:

请我喝杯咖啡吧~

支付宝
微信