# Benders Decomposition (Lazy Constraints)

Originally Contributed by: Mathieu Besançon

This notebook describes how to implement the Benders decomposition in JuMP using lazy constraints. We keep the same notation and problem form as the first notebook Benders decomposition.

\begin{aligned} & \text{maximize} \quad &&c_1^T x+c_2^T v \\ & \text{subject to} \quad &&A_1 x+ A_2 v \preceq b \\ & &&x \succeq 0, x \in \mathbb{Z}^n \\ & &&v \succeq 0, v \in \mathbb{R}^p \\ \end{aligned}

where $b \in \mathbb{R}^m$, $A_1 \in \mathbb{R}^{m \times n}$, $A_2 \in \mathbb{R}^{m \times p}$ and $\mathbb{Z}$ is the set of integers.

Here the symbol $\succeq$ ($\preceq$) stands for element-wise greater (less) than or equal to. Any mixed integer programming problem can be written in the form above.

For a detailed explanation on the Benders decomposition algorithm, see the introduction notebook.

### Lazy constraints

Some optimization solvers allow users to interact with them during the solution process by providing user-defined functions which are triggered under certain conditions.

The generic term for these functions is callback.

In integer optimization, the main callback types are lazy constraints, user-cuts and heuristic solutions. See the Callbacks section for an introduction on how to use them.

Some callbacks define a new constraint which is only activated when necessary, i.e., when a current solution does not respect them. It can avoid building an optimization model with too many constraints up-front.

This is the case for Benders decomposition, since the sub-problem defines an exponential number of primal vertices and therefore dual cuts.

A detailed explanation on the distinction between user-cuts and lazy constraints is also available on Paul Rubin's blog. He also describes this approach to Benders Decomposition.

We use the data from the original notebook and change the solution algorithm to leverage lazy constraints:

• Step 1 (Initialization)
• Step 2 (defining the subproblem model)
• Step 3 (registering the lazy constraint of the subproblem)

### Data for the problem

The input data is from page 139, Integer programming by Garfinkel and Nemhauser.

c1 = [-1; -4]
c2 = [-2; -3]

dim_x = length(c1)
dim_u = length(c2)

b = [-2; -3]

A1 = [1 -3;
-1 -3]
A2 = [1 -2;
-1 -1]

M = 1000;

using JuMP
import GLPK
import Test

Subproblem creation

function build_subproblem()
sub_problem_model = Model(GLPK.Optimizer)
@variable(sub_problem_model, u[1:dim_u] >= 0)
@constraint(
sub_problem_model,
constr_ref_subproblem[j = 1:size(A2, 2)],
A2[:, j]' * u >= c2[j],
)
return (sub_problem_model, u)
end
build_subproblem (generic function with 1 method)

### Master Problem Description

master_problem_model = Model(GLPK.Optimizer);

(sub_problem_model, u) = build_subproblem();

Variable Definition

@variable(master_problem_model, 0 <= x[1:dim_x] <= 1e6, Int)
@variable(master_problem_model, t <= 1e6)
$$$t$$$

Objective Setting

@objective(master_problem_model, Max, t)

print(master_problem_model)
Max t
Subject to
x ≥ 0.0
x ≥ 0.0
x ≤ 1.0e6
x ≤ 1.0e6
t ≤ 1.0e6
x integer
x integer

Track the calls to the callback

iter_num = 0
0

Define lazy constraints

function benders_lazy_constraint_callback(cb_data)
global iter_num
iter_num += 1
println("Iteration number = ", iter_num)

x_current = callback_value.(Ref(cb_data), x)
fm_current = callback_value(cb_data, t)

c_sub = b - A1 * x_current
@objective(sub_problem_model, Min, c1' * x_current + c_sub' * u)
optimize!(sub_problem_model)

print("\nThe current subproblem model is \n")
print(sub_problem_model)

t_status_sub = termination_status(sub_problem_model)
p_status_sub = primal_status(sub_problem_model)

fs_x_current = objective_value(sub_problem_model)

u_current = value.(u)

γ = b' * u_current

