Warm-Starting with Initial Solutions

Warm-Starting with Initial Solutions#

When you already have a good feasible solution (from a previous solve, a heuristic, domain knowledge, or a related problem), you can pass it to m.solve() as an initial_solution. This provides two benefits in the Branch and Bound algorithm [Belotti et al., 2013, Land and Doig, 1960]:

  1. The solution is injected as the initial incumbent, giving the solver an upper bound from the start. This allows it to prune suboptimal branches immediately rather than exploring them first.

  2. For continuous (NLP) subproblems, the solution is used as the starting point for the nonlinear solver [Nocedal and Wright, 2006], which can improve convergence.

This is analogous to “MIP-start” functionality in commercial solvers [Koch et al., 2011].

import discopt.modeling as dm

Basic Usage#

The initial_solution parameter accepts a dictionary mapping Variable objects to their values. Scalars, lists, and numpy arrays are all supported.

The cell below prints Nodes: 0, and that is not the warm start doing its work: this model is a convex MIQP, so it is routed to outer approximation and certified at the root whether or not a starting point is given. Measured on this model, cold and warm both open 0 nodes on the default route and 1 node under solver="bb". Benefit (1) above — pruning with an early incumbent — needs a tree to prune, so the cell reports the route alongside the node count rather than letting 0 read as a result.

# A simple MINLP: minimize (x - 2.5)^2 + (y - 1)^2 with y binary
m = dm.Model("warm_start_demo")
x = m.continuous("x", lb=0, ub=5)
y = m.binary("y")

m.minimize((x - 2.5) ** 2 + (y - 1) ** 2)
m.subject_to(x <= 4 * y + 1, name="linking")

# Solve with a known good initial solution
result = m.solve(initial_solution={x: 2.5, y: 1.0})

print(f"Status:    {result.status}")
print(f"Objective: {result.objective:.6f}")
print(f"x = {result.x['x']:.4f}")
print(f"y = {result.x['y']:.4f}")
print(f"Nodes:     {result.node_count}")
print(f"Route:     {result.algorithm_route}")

# The same model without the warm start: identical node count, because the route
# certifies at the root either way. Pin that, so "Nodes: 0" is never read as the
# warm start having pruned anything.
cold = m.solve()
print(f"\nWithout initial_solution: nodes={cold.node_count}, obj={cold.objective:.6f}")
assert cold.node_count == result.node_count
assert abs(cold.objective - result.objective) < 1e-9

OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)

OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)

Status:    optimal
Objective: 0.000000
x = 2.5000
y = 1.0000
Nodes:     0
Route:     mip-nlp/oa: miqp certified convex at the root (DISCOPT_CONVEX_MINLP_ROUTE)

Without initial_solution: nodes=0, obj=0.000000

Warm-Starting a Continuous NLP#

For purely continuous problems, the initial solution is used as the starting point for the interior-point solver. When the starting point is close to the optimum, this can reduce the number of iterations needed.

# Continuous NLP: minimize (x-3)^2 + (y-4)^2
m = dm.Model("nlp_warm")
x = m.continuous("x", lb=-10, ub=10)
y = m.continuous("y", lb=-10, ub=10)
m.minimize((x - 3) ** 2 + (y - 4) ** 2)

# Cold start (default midpoint)
result_cold = m.solve()
print(f"Cold start: obj={result_cold.objective:.6f}, time={result_cold.wall_time:.4f}s")

# Warm start near the optimum
result_warm = m.solve(initial_solution={x: 2.9, y: 3.9})
print(f"Warm start: obj={result_warm.objective:.6f}, time={result_warm.wall_time:.4f}s")

Cold start: obj=0.000000, time=0.0061s
Warm start: obj=0.000000, time=0.0024s

Primal-dual warm start: SolveResult.kkt#

