army planners with coordination issues on a scale that had not often been imagined. Plane, personnel, gas, gear, coaching, and provides needed to be allotted throughout an immense system, usually below extreme time strain and uncertainty. This setting accelerated the event of operations analysis, a self-discipline based on the concept that advanced choices could possibly be studied systematically by knowledge, scientific reasoning, and arithmetic slightly than addressed solely by expertise and instinct [1].
George Dantzig was among the many researchers immersed on this new world of large-scale planning. Throughout the struggle, he labored for the US Military Air Forces, the place he contributed to the evaluation and planning of army operations. After the struggle, the Air Power challenged him to assist mechanize planning processes that had beforehand been carried out largely by hand. In 1947, Dantzig formulated the overall linear programming downside and developed the simplex methodology, offering a scientific process for locating the perfect allocation of restricted assets amongst competing actions [2, 3]. This analysis continued by Venture SCOOP, the Scientific Computation of Optimum Packages, and later on the RAND Company, which turned one of many principal facilities for the early growth of linear programming and mathematical optimization [2].
This mathematical revolution ultimately reached the Shell Laboratories in Amsterdam, the place Jacques Benders started working in 1955. Collectively along with his workplace mate, Guus Zoutendijk, he encountered troublesome refinery-planning issues involving manufacturing choices, materials flows, processing assets, and prices. After visiting Shell’s Californian department and studying in regards to the potential of linear programming, they returned to the Netherlands and started implementing the simplex methodology for manufacturing planning within the refinery and petrochemical industries. Mathematical optimization promised to assist coordinate industrial techniques whose many interconnected choices might simply overwhelm unaided human reasoning [4].
Nevertheless, the promise of linear programming got here with a irritating limitation. It labored fantastically when choices might fluctuate repeatedly, corresponding to how a lot materials to course of, produce, or transport. Many necessary industrial choices weren’t steady. They required selecting an possibility or rejecting it, activating another or leaving it unused, and answering sure or no. Introducing these discrete variables reworked an in any other case manageable linear program right into a a lot more durable mixed-variable optimization downside. Benders targeted exactly on this complication. In 1959, he and two Shell colleagues introduced an algorithm for linear packages containing binary variables at a RAND symposium, reporting computational experiments carried out on a Ferranti Mark I laptop whose addition time was roughly 1.2 milliseconds [4].
Benders might have searched just for a sooner method to resolve the whole mannequin, however as a substitute he questioned whether or not the whole mannequin wanted to be solved suddenly. The discrete variables created a lot of the problem, but as soon as their values have been fastened, the remaining choices usually shaped an peculiar linear program. This urged a unique technique: isolate the complicating choices, permit one downside to decide on them, and let a second downside decide their operational penalties. As an alternative of treating the optimization mannequin as an indivisible object, Benders divided it into elements that would talk with and progressively enhance each other [5, 6].
This perception produced one thing extra refined than merely fixing two smaller issues. The primary downside, which we now name the grasp downside, proposes a high-level determination whereas initially possessing solely an incomplete understanding of its true penalties. The subproblem evaluates that proposal and returns mathematical data exhibiting the place the grasp has been too optimistic. This data is transformed into a brand new constraint, referred to as a Benders minimize, and added to the grasp earlier than it tries once more. With every iteration, the grasp develops a extra correct illustration of the operational actuality hidden behind its choices.
The concept turned the muse of Benders’s doctoral thesis, Partitioning in Mathematical Programming, accomplished at Utrecht College in 1960 below the supervision of Hans Freudenthal [4, 5]. In 1962, he printed the tactic within the now-classic paper “Partitioning Procedures for Fixing Blended-Variables Programming Issues” [6]. Greater than six a long time later, Benders decomposition stays probably the most influential decomposition strategies in mathematical optimization. It has been prolonged, strengthened, and efficiently utilized to facility location, transportation, supply-chain design, vitality techniques, stochastic programming, scheduling, and lots of different troublesome optimization issues [7].
Benders decomposition has additionally turn into a staple in my very own operations analysis toolkit, however I need to admit that it took me significantly longer than it ought to have to grasp it correctly. Its central perception is elegant and intuitive, but many studying supplies introduce it by dense notation, excessive factors, twin rays, and prolonged derivations earlier than clearly explaining what the algorithm is making an attempt to perform. In consequence, a method constructed round a remarkably easy dialog between two optimization issues can initially seem much more mysterious than it truly is.
This sequence of three articles is my try to supply the Benders decomposition tutorial that I might have cherished to have throughout the first years of my PhD. My goal is to make the tactic accessible to anybody fascinated by operations analysis and mathematical optimization, with out sacrificing the mathematical reasoning required to grasp why it really works. We are going to progress from classical optimality cuts, to infeasible subproblems and feasibility cuts, and eventually to logic-based Benders decomposition for combinatorial issues.
On this first article, we start with essentially the most accessible attainable setting: an issue wherein acquiring a possible resolution is easy, however figuring out the optimum one stays difficult. We are going to first clarify the overall mechanics of Benders decomposition and the roles performed by the grasp downside, the subproblem, and the optimality cuts. We are going to then work by a small toy downside by hand, deriving every minimize step-by-step till the decrease and higher bounds converge to the optimum resolution. Lastly, we are going to apply the identical equipment to the uncapacitated facility location downside and implement the whole algorithm in Python utilizing Pyomo and the open-source HiGHS solver.
1. The fundamental concept behind Benders decomposition
Many optimization issues mix two basically various kinds of choices. The primary are strategic choices that decide the construction of the answer and are sometimes binary or integer. The second are operational choices that turn into a lot simpler to optimize as soon as the strategic selections are recognized.
Facility location supplies a pure instance. Deciding which services to open is a strategic determination, whereas assigning clients to the chosen services is an operational determination. Benders decomposition separates these two layers as a substitute of fixing them concurrently.
Think about an optimization downside of the next type:
topic to
The variables symbolize the troublesome strategic choices, whereas accommodates the operational choices. If we quickly repair the strategic variables at some worth (), the remaining optimization downside turns into:
Topic to
The operate represents the perfect operational price related to a strategic determination (). Subsequently, the unique downside might be seen as:
The issue is that the grasp downside doesn’t initially know the operate . Benders decomposition introduces a variable () to approximate it:
The grasp and subproblem are then solved iteratively, with every subproblem resolution offering new details about . The subsequent query is the place this data comes from and the way it may be expressed as a sound constraint.
2. The place do Benders cuts come from?
For a set strategic determination (), the operational subproblem is:
topic to
As a result of it is a linear program, it has the next twin:
topic to
Each possible twin resolution () supplies a decrease sure on the operational price:
Since () represents the grasp downside’s approximation of , we are able to impose:
This inequality is a Benders optimality minimize. Its key benefit is that it doesn’t solely describe the operational price of the present grasp resolution. It stays legitimate for different values of (), permitting one subproblem resolve to enhance the grasp’s understanding of a number of attainable strategic choices.
3. A small instance, solved step-by-step
The earlier sections launched the formal construction of Benders decomposition and defined the place optimality cuts come from. Nevertheless, these concepts can nonetheless really feel summary when they’re introduced solely by basic notation. To make the mechanism extra concrete, we are going to now work by a really small instance by hand and comply with the algorithm one iteration at a time.
The aim of this instance is to not resolve a practical optimization downside. It’s to make the interplay between the grasp downside and the subproblem utterly seen. We are going to see how the grasp begins with a very optimistic estimate, how the subproblem reveals the true operational price of that call, and the way the twin resolution generates a minimize that makes the grasp barely wiser within the subsequent iteration.
3.1 The toy optimization downside
So, contemplate the next small optimization downside:
topic to
We are able to interpret () and () as two strategic funding choices. Activating the primary possibility prices 6 models, whereas activating the second prices 4 models. The variables () and () symbolize the remaining operational necessities after these strategic choices have been made, with every unit carrying a value of 1.
3.2 The preliminary grasp downside
Benders decomposition begins by putting the strategic variables () and () within the grasp downside. The operational variables () and () are quickly eliminated and changed by the variable (), which represents an estimate of their mixed price.
The preliminary grasp downside is subsequently:
topic to
At this stage, the grasp has not but acquired any Benders cuts. It solely is aware of that the operational price can’t be adverse. Consequently, its optimum resolution is:
The grasp subsequently proposes making no strategic funding and assumes that the operational price will even be zero. Its goal worth is
This worth supplies a decrease sure on the optimum goal as a result of the grasp is presently underestimating the operational price. The answer shouldn’t be but legitimate for the unique downside, for the reason that operational necessities represented by () and () haven’t been evaluated. The subsequent step is subsequently to ship the choice () to the subproblem and decide its true operational penalties.
3.3 Evaluating the grasp determination
The grasp downside has proposed the strategic determination:
The subproblem is obtained by fixing the strategic variables () and () at these values and optimizing solely the operational variables () and ().
Recall the target operate of the unique downside:
The phrases () symbolize the strategic price and are already accounted for within the grasp downside. Subsequently, the subproblem solely wants to attenuate the operational price:
We now substitute the grasp resolution () and () into every operational constraint. The primary constraint turns into:
The second constraint turns into:
The entire subproblem is subsequently:
topic to
As a result of each variables have constructive coefficients within the goal operate, the subproblem chooses the smallest possible values:
The operational price related to the grasp determination is subsequently:
The strategic price of selecting ) is zero, so the entire price of this possible resolution to the unique downside is:
This offers us an higher sure of (16), whereas the grasp goal supplies a decrease sure of (0):
3.4 Developing the twin subproblem
The subproblem tells us the operational price related to the present grasp determination. Nevertheless, to generate a Benders minimize, we want greater than the worth (). We’d like an expression that is still legitimate when the grasp later considers completely different values of () and (). This expression comes from the twin of the subproblem.
Earlier than substituting the present grasp resolution, allow us to write the subproblem utilizing the generic strategic variables () and ():
topic to
To assemble the twin, we use two fundamental guidelines:
- Every constraint within the primal subproblem creates one twin variable.
- Every variable within the primal subproblem creates one constraint within the twin.
The primary primal constraint receives the twin variable (), whereas the second receives ():
As a result of the primal is a minimization downside and its constraints use the greater-than-or-equal-to route, each twin variables have to be nonnegative:
The twin goal is constructed by multiplying the right-hand aspect of every primal constraint by its corresponding twin variable. Subsequently, the twin goal is:
We now assemble the twin constraints. There’s one twin constraint for every primal variable.
The variable () seems solely within the first primal constraint, with a coefficient of 1. Its coefficient within the primal goal can also be 1. Subsequently, the corresponding twin constraint is:
Equally, () seems solely within the second primal constraint, additionally with a coefficient of 1, and its goal coefficient is 1. Its twin constraint is subsequently:
Combining these parts, the twin of the overall subproblem is:
topic to
We are able to now consider this twin on the present grasp resolution () and (). Substituting these values into the target offers:
topic to
Each coefficients within the goal are constructive, so the target is maximized by setting each twin variables to their largest possible values:
The optimum twin goal worth is subsequently:
This matches the optimum worth of the primal subproblem, Q(0,0)=16. This equality is a consequence of sturdy duality. As a result of the primal and twin subproblems are each possible linear packages, their optimum goal values are equal. The twin resolution () will now permit us to remodel the operational data obtained at () right into a Benders minimize that is still legitimate for different strategic choices.
3.5 Producing the primary Benders minimize
The optimum twin resolution obtained for the present grasp determination is, Recall that the target of the overall twin subproblem is:
To generate the Benders minimize, we substitute the optimum twin values into this expression:
After accumulating phrases, we receive:
This inequality is the primary Benders optimality minimize. On the present resolution (), it forces () to be no less than 16, which is strictly the operational price discovered by the subproblem. Nevertheless, the minimize additionally describes how this decrease sure adjustments when both strategic possibility is chosen.
We now add the minimize to the grasp downside:
topic to
The grasp can not select () whereas pretending that the operational price is zero. If it makes no funding, the brand new minimize requires (), giving a complete estimated price of 16.
If it selects solely the second possibility, the minimize offers:
The corresponding grasp goal worth is:
The up to date grasp subsequently strikes away from the unique determination and proposes:
Its new goal worth is:
The grasp has now discovered that avoiding strategic funding creates a considerable operational burden. The subsequent step is to ship the brand new determination () to the subproblem and confirm whether or not the estimated operational price () is correct.
3.6 The second iteration and convergence
After including the primary optimality minimize, the grasp proposes the answer, with an goal worth of (10), which turns into the brand new decrease sure
It is very important spotlight, that the up to date grasp is definitely detached between ( with ()) and ( with ()), since each produce an goal worth of (10). The answer chosen might subsequently depend upon the solver’s tie-breaking guidelines. For this rationalization, we are going to proceed with ().
To guage this determination, we as soon as once more repair the strategic variables within the authentic operational constraints of the subproblem. The primary constraint turns into:
whereas the second turns into:
The ensuing subproblem is then:
topic to
The optimum operational resolution is subsequently:
with an operational price of, That is precisely the worth predicted by the grasp by (). Including the strategic price offers a complete possible price of:
This resolution improves the earlier higher sure of (16), so, We now have:
The decrease and higher bounds coincide, which proves that the present resolution () is perfect. The grasp can not discover a resolution costing lower than (10), whereas the subproblem has confirmed {that a} possible resolution with price (10) exists. Benders decomposition has subsequently converged.
For the sake of completeness, the twin subproblem at () is:
topic to
One optimum twin resolution is () and (), which might generate the extra minimize:
Nevertheless, this minimize shouldn’t be required to show optimality as a result of the decrease and higher bounds have already converged.
4. Benders decomposition for the uncapacitated facility location downside
Now that we perceive how Benders optimality cuts are generated, we are able to apply the equipment to one of many quintessential issues in operations analysis: the uncapacitated facility location downside, or UFL.
Facility location issues seem each time a corporation should resolve the place to ascertain assets that may serve geographically distributed demand. Typical functions embrace finding warehouses, distribution facilities, hospitals, emergency providers, knowledge facilities, and charging stations. The strategic determination is the place to open the services, whereas the operational determination is how clients ought to be assigned to the chosen areas.
The uncapacitated model assumes that any open facility can serve any variety of clients. This makes it an particularly handy introduction to Benders decomposition as a result of, supplied that no less than one facility is opened and each buyer might be served by each facility, the task subproblem is at all times possible.
4.1 The uncapacitated facility location mannequin
Let () denote the set of candidate services and () the set of consumers. Opening facility () incurs a set price (), whereas assigning buyer () to facility () incurs a service price ().
We outline the strategic variables:
The operational variables are:
The entire, or monolithic, formulation is then:
topic to
The primary constraint requires each buyer to be totally assigned. The second prevents clients from being assigned to closed services, the third constraints, forces to open no less than one facility. Though the task variables are steady, the construction of the issue ensures that an optimum resolution assigns every buyer to one of many least expensive open services.
4.2 The grasp downside
Within the Benders decomposition, the facility-opening variables stay within the grasp downside. The task variables are eliminated and their price is represented by ():
topic to
Please discover that the constraints above, can be joined along with the Benders cuts generated throughout the completely different iterations of the Benders algorithm. The grasp subsequently decides which services to open whereas regularly studying the ensuing task price.
4.3 The task subproblem
As soon as the grasp produces a facility-opening determination (), we repair these values and resolve the operational task downside:
topic to
4.4 The twin subproblem and optimality minimize
To generate a Benders minimize, we assemble the twin of the task subproblem. Let () be the twin variable related to the task equality for buyer (). Because the corresponding primal constraint is an equality, () is unrestricted in signal.
We rewrite the linking constraint as:
and affiliate it with a nonnegative twin variable ().
The twin subproblem is then:
topic to
Suppose that fixing the twin produces the optimum values () and (). The corresponding Benders optimality minimize is then:
This minimize offers the grasp a brand new decrease sure on the task price. The grasp is then solved once more, a brand new facility configuration is evaluated, and extra cuts are generated till the decrease and higher bounds converge.
5. Implementing Benders decomposition with Pyomo and HiGHS
We now have all of the mathematical equipment required to implement the Benders decomposition algorithm mentioned within the earlier sections. The grasp downside will choose the services to open, the twin subproblem will consider the ensuing customer-assignment price, and the generated optimality cuts will progressively enhance the grasp downside.
As an example the whole process, we are going to resolve a small uncapacitated facility location occasion containing 5 candidate services and twenty clients. The occasion is saved in JSON format in order that the optimization mannequin stays separate from the info and might simply be reused or modified. Each the occasion and the whole pocket book might be downloaded from the accompanying GitHub repository.
The implementation makes use of Pyomo to formulate the optimization fashions and the open-source HiGHS solver to unravel the grasp and subproblems. We are going to first load and put together the info, then assemble every element of the decomposition algorithm earlier than combining them into the whole iterative process.
5.1 Loading and getting ready the occasion
We start by putting in the required packages and importing the libraries used all through the implementation. The JSON occasion is then downloaded straight from the GitHub repository and transformed into the units and dictionaries required by Pyomo.
!pip -q set up pyomo highspy
import requests
import pandas as pd
import matplotlib.pyplot as plt
import pyomo.environ as pyo
from pyomo.decide import SolverFactory, TerminationCondition
# Load the occasion from GitHub
JSON_URL = (
"https://uncooked.githubusercontent.com/ceche1212/"
"Benders_Tutorials_TDS/refs/heads/foremost/Knowledge/"
"ufl_benders_5x20_instance.json"
)
response = requests.get(JSON_URL, timeout=30)
response.raise_for_status()
occasion = response.json()
# Units
services = [facility["id"] for facility in occasion["facilities"]]
clients = [customer["id"] for buyer in occasion["customers"]]
# Facility parameters
fixed_cost = {
facility["id"]: facility["fixed_cost"]
for facility in occasion["facilities"]
}
facility_coordinates = {
facility["id"]: (facility["x"], facility["y"])
for facility in occasion["facilities"]
}
# Buyer parameters
demand = {
buyer["id"]: buyer["demand"]
for buyer in occasion["customers"]
}
customer_coordinates = {
buyer["id"]: (buyer["x"], buyer["y"])
for buyer in occasion["customers"]
}
# Project prices
assignment_cost = {
(file["facility"], file["customer"]): file["cost"]
for file in occasion["assignment_costs"]
}
Earlier than constructing the optimization fashions, we carry out a small validation to verify that every one identifiers are distinctive and that each facility-customer mixture has an task price.
if len(services) != len(set(services)):
increase ValueError("Duplicate facility IDs have been discovered.")
if len(clients) != len(set(clients)):
increase ValueError("Duplicate buyer IDs have been discovered.")
expected_pairs = {
(facility, buyer)
for facility in services
for buyer in clients
}
available_pairs = set(assignment_cost)
if expected_pairs != available_pairs:
missing_pairs = expected_pairs - available_pairs
extra_pairs = available_pairs - expected_pairs
increase ValueError(
f"Lacking task prices: {sorted(missing_pairs)}n"
f"Surprising task prices: {sorted(extra_pairs)}"
)
print("Occasion loaded efficiently")
print(f"Identify: {occasion['instance_name']}")
print(f"Services: {len(services)}")
print(f"Clients: {len(clients)}")
print(f"Project prices: {len(assignment_cost)}")The JSON knowledge is reworked into peculiar Python lists and dictionaries. The ability and buyer coordinates will later be used to visualise the ultimate resolution, whereas fixed_cost and assignment_cost include the parameters required by the grasp and subproblems.
5.2 Constructing the preliminary grasp downside
We first assemble the grasp downside, which accommodates solely the strategic facility-opening choices and the variable (), representing the present approximation of the customer-assignment price. The ConstraintList is initially empty as a result of the Benders optimality cuts can be generated and added throughout the iterative algorithm.
# ============================================================
# Construct the preliminary grasp downside
# ============================================================
grasp = pyo.ConcreteModel()
grasp.FACILITIES = pyo.Set(
initialize=services,
ordered=True
)
grasp.fixed_cost = pyo.Param(
grasp.FACILITIES,
initialize=fixed_cost
)
# Facility-opening choices and assignment-cost approximation
grasp.x = pyo.Var(
grasp.FACILITIES,
area=pyo.Binary
)
grasp.theta = pyo.Var(
area=pyo.NonNegativeReals
)
# Not less than one facility have to be opened
grasp.OpenAtLeastOne = pyo.Constraint(
expr=sum(grasp.x[i] for i in grasp.FACILITIES) >= 1
)
# Optimality cuts can be added throughout the algorithm
grasp.BendersCuts = pyo.ConstraintList()
grasp.TotalCost = pyo.Goal(
expr=(
sum(
grasp.fixed_cost[i] * grasp.x[i]
for i in grasp.FACILITIES
)
+ grasp.theta
),
sense=pyo.decrease
)
# ============================================================
# Configure and resolve with HiGHS
# ============================================================
solver = SolverFactory("appsi_highs")
solver.highs_options.replace({
"output_flag": True,
"mip_rel_gap": 0.0,
"time_limit": 360
})
outcomes = solver.resolve(grasp, tee=True)
if outcomes.solver.termination_condition != TerminationCondition.optimum:
increase RuntimeError(
"The preliminary grasp downside was not solved to optimality."
)
# ============================================================
# Extract the preliminary resolution
# ============================================================
x_solution = {
i: int(spherical(pyo.worth(grasp.x[i])))
for i in grasp.FACILITIES
}
theta_value = pyo.worth(grasp.theta)
fixed_cost_value = sum(
fixed_cost[i] * x_solution[i]
for i in grasp.FACILITIES
)
master_objective_value = pyo.worth(grasp.TotalCost)
print("nInitial grasp resolution")
for i, worth in x_solution.gadgets():
standing = "Open" if worth else "Closed"
print(f"Facility {i}: x = {worth} ({standing})")
print(f"nFixed opening price: {fixed_cost_value:,.2f}")
print(f"Estimated task price: {theta_value:,.2f}")
print(f"Grasp goal worth: {master_objective_value:,.2f}")Preliminary grasp resolution
Facility F1: x = 0 (Closed)
Facility F2: x = 0 (Closed)
Facility F3: x = 0 (Closed)
Facility F4: x = 1 (Open)
Facility F5: x = 0 (Closed)
Fastened opening price: 2,100.00
Estimated task price: 0.00
Grasp goal worth: 2,100.00The preliminary grasp has no details about the precise task price past its nonnegativity, so it units (). It subsequently opens solely the power with the bottom fastened price, producing an optimistic decrease sure for the unique downside. The twin subproblem will now consider this facility configuration and generate the primary optimality minimize.
5.3 Evaluating the preliminary determination and producing the primary minimize
The preliminary grasp resolution opens the power with the bottom fastened price whereas assuming that buyer assignments price nothing. We now consider this determination utilizing the twin task subproblem derived earlier. The operate under receives the present facility-opening vector (), solves the twin downside, and returns the knowledge wanted to assemble an optimality minimize.
# ============================================================
# Construct the twin task subproblem
# ============================================================
def build_dual_subproblem(x_bar):
twin = pyo.ConcreteModel()
twin.FACILITIES = pyo.Set(
initialize=services,
ordered=True
)
twin.CUSTOMERS = pyo.Set(
initialize=clients,
ordered=True
)
twin.assignment_cost = pyo.Param(
twin.FACILITIES,
twin.CUSTOMERS,
initialize=assignment_cost
)
twin.x_bar = pyo.Param(
twin.FACILITIES,
initialize=x_bar
)
# alpha_j is unrestricted; beta_ij is nonnegative
twin.alpha = pyo.Var(
twin.CUSTOMERS,
area=pyo.Reals
)
twin.beta = pyo.Var(
twin.FACILITIES,
twin.CUSTOMERS,
area=pyo.NonNegativeReals
)
twin.DualFeasibility = pyo.Constraint(
twin.FACILITIES,
twin.CUSTOMERS,
rule=lambda twin, i, j: (
twin.alpha[j] - twin.beta[i, j]
<= twin.assignment_cost[i, j]
)
)
twin.AssignmentCost = pyo.Goal(
expr=(
sum(
twin.alpha[j]
for j in twin.CUSTOMERS
)
-
sum(
twin.x_bar[i] * twin.beta[i, j]
for i in twin.FACILITIES
for j in twin.CUSTOMERS
)
),
sense=pyo.maximize
)
return twin# ============================================================
# Resolve the twin for the preliminary grasp resolution
# ============================================================
twin = build_dual_subproblem(x_solution)
dual_results = solver.resolve(twin, tee=False)
if dual_results.solver.termination_condition != TerminationCondition.optimum:
increase RuntimeError(
"The twin subproblem was not solved to optimality."
)
alpha_solution = {
j: pyo.worth(twin.alpha[j])
for j in twin.CUSTOMERS
}
beta_solution = {
(i, j): pyo.worth(twin.beta[i, j])
for i in twin.FACILITIES
for j in twin.CUSTOMERS
}
assignment_cost_value = pyo.worth(
twin.AssignmentCost
)
print("nDual subproblem resolution")
print(f"Project price: {assignment_cost_value:,.2f}")
print("nAlpha values:")
for j, worth in alpha_solution.gadgets():
print(f" {j}: {worth:,.2f}")# ============================================================
# Add the primary Benders optimality minimize
# ============================================================
cut_expression = (
sum(
alpha_solution[j]
for j in clients
)
-
sum(
beta_solution[i, j] * grasp.x[i]
for i in services
for j in clients
)
)
grasp.BendersCuts.add(
grasp.theta >= cut_expression
)
print("First Benders optimality minimize added efficiently.")The parameters x_bar repair the power choices proposed by the grasp. The variables () correspond to the customer-assignment equalities, whereas () correspond to the constraints stopping assignments to closed services. Fixing the twin supplies each the true task price of the present configuration and the coefficients required for the minimize. The minimize is then added to grasp.BendersCuts. When the grasp is solved once more, it might probably not assign an unrealistically low worth to () for a similar facility configuration.
5.4 Automating the Benders iterations
The earlier cells carried out the primary Benders iteration explicitly. We are able to now automate the identical sequence: resolve the grasp downside, consider its facility configuration with the twin subproblem, replace the decrease and higher bounds, and add a brand new optimality minimize. The process stops when each bounds coincide, proving that the present resolution is perfect.
# ============================================================
# Benders decomposition algorithm
# ============================================================
maximum_iterations = 100
tolerance = 1e-6
lower_bound = -float("inf")
upper_bound = float("inf")
best_x_solution = None
best_fixed_cost = None
best_assignment_cost = None
best_total_cost = None
iteration_results = []
# Suppress the HiGHS log throughout the iterative process
solver.highs_options["output_flag"] = False
solver.highs_options["log_to_console"] = False
for iteration in vary(1, maximum_iterations + 1):
# --------------------------------------------------------
# 1. Resolve the grasp downside
# --------------------------------------------------------
master_results = solver.resolve(grasp, tee=False)
if master_results.solver.termination_condition != TerminationCondition.optimum:
increase RuntimeError(
f"Grasp downside failed at iteration {iteration}."
)
x_bar = {
i: int(spherical(pyo.worth(grasp.x[i])))
for i in grasp.FACILITIES
}
theta_value = pyo.worth(grasp.theta)
fixed_cost_value = sum(
fixed_cost[i] * x_bar[i]
for i in services
)
lower_bound = pyo.worth(grasp.TotalCost)
# --------------------------------------------------------
# 2. Resolve the twin subproblem
# --------------------------------------------------------
twin = build_dual_subproblem(x_bar)
dual_results = solver.resolve(twin, tee=False)
if dual_results.solver.termination_condition != TerminationCondition.optimum:
increase RuntimeError(
f"Twin subproblem failed at iteration {iteration}."
)
assignment_cost_value = pyo.worth(twin.AssignmentCost)
total_cost = fixed_cost_value + assignment_cost_value
# --------------------------------------------------------
# 3. Replace the perfect possible resolution
# --------------------------------------------------------
if total_cost < upper_bound:
upper_bound = total_cost
best_x_solution = dict(x_bar)
best_fixed_cost = fixed_cost_value
best_assignment_cost = assignment_cost_value
best_total_cost = total_cost
absolute_gap = upper_bound - lower_bound
relative_gap = absolute_gap / max(1.0, abs(upper_bound))
open_facilities = [
i for i in facilities
if x_bar[i] == 1
]
iteration_results.append({
"iteration": iteration,
"open services": ", ".be a part of(open_facilities),
"fastened price": fixed_cost_value,
"theta": theta_value,
"task price": assignment_cost_value,
"decrease sure": lower_bound,
"higher sure": upper_bound,
"absolute hole": absolute_gap,
"relative hole": relative_gap
})
print("-" * 70)
print(f"Iteration: {iteration}")
print(f"Open services: {open_facilities}")
print(f"Fastened price: {fixed_cost_value:,.2f}")
print(f"Theta: {theta_value:,.2f}")
print(f"Precise task price: {assignment_cost_value:,.2f}")
print(f"Decrease sure: {lower_bound:,.2f}")
print(f"Higher sure: {upper_bound:,.2f}")
print(f"Relative hole: {100 * relative_gap:.6f}%")
# --------------------------------------------------------
# 4. Test convergence
# --------------------------------------------------------
if absolute_gap <= tolerance:
print("nBenders decomposition converged.")
break
# --------------------------------------------------------
# 5. Generate and add a brand new optimality minimize
# --------------------------------------------------------
alpha_solution = {
j: pyo.worth(twin.alpha[j])
for j in twin.CUSTOMERS
}
beta_solution = {
(i, j): pyo.worth(twin.beta[i, j])
for i in twin.FACILITIES
for j in twin.CUSTOMERS
}
cut_expression = (
sum(alpha_solution[j] for j in clients)
-
sum(
beta_solution[i, j] * grasp.x[i]
for i in services
for j in clients
)
)
grasp.BendersCuts.add(
grasp.theta >= cut_expression
)
else:
increase RuntimeError(
"Most variety of iterations reached with out convergence."
)After 9 iterations, the Benders algorithm converges to the optimum resolution. At that time, the decrease sure supplied by the grasp downside and the higher sure obtained from the evaluated subproblem attain the identical worth, which proves optimality.
----------------------------------------------------------------------
Iteration: 1
Open services: ['F2', 'F3', 'F5']
Fastened price: 7,000.00
Theta: 56.00
Precise task price: 8,955.00
Decrease sure: 7,056.00
Higher sure: 15,955.00
Relative hole: 55.775619%
----------------------------------------------------------------------
Iteration: 2
Open services: ['F1', 'F2']
Fastened price: 4,600.00
Theta: 7,306.00
Precise task price: 13,781.00
Decrease sure: 11,906.00
Higher sure: 15,955.00
Relative hole: 25.377625%
----------------------------------------------------------------------
Iteration: 3
Open services: ['F2', 'F4', 'F5']
Fastened price: 6,600.00
Theta: 6,403.00
Precise task price: 8,143.00
Decrease sure: 13,003.00
Higher sure: 14,743.00
Relative hole: 11.802211%
----------------------------------------------------------------------
Iteration: 4
Open services: ['F1', 'F3', 'F4', 'F5']
Fastened price: 9,300.00
Theta: 4,206.00
Precise task price: 5,455.00
Decrease sure: 13,506.00
Higher sure: 14,743.00
Relative hole: 8.390423%
----------------------------------------------------------------------
Iteration: 5
Open services: ['F3', 'F5']
Fastened price: 4,800.00
Theta: 8,955.00
Precise task price: 13,792.00
Decrease sure: 13,755.00
Higher sure: 14,743.00
Relative hole: 6.701485%
----------------------------------------------------------------------
Iteration: 6
Open services: ['F2', 'F3', 'F4']
Fastened price: 6,800.00
Theta: 7,057.00
Precise task price: 8,750.00
Decrease sure: 13,857.00
Higher sure: 14,743.00
Relative hole: 6.009632%
----------------------------------------------------------------------
Iteration: 7
Open services: ['F2', 'F5']
Fastened price: 4,500.00
Theta: 9,594.00
Precise task price: 11,243.00
Decrease sure: 14,094.00
Higher sure: 14,743.00
Relative hole: 4.402089%
----------------------------------------------------------------------
Iteration: 8
Open services: ['F1', 'F5']
Fastened price: 4,700.00
Theta: 9,594.00
Precise task price: 12,241.00
Decrease sure: 14,294.00
Higher sure: 14,743.00
Relative hole: 3.045513%
----------------------------------------------------------------------
Iteration: 9
Open services: ['F2', 'F4', 'F5']
Fastened price: 6,600.00
Theta: 8,143.00
Precise task price: 8,143.00
Decrease sure: 14,743.00
Higher sure: 14,743.00
Relative hole: 0.000000%
Benders decomposition converged.
5.5 Recovering the ultimate facility and buyer choices
As soon as the bounds converge, we are able to report the optimum facility configuration and its price. The twin subproblem supplies the task price and the coefficients required to generate cuts, however it doesn’t straight return the primal task variables (). We subsequently resolve the primal task subproblem one ultimate time, fixing the services at the perfect resolution discovered by Benders.
# ============================================================
# Show the ultimate Benders resolution
# ============================================================
print("nFinal Benders resolution")
print("=" * 50)
for i in services:
standing = "Open" if best_x_solution[i] else "Closed"
print(f"Facility {i}: {standing}")
print(f"nFixed opening price: {best_fixed_cost:,.2f}")
print(f"Project price: {best_assignment_cost:,.2f}")
print(f"Complete price: {best_total_cost:,.2f}")
print(f"Iterations: {iteration}")# ============================================================
# Construct the ultimate primal task subproblem
# ============================================================
assignment_model = pyo.ConcreteModel()
assignment_model.FACILITIES = pyo.Set(
initialize=services,
ordered=True
)
assignment_model.CUSTOMERS = pyo.Set(
initialize=clients,
ordered=True
)
assignment_model.assignment_cost = pyo.Param(
assignment_model.FACILITIES,
assignment_model.CUSTOMERS,
initialize=assignment_cost
)
assignment_model.x_bar = pyo.Param(
assignment_model.FACILITIES,
initialize=best_x_solution
)
assignment_model.y = pyo.Var(
assignment_model.FACILITIES,
assignment_model.CUSTOMERS,
area=pyo.NonNegativeReals
)
assignment_model.AssignEachCustomer = pyo.Constraint(
assignment_model.CUSTOMERS,
rule=lambda mannequin, j: sum(
mannequin.y[i, j]
for i in mannequin.FACILITIES
) == 1
)
assignment_model.UseOnlyOpenFacilities = pyo.Constraint(
assignment_model.FACILITIES,
assignment_model.CUSTOMERS,
rule=lambda mannequin, i, j: (
mannequin.y[i, j] <= mannequin.x_bar[i]
)
)
assignment_model.TotalAssignmentCost = pyo.Goal(
expr=sum(
assignment_model.assignment_cost[i, j]
* assignment_model.y[i, j]
for i in assignment_model.FACILITIES
for j in assignment_model.CUSTOMERS
),
sense=pyo.decrease
)
assignment_results = solver.resolve(
assignment_model,
tee=False
)
if (
assignment_results.solver.termination_condition
!= TerminationCondition.optimum
):
increase RuntimeError(
"The ultimate task downside was not solved to optimality."
)# ============================================================
# Extract and show the shopper assignments
# ============================================================
assignment_rows = [
{
"customer": j,
"facility": i,
"assignment": pyo.value(assignment_model.y[i, j]),
"demand": demand[j],
"task price": assignment_cost[i, j]
}
for j in clients
for i in services
if pyo.worth(assignment_model.y[i, j]) > 1e-6
]
assignment_table = pd.DataFrame(assignment_rows)
show(
assignment_table.fashion.format({
"task": "{:.2f}",
"demand": "{:,.0f}",
"task price": "{:,.2f}"
})
)The ultimate output identifies the services chosen by the grasp downside and separates the entire goal into fastened opening and customer-assignment prices. We then repair these facility choices within the primal subproblem and get better the corresponding () values. The ensuing desk exhibits which open facility serves every buyer, along with the related demand and task price.
| buyer | facility | task | demand | task price | |
|---|---|---|---|---|---|
| 0 | C1 | F2 | 1.00 | 2 | 537.00 |
| 1 | C2 | F2 | 1.00 | 3 | 562.00 |
| 2 | C3 | F2 | 1.00 | 1 | 215.00 |
| 3 | C4 | F2 | 1.00 | 4 | 1,024.00 |
| 4 | C5 | F2 | 1.00 | 2 | 113.00 |
| 5 | C6 | F2 | 1.00 | 3 | 129.00 |
| 6 | C7 | F2 | 1.00 | 2 | 195.00 |
| 7 | C8 | F2 | 1.00 | 1 | 115.00 |
| 8 | C9 | F2 | 1.00 | 3 | 840.00 |
| 9 | C10 | F2 | 1.00 | 2 | 691.00 |
| 10 | C11 | F2 | 1.00 | 4 | 1,152.00 |
| 11 | C12 | F2 | 1.00 | 1 | 437.00 |
| 12 | C13 | F4 | 1.00 | 3 | 206.00 |
| 13 | C14 | F4 | 1.00 | 2 | 113.00 |
| 14 | C15 | F4 | 1.00 | 4 | 525.00 |
| 15 | C16 | F4 | 1.00 | 1 | 131.00 |
| 16 | C17 | F5 | 1.00 | 2 | 138.00 |
| 17 | C18 | F5 | 1.00 | 3 | 170.00 |
| 18 | C19 | F5 | 1.00 | 2 | 262.00 |
| 19 | C20 | F5 | 1.00 | 4 | 588.00 |
The task desk is helpful for verification, however a determine supplies a a lot clearer overview of the answer. Within the subsequent cell, we create a easy plot exhibiting which services are open, that are closed, and the way clients are assigned to the chosen services.
# ============================================================
# Visualize the ultimate facility-location resolution
# ============================================================
fig, ax = plt.subplots(figsize=(12, 8))
open_facilities = [i for i in facilities if best_x_solution[i] == 1]
closed_facilities = [i for i in facilities if best_x_solution[i] == 0]
# Plot task connections
for _, row in assignment_table.iterrows():
facility, buyer = row["facility"], row["customer"]
fx, fy = facility_coordinates[facility]
cx, cy = customer_coordinates[customer]
ax.plot(
[fx, cx], [fy, cy],
linewidth=1.2,
alpha=0.55 * row["assignment"],
zorder=1
)
# Plot clients
ax.scatter(
[customer_coordinates[j][0] for j in clients],
[customer_coordinates[j][1] for j in clients],
marker="o",
s=70,
label="Clients",
zorder=3
)
# Plot open services
ax.scatter(
[facility_coordinates[i][0] for i in open_facilities],
[facility_coordinates[i][1] for i in open_facilities],
marker="*",
s=350,
label="Open services",
zorder=4
)
# Plot closed services
ax.scatter(
[facility_coordinates[i][0] for i in closed_facilities],
[facility_coordinates[i][1] for i in closed_facilities],
marker="X",
s=150,
alpha=0.45,
label="Closed services",
zorder=2
)
# Add buyer labels
for buyer in clients:
x, y = customer_coordinates[customer]
ax.annotate(
buyer, (x, y),
xytext=(5, 5),
textcoords="offset factors",
fontsize=9
)
# Add facility labels
for facility in services:
x, y = facility_coordinates[facility]
standing = "Open" if best_x_solution[facility] else "Closed"
ax.annotate(
f"{facility}n{standing}",
(x, y),
xytext=(7, 7),
textcoords="offset factors",
fontsize=10,
fontweight="daring"
)
ax.set_title("Optimum Uncapacitated Facility Location Resolution", fontsize=15)
ax.set_xlabel("X coordinate")
ax.set_ylabel("Y coordinate")
ax.set_aspect("equal", adjustable="field")
ax.grid(alpha=0.25)
ax.legend()
plt.tight_layout()
plt.present()The plot supplies a hen’s-eye view of the optimum resolution. Open services are proven with star markers, closed services with X markers, and clients with circles. The connecting segments point out which facility serves every buyer, making it straightforward to interpret the geographic construction of the answer at a look.
Conclusions
On this article, we explored the central mechanism of classical Benders decomposition. The grasp downside chosen which services to open, the subproblem evaluated the ensuing task price, and the twin resolution generated optimality cuts that regularly corrected the grasp’s initially optimistic estimate. We then applied the whole algorithm in Python utilizing Pyomo and the open-source HiGHS solver.
Benders decomposition is an awfully highly effective approach, however it’s not assured to outperform the equal monolithic formulation. For a small occasion, a contemporary solver might resolve the whole mannequin virtually instantly. In that state of affairs, developing separate grasp and subproblems, transferring options between them, producing cuts, and repeatedly fixing the fashions might introduce extra computational overhead than profit. Decomposition ought to subsequently not be adopted just because it seems extra subtle.
It’s not a silver bullet both. Benders decomposition requires an issue with an exploitable construction, the place fixing a set of complicating variables leaves a subproblem that’s considerably simpler to unravel. Even when this construction exists, efficiency relies upon closely on the standard of the generated cuts. Weak cuts might present little or no data, forcing the algorithm to finish many iterations earlier than the grasp obtains a sufficiently correct illustration of the subproblem. In follow, a lot of the artwork of designing an efficient Benders algorithm lies in producing stronger cuts, eradicating redundant ones, producing a number of cuts, or exploiting further problem-specific information.
The implementation introduced on this article is usually referred to as an outer-loop Benders algorithm. Python explicitly solves the grasp, solves the subproblem, provides a minimize, after which solves the grasp once more. This method is extraordinarily helpful for studying as a result of each step of the algorithm stays seen. Nevertheless, the grasp is repeatedly handed again to the solver as a substitute of permitting the whole process to function straight inside a single branch-and-cut search.
Industrial solvers corresponding to Gurobi and CPLEX assist callback mechanisms by which cuts might be generated whereas the solver is exploring the branch-and-cut tree. Moderately than repeatedly terminating and restarting the grasp optimization course of, a callback can examine candidate or leisure options and add lazy constraints or person cuts throughout the search. Relying on the issue and the energy of the cuts, this could produce a considerable enchancment in computational efficiency.
Nonetheless, the worth of Benders decomposition shouldn’t be restricted to fixing an issue sooner. Generally decomposition makes it attainable to unravel an issue that can’t even be represented as a monolithic mannequin on the accessible {hardware}.
I’m presently encountering exactly this case in certainly one of my analysis initiatives involving a variation of the power location downside. On our analysis machine, which has 128 GB of RAM, making an attempt to assemble the whole monolithic mannequin exhausts the accessible reminiscence earlier than the solver may even start optimizing it. The machine merely reviews that there’s inadequate reminiscence to construct the formulation. With Benders decomposition, nevertheless, we don’t have to create and preserve all grasp and operational variables concurrently. By separating the issue and producing solely the knowledge required by the grasp, we are able to deal with cases that might in any other case stay utterly inaccessible.
That is one other aspect of the ability of decomposition. Generally Benders helps us resolve a mannequin extra effectively. Generally it permits us to unravel a mannequin that can’t slot in reminiscence within the first place.
To maintain this primary introduction as pleasant as attainable, we intentionally thought of a setting wherein the subproblem was at all times possible. By requiring no less than one facility to open and assuming that each facility might serve each buyer, every grasp determination might at all times be translated into a sound task. The one query was how costly that task can be.
However Benders decomposition shouldn’t be at all times this properly behaved. What occurs when the grasp proposes a strategic determination for which no possible operational plan exists? In that state of affairs, the subproblem can not return an optimum price as a result of there isn’t any possible resolution to judge. As an alternative, the algorithm requires a unique sort of suggestions: a feasibility minimize that blocks the present determination and, ideally, a whole household of equally inconceivable choices.
That would be the topic of How Benders Decomposition Works, Half II: Feasibility Cuts. We are going to examine the capacitated facility location downside, the place every facility can serve solely a restricted quantity of demand. The grasp might subsequently open inadequate capability and render the task subproblem infeasible. We are going to look at how feasibility cuts permit Benders decomposition to be taught from these failures and proceed looking for an optimum resolution.
I sincerely hope you discovered this text helpful and that it made Benders decomposition really feel rather less mysterious. Be at liberty to go away a remark or present your appreciation with a clap 👏.
You may also comply with the most recent work from Sávila Training and join with me on LinkedIn. All of the code and knowledge used on this article might be discovered within the accompanying GitHub repository.
Thanks for taking the time to learn. See you in Half II.
References
[1] McCloskey, J. F. (1987). U.S. operations analysis in World Struggle II. Operations Analysis, 35(6), 910-925. DOI: 10.1287/opre.35.6.910.
[2] Dantzig, G. B. (1982). Reminiscences in regards to the origins of linear programming. Operations Analysis Letters, 1(2), 43-48. DOI: 10.1016/0167-6377(82)90043-8.
[3] Dantzig, G. B. (1951). Maximization of a linear operate of variables topic to linear inequalities. In T. C. Koopmans (Ed.), Exercise evaluation of manufacturing and allocation (pp. 339-347). John Wiley & Sons.
[4] Aardal, Ok. I., Hurkens, C. A. J., & Lenstra, J. Ok. (2025). Jacques Benders and his decomposition algorithm. Operations Analysis Letters, 63, Article 107361. DOI: 10.1016/j.orl.2025.107361.
[5] Benders, J. F. (1960). Partitioning in mathematical programming [Doctoral dissertation, Utrecht University].
[6] Benders, J. F. (1962). Partitioning procedures for fixing mixed-variables programming issues. Numerische Mathematik, 4, 238-252. DOI: 10.1007/BF01386316.
[7] Rahmaniani, R., Crainic, T. G., Gendreau, M., & Rei, W. (2017). The Benders decomposition algorithm: A literature overview. European Journal of Operational Analysis, 259(3), 801-817. DOI: 10.1016/j.ejor.2016.12.005.