if p_status_sub == MOI.FEASIBLE_POINT && fs_x_current  ≈  fm_current # we are done
@info("No additional constraint from the subproblem")
end

if p_status_sub == MOI.FEASIBLE_POINT && fs_x_current < fm_current
println("\nThere is a suboptimal vertex, add the corresponding constraint")
cv = A1' * u_current - c1
new_optimality_cons = @build_constraint(t + cv' * x <= γ)
MOI.submit(
master_problem_model,
MOI.LazyConstraint(cb_data),
new_optimality_cons,
)
end

if t_status_sub == MOI.INFEASIBLE_OR_UNBOUNDED
println("\nThere is an  extreme ray, adding the corresponding constraint")
ce = A1' * u_current
new_feasibility_cons = @build_constraint(dot(ce, x) <= γ)
MOI.submit(
master_problem_model,
MOI.LazyConstraint(cb_data),
new_feasibility_cons,
)
end
end

MOI.set(
master_problem_model,
MOI.LazyConstraintCallback(),
benders_lazy_constraint_callback,
)

optimize!(master_problem_model)

t_status = termination_status(master_problem_model)
p_status = primal_status(master_problem_model)

if p_status == MOI.INFEASIBLE_POINT
println("The problem is infeasible :-(")
end

if t_status == MOI.INFEASIBLE_OR_UNBOUNDED
fm_current = M
x_current = M * ones(dim_x)
end

if p_status == MOI.FEASIBLE_POINT
fm_current = value(t)
x_current = value.(x)
end

println("Status of the master problem is ", t_status,
"\nwith fm_current = ", fm_current,
"\nx_current = ", x_current)
Iteration number = 1

The current subproblem model is
Min -2 u - 3 u
Subject to
constr_ref_subproblem : u - u ≥ -2.0
constr_ref_subproblem : -2 u - u ≥ -3.0
u ≥ 0.0
u ≥ 0.0

There is a suboptimal vertex, add the corresponding constraint
Iteration number = 2

The current subproblem model is
Min 750003.75 u + 750002.75 u - 1.0000076666666666e6
Subject to
constr_ref_subproblem : u - u ≥ -2.0
constr_ref_subproblem : -2 u - u ≥ -3.0
u ≥ 0.0
u ≥ 0.0

There is a suboptimal vertex, add the corresponding constraint
Iteration number = 3

The current subproblem model is
Min 0.875 u - 0.125 u - 3.8333333333333335
Subject to
constr_ref_subproblem : u - u ≥ -2.0
constr_ref_subproblem : -2 u - u ≥ -3.0
u ≥ 0.0
u ≥ 0.0

There is a suboptimal vertex, add the corresponding constraint
Iteration number = 4

The current subproblem model is
Min -3.8333333333333335
Subject to
constr_ref_subproblem : u - u ≥ -2.0
constr_ref_subproblem : -2 u - u ≥ -3.0
u ≥ 0.0
u ≥ 0.0
[ Info: No additional constraint from the subproblem
Iteration number = 5

The current subproblem model is
Min -0.875 u + 0.125 u - 3.8333333333333335
Subject to
constr_ref_subproblem : u - u ≥ -2.0
constr_ref_subproblem : -2 u - u ≥ -3.0
u ≥ 0.0
u ≥ 0.0

There is a suboptimal vertex, add the corresponding constraint
Iteration number = 6

The current subproblem model is
Min u - 5
Subject to
constr_ref_subproblem : u - u ≥ -2.0
constr_ref_subproblem : -2 u - u ≥ -3.0
u ≥ 0.0
u ≥ 0.0
[ Info: No additional constraint from the subproblem
Iteration number = 7

The current subproblem model is
Min u - 4
Subject to
constr_ref_subproblem : u - u ≥ -2.0
constr_ref_subproblem : -2 u - u ≥ -3.0
u ≥ 0.0
u ≥ 0.0
[ Info: No additional constraint from the subproblem
Status of the master problem is OPTIMAL
with fm_current = -4.0
x_current = [0.0, 1.0]
1. Garfinkel, R. & Nemhauser, G. L. Integer programming. (Wiley, 1972).