The cells above warm-start the primal point. On the single-NLP routes whose backend reports them, a solve also hands back its terminal KKT residuals, and solve(warm_start=...) uses them to seed the next solve’s interior point — not just its x (#1247):

first  = model.solve()
second = model.solve(warm_start=first)   # reuses first.kkt["barrier_parameter"]

result.kkt is a dict[str, float], or None. It is None on every branch-and-bound route — there is no single NLP there whose residuals these would be — so read it with .get() rather than assuming it is populated:

key

meaning

primal_infeasibility, dual_infeasibility, complementarity, kkt_error

POUNCE’s final_*, measured on its internally scaled problem — what its own convergence test runs on

the same four suffixed _unscaled

the residuals in the model’s own units

barrier_parameter

terminal interior-point mu, forwarded to seed the next solve

 
res = model.solve()
kkt = res.kkt or {}
print(kkt.get("kkt_error_unscaled"), kkt.get("barrier_parameter"))

Which pair to read matters. A certificate stated in problem units must be built from the _unscaled residuals; the scaled ones answer “did the solver converge”, which is a different question and a different number.

This field is diagnostic and never load-bearing: nothing in the solver reads it back, so a backend that omits a key cannot change a verdict.

Partial Solutions#

You do not need to provide values for every variable. Any variables not included in initial_solution are initialized to the midpoint of their bounds.

# Only provide the binary variable; let the solver find x
m = dm.Model("partial")
x = m.continuous("x", lb=0, ub=5)
y = m.binary("y")
m.minimize((x - 2.5) ** 2 + (y - 1) ** 2)
m.subject_to(x <= 4 * y + 1)

result = m.solve(initial_solution={y: 1.0})
print(f"Status: {result.status}, obj={result.objective:.6f}")
print(f"x={result.x['x']:.4f}, y={result.x['y']:.4f}")

OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)

Status: optimal, obj=0.000000
x=2.5000, y=1.0000

Validation and Safety#

The initial_solution is validated before solving. discopt checks variable bounds and integrality, producing warnings and automatically correcting violations rather than failing silently.

import warnings

m = dm.Model("validation_demo")
x = m.continuous("x", lb=0, ub=5)
y = m.binary("y")
m.minimize((x - 2.5) ** 2 + (y - 1) ** 2)
m.subject_to(x <= 4 * y + 1)

# Out-of-bounds value: x=10 exceeds ub=5, gets clamped to 5
# Non-integer binary: y=0.7 gets rounded to 1
with warnings.catch_warnings(record=True) as w:
    warnings.simplefilter("always")
    result = m.solve(initial_solution={x: 10.0, y: 0.7})
    for warning in w:
        print(f"Warning: {warning.message}")

print(f"\nStatus: {result.status}, obj={result.objective:.6f}")

OA: generating OA cuts for 0 of 1 rows (1 convex; the rest are either nonconvex or already exact in the master)

Warning: Variable 'x': 1 value(s) outside bounds. Clamping to [lb, ub].
Warning: Variable 'y': 1 value(s) are not integer-valued. Rounding to nearest integer.

Status: optimal, obj=0.000000

Programmatic Feasibility Checking#

You can also check feasibility of a candidate solution without solving, using discopt.warm_start.check_feasibility.

from discopt.warm_start import check_feasibility, validate_initial_solution

m = dm.Model("feasibility_check")
x = m.continuous("x", lb=0, ub=5)
y = m.binary("y")
m.minimize((x - 2.5) ** 2 + (y - 1) ** 2)
m.subject_to(x <= 4 * y + 1, name="linking")

# Check a feasible solution
x_flat = validate_initial_solution(m, {x: 2.5, y: 1.0})
is_feasible, violations = check_feasibility(m, x_flat)
print(f"Solution {{x=2.5, y=1}}: feasible={is_feasible}")

# Check an infeasible solution (x=5 but y=0 means x <= 1)
x_flat = validate_initial_solution(m, {x: 5.0, y: 0.0})
is_feasible, violations = check_feasibility(m, x_flat)
print(f"Solution {{x=5, y=0}}: feasible={is_feasible}")
for v in violations:
    print(f"  - {v}")

Solution {x=2.5, y=1}: feasible=True
Solution {x=5, y=0}: feasible=False
  - Constraint 'linking': value 4 > 0 (sense <=)