Callbacks and User-Defined Cuts#

discopt’s Branch and Bound solver supports three callback types that let you interact with the search process: node callbacks for monitoring progress, lazy constraint callbacks for adding cutting planes on the fly, and incumbent callbacks for filtering candidate solutions.

Callbacks are a standard mechanism in mixed-integer solvers [Quesada and Grossmann, 1992] for implementing problem-specific logic without modifying the solver itself. Lazy constraints, in particular, are essential for problems where the full constraint set is exponentially large (such as subtour elimination in the TSP) and constraints are only generated as needed [Sawaya and Grossmann, 2005].

Note

Callbacks never affect the correctness of the solver’s internal math. The solver catches and logs any exceptions raised by user callbacks, so a buggy callback will not crash the solve.

import os

os.environ["JAX_PLATFORMS"] = "cpu"
os.environ["JAX_ENABLE_X64"] = "1"

import discopt
import numpy as np
from discopt.callbacks import CallbackContext, CutResult

Node Callback: Logging B&B Progress#

The simplest callback type is the node callback, which is called after each batch of B&B nodes is processed. It receives a CallbackContext with the current state of the search (node count, incumbent objective, best bound, gap, and elapsed time). This is useful for custom logging, progress bars, or early termination logic.

# Build a small MILP
m = discopt.Model("progress_demo")
x = m.binary("x", shape=(5,))
costs = np.array([3.0, 7.0, 2.0, 5.0, 8.0])
m.minimize(sum(costs[i] * x[i] for i in range(5)))
m.subject_to(x[0] + x[1] + x[2] >= 2, name="cover_1")
m.subject_to(x[2] + x[3] + x[4] >= 2, name="cover_2")

# Define a node callback that logs progress
log_entries = []


def progress_logger(ctx: CallbackContext, model: discopt.Model) -> None:
    entry = {
        "nodes": ctx.node_count,
        "incumbent": ctx.incumbent_obj,
        "bound": ctx.best_bound,
        "gap": ctx.gap,
        "time": f"{ctx.elapsed_time:.3f}s",
    }
    log_entries.append(entry)
    print(
        f"  Nodes: {ctx.node_count:>5}  "
        f"Incumbent: {ctx.incumbent_obj}  "
        f"Bound: {ctx.best_bound:.4f}  "
        f"Gap: {ctx.gap}"
    )


result = m.solve(node_callback=progress_logger, time_limit=30)
print(f"\nFinal: status={result.status}, obj={result.objective}")
Final: status=optimal, obj=10.0

Lazy Constraints: Subtour Elimination for TSP#

Lazy constraints are the most powerful callback type. They allow you to add linear cutting planes during the solve, which is essential for problems with an exponential number of constraints.

The classic example is the Travelling Salesman Problem (TSP), where subtour elimination constraints (SECs) prevent solutions that form disconnected cycles. Rather than adding all \(O(2^n)\) SECs upfront, we add them lazily: whenever the solver finds an integer-feasible solution that contains a subtour, we add a cut that eliminates it [Westerlund and Pettersson, 1995].

The callback receives a CallbackContext and returns a list of CutResult objects. Each CutResult specifies the cut as a list of (variable, coefficient) pairs, a sense ("<=", ">=", or "=="), and a right-hand side value.

# Small 4-city TSP with lazy subtour elimination
# Cities: 0, 1, 2, 3
# Decision variables: x[i,j] = 1 if edge (i,j) is in the tour
n = 4
# Distance matrix (symmetric)
dist = np.array(
    [
        [0, 10, 15, 20],
        [10, 0, 35, 25],
        [15, 35, 0, 30],
        [20, 25, 30, 0],
    ],
    dtype=float,
)

m = discopt.Model("tsp_4city")

# Binary variables for edges (upper triangle only, i < j)
edges = {}
edge_vars = []
for i in range(n):
    for j in range(i + 1, n):
        var = m.binary(f"x_{i}_{j}")
        edges[(i, j)] = var
        edge_vars.append((i, j, var))

# Minimize total distance
m.minimize(sum(dist[i, j] * var for i, j, var in edge_vars))

# Degree constraints: each city has exactly 2 edges
for city in range(n):
    incident = []
    for i, j, var in edge_vars:
        if i == city or j == city:
            incident.append(var)
    m.subject_to(sum(incident) == 2, name=f"degree_{city}")


