Skip to content

Callback: User Heuristic and Solution Verification

Consider again the Vehicle Routing Problem with Time Windows and let's add a primal heuristic for finding good upper bounds.

In Where.PathMIPHeuristic one get's access to the current liner relaxation and can use it as a starting point for a heuristic. The example heuristic is not too good, we will add the one-customer route solution every time.

Using the the Where.PathMIPSolution the user can verify feasibility of solutions found by Flowty. For this scenario the solutions are never infeasible.

# Vehicle Routing Problem with Time Windows

from flowty import Model, xsum, CallbackModel, Where
from flowty.datasets import vrp_rep

bunch = vrp_rep.fetch_vrp_rep("solomon-1987-r1", instance="R101_025")
name, n, es, c, d, Q, t, a, b, x, y = bunch["instance"]

m = Model()

# one graph, it is identical for all vehicles
g = m.addGraph(obj=c, edges=es, source=0, sink=n - 1, L=1, U=n - 2, type="B")

# adds resources variables to the graph.
# demand and capacity
    graph=g, consumptionType="V", weight=d, boundsType="V", lb=0, ub=Q, name="d"

# travel time and customer time windows
    graph=g, consumptionType="E", weight=t, boundsType="V", lb=a, ub=b, name="t"

# Add a solution consisting of all 1-customer routes
# not a clever heuristic :)
def callback(cb: CallbackModel, where: Where):
    # Heuristic
    if where == Where.PathMIPHeuristic:
        x = cb.x  # lp relaxation

        # add all 1-customer routes
        cost = 0
        xEdges = [0] * len(x)
        for i in range(n)[1:-1]:
            index = es.index((0, i))
            xEdges[index] = 1
            cost += c[index]
            index = es.index((i, n - 1))
            xEdges[index] = 1
            cost += c[index]

        cb.addSolution(cost, xEdges)

    # Verify solution
    if where == Where.PathMIPSolution:
        x = cb.x  # candidate solution

        # check is solution is infeasible and skip it if so
        isInfeasible = False
        if isInfeasible:


# set partitioning constraints
for i in range(n)[1:-1]:
    m.addConstr(xsum(1 * x for x in g.vars if i == x.source) == 1)
    m.addPackingSet([x for x in g.vars if i == x.source])

status = m.optimize()
print(f"ObjectiveValue {round(m.objectiveValue, 1)}")

# get the variable values
for var in m.vars:
    if var.x > 0:
        print(f"{} = {round(var.x, 1)}")