def find_subtours(sol_flat, edge_list, n_cities):
    """Find connected components (subtours) in the solution."""
    # Build adjacency from solution
    adj = {i: [] for i in range(n_cities)}
    for idx, (i, j, _) in enumerate(edge_list):
        if sol_flat[idx] > 0.5:
            adj[i].append(j)
            adj[j].append(i)

    visited = set()
    components = []
    for start in range(n_cities):
        if start in visited:
            continue
        component = set()
        stack = [start]
        while stack:
            node = stack.pop()
            if node in component:
                continue
            component.add(node)
            for nb in adj[node]:
                if nb not in component:
                    stack.append(nb)
        visited |= component
        components.append(component)
    return components


def subtour_callback(ctx: CallbackContext, model: discopt.Model):
    """Add subtour elimination cuts for any disconnected component."""
    components = find_subtours(ctx.x_relaxation, edge_vars, n)

    if len(components) <= 1:
        return []  # Single tour, no subtours

    cuts = []
    for comp in components:
        if len(comp) == n:
            continue  # Full tour
        # SEC: sum of edges within the subtour <= |S| - 1
        terms = []
        for i, j, var in edge_vars:
            if i in comp and j in comp:
                terms.append((var, 1.0))
        if terms:
            cuts.append(
                CutResult(
                    terms=terms,
                    sense="<=",
                    rhs=float(len(comp) - 1),
                )
            )
    return cuts


result = m.solve(lazy_constraints=subtour_callback, time_limit=60)
print(f"Status: {result.status}")
print(f"Optimal tour cost: {result.objective}")
if result.x is not None:
    print("Edges in tour:")
    for i, j, var in edge_vars:
        val = result.x[f"x_{i}_{j}"]
        if np.all(val > 0.5):
            print(f"  {i} -- {j}  (dist={dist[i, j]:.0f})")
Status: optimal
Optimal tour cost: 80.0
Edges in tour:
  0 -- 1  (dist=10)
  0 -- 2  (dist=15)
  1 -- 3  (dist=25)
  2 -- 3  (dist=30)

Incumbent Callback: Filtering Solutions#

The incumbent callback is called whenever the solver finds a new best integer-feasible solution. You can inspect the solution and return False to reject it. This is useful for enforcing problem-specific feasibility conditions that are difficult to express as algebraic constraints.

For example, you might reject solutions that violate a complex business rule, or solutions that fail an external simulation check.

# Example: reject solutions where fewer than 2 items are selected
m = discopt.Model("filter_demo")
x = m.binary("x", shape=(4,))
profits = np.array([10.0, 6.0, 4.0, 2.0])
m.maximize(sum(profits[i] * x[i] for i in range(4)))
m.subject_to(sum(3 * x[i] for i in range(4)) <= 7, name="capacity")


def require_min_items(ctx, model, solution):
    """Reject solutions with fewer than 2 selected items."""
    x_val = solution["x"]
    n_selected = int(np.sum(x_val > 0.5))
    if n_selected < 2:
        print(f"  Rejected: only {n_selected} item(s) selected")
        return False
    print(f"  Accepted: {n_selected} items selected, obj={ctx.incumbent_obj}")
    return True


result = m.solve(incumbent_callback=require_min_items, time_limit=30)
print(f"\nFinal: status={result.status}, obj={result.objective}")
if result.x is not None:
    print(f"Selected items: {np.where(result.x['x'] > 0.5)[0]}")
Final: status=optimal, obj=16.0
Selected items: [0 1]

API Summary#

The callback API consists of three components in discopt.callbacks:

CallbackContext is a dataclass passed to all callbacks containing the current B&B state: node_count, incumbent_obj (or None), best_bound, gap (or None), elapsed_time, x_relaxation (the current node’s solution as a flat numpy array), and node_bound.

CutResult specifies a linear cut via terms (a list of (Variable, coefficient) tuples), sense ("<=", ">=", or "=="), and rhs (the right-hand side value). Cuts are converted to dense constraint vectors internally and added to the solver’s cut pool.

The three callback types are passed as keyword arguments to Model.solve():

  • node_callback(ctx, model) -> None: called after each batch of nodes.

  • lazy_constraints(ctx, model) -> list[CutResult]: called at integer-feasible nodes, returns cuts to add.

  • incumbent_callback(ctx, model, solution) -> bool: called when a new incumbent is found, returns False to reject.

All callbacks are optional and can be combined freely [Smith and Pantelides, 1999].