A Gentle Introduction to Branch & Bound
The most fundamental integer and mixed-integer programming algorithm explained with Python
Numerical optimization problems are a fundamental tool in quantitative decision-making processes. Suppose a system can be described by a set of mathematical equations that adequately encompasses the impact of decision variables in objectives and constraints. In that case, one might search for the values of those decision variables that produce the best possible outcome using optimization algorithms.
The presence or not of nonlinear functions to describe how the decision variables influence the objectives and constraints is used to divide optimization problems into two major categories: Linear and Nonlinear Programming. Management sciences and operations research make extensive use of linear models, whereas nonlinear programming problems tend to arise naturally in the physical sciences and engineering (Nocedal & Wright, 2006). Throughout this article, the focus will be on linear models.
In some problems, it makes sense that decision variables assume only integer values. For instance, in a personnel scheduling problem, it does not make sense to assign 2.57 employees to a shift; or in a vehicle routing problem, one can not choose to take half a route to go from a to b. These are denoted as either integer or mixed-integer problems.
In practice, most of these problems are solved by the Branch & Bound algorithm combined with real-valued algorithms and other strategies, such as cutting planes, pricing, and custom heuristics. Therefore, understanding how Branch & Bound works might help us to gain insight on how to combine it with other strategies and formulate problems better when solving complex problems.
In the first part of this article, we will see how to formulate a linear programming problem with a visual example of a two-variable problem. In the following part, a step-by-step branch and bound algorithm will be applied to the same problem to obtain the optimal integer solution. All the code used throughout the article is available in this example notebook.
It might be helpful if, before reading this article, the reader is already familiar with Linear Programming, although its main concepts will be briefly presented here too. You can find an interesting overview in this other article.
Definition of a linear programming problem
When formulating an optimization problem, one must define an objective that is a function of a vector decision variables x and might be subject to some equality and inequality constraints, which are functions of x as well. The objective can be defined either in a minimization or maximization sense although the former is the most usual. Notice that by simply multiplying the coefficients cᵢ of a maximization objective, it can be re-formulated in a minimization sense.
In linear problems, the decision variable space must be limited somehow. Otherwise, decision variables would converge to positive or negative infinity values according to the objective function. Constraints might be formulated by either equality or inequality relationships. They are usually stated as rows in the problem matrix. Notice inequality constraints might be formulated as equality constraints by adding non-negative slack variables. However, throughout this article, let us distinguish the matrix of equality constraints A_eq from the matrix of inequality constraints A_ub.
Lower and upper boundaries for each component of x might be explicit in the formulation, which reduces the search space. As a convention, due to solution techniques, usually lower bounds for decision variables are by default equal to zero. This leads to a general problem formulation such as the following.
Let us consider a numerical example to illustrate these concepts.
The feasible space is a region in the decision variable space in which all constraints are valid. Notice that in some situations the constraints might be inconsistent among themselves. Suppose we had one constraint x₁ + x₂ ≤ 2 and another x₁+ x₂ ≥ 3. In such a case, the problem would be considered infeasible. However, in the numerical example previously defined, we can visually represent the feasible space as follows.
As we are considering a maximization problem, it is expected that the best solution occurs within the limits established by the constraints. Let us visualize how the objective function is evaluated in the decision space by a contour plot.
As we can see, the maximum value of the objective is expected to occur in the intersection between the two constraints. In this case, both are thus considered active constraints.
To solve this problem in Python, we will use the linprog function from scipy.optimize. Notice, as in scipy standard form, it will be stated as a minimization problem, multiplying objective coefficients by -1.
import numpy as np
from scipy.optimize import linprog
c = np.array([-5.0, -4.0])
A_ub = np.array(
[[2.0, 3.0],
[2.0, 1.0]]
)
b_ub = np.array([12.0, 6.0])
sol_relaxed = linprog(c, A_ub=A_ub, b_ub=b_ub)This should return a solution object with the properties:
- x: [1.5, 3.0]
- fun: -19.5
In this solution, x₁ has a fraction value, which would not be feasible in some practical situations. In case we needed all variables to assume integer values, we might have stated this condition in the linprog function (at least since scipy version 1.9).
sol_int = linprog(c, A_ub=A_ub, b_ub=b_ub, integrality=np.ones(2))And this should return a solution object with the properties:
- x: [2.0, 2.0]
- fun: -18.0
Notice that the objective value is now worse than it was in the relaxed formulation, as the feasible space has been reduced due to integrality constraints. We can state that the solution of the problem in its relaxed form is necessarily an optimality limit for the integer solution. To understand the Branch & Bound algorithm in the next section let us first establish a simple yet meaningful rule.
If you solve the LP (linear program) relaxation of a pure IP (integer program) and obtain a solution in which all variables are integers, then the optimal solution to the LP relaxation is also the optimal solution to the IP (Winston & Goldberg, 2004).
Now let us dive into integrality constraints.
Integrality and the Branch & Bound algorithm
As previously defined, if in the solution of a linear problem in its relaxed form all integrality constraints are valid, the solution is also optimal to the integer or mixed integer problem. Let us then consider that, in the same example from the previous section, both x₁ and x₂ needed to be integers. In such a situation, we could apply the Branch & Bound algorithm to find its optimal integer solution.
The first step to find the optimal integer solution was performed in the previous section, as we found the solution to the problem in its relaxed form. To apply the Branch & Bound algorithm, we must now partition the feasible space in the fractional variable of the relaxed solution x₁.
Notice that in most real-world problems it is likely that one faces steps in which more than one integer variable is fractional in an intermediate solution. Several heuristics can be then used to define which variable to branch on. For instance, branching on the fractional-valued variable that has the greatest economic importance is often the best strategy (Winston & Goldberg, 2004). In our simple example, this decision is straightforward as only one variable is fractional.
Two subproblems are then created from the previous problem, such that, in neither of them, 1 < x₁ < 2 is a valid solution.
- P1: Original problem with new constraint x₁ ≤ 1.
- P2: Original problem with new constraint x₁ ≥ 2.
A display of all subproblems that have been created is called a tree. Each subproblem is referred to as a node of the tree, and each line connecting two nodes of the tree is called an arc. The constraints associated with any node of the tree are the constraints for the LP relaxation plus the constraints associated with the arcs leading from subproblem 1 (the original problem) to the node. (Winston & Goldberg, 2004). In this structure, the original problem with no integrality restrictions is usually referred to as the root node. Here, I adopted the index 0 to identify it.
In this example, we have so far:
In Python, it can be formulated as the following.
# P1:
A_ub_p1 = np.vstack((A_ub, np.atleast_2d([1.0, 0.0])))
b_ub_p1 = np.append(b_ub, np.floor(sol_relaxed.x[0]))
# P2:
A_ub_p2 = np.vstack((A_ub, np.atleast_2d([-1.0, 0.0])))
b_ub_p2 = np.append(b_ub, -np.ceil(sol_relaxed.x[0]))There are also several strategies to choose which subproblem to solve first. One alternative is to use the Last In First Out strategy, which states that the most recently created subproblems should be solved in advance. However, it still does not define how to define priorities between two subproblems simultaneously created. In this example, for illustration purposes, we will arbitrarily solve P1 before P2, but P2 before subproblems generated from P1.
sol_p1 = linprog(c, A_ub=A_ub_p1, b_ub=b_ub_p1)And now we have another fractional solution [1.0, 3.3] with an objective value of 18.33, which is evidently worse than the relaxed solution. As P2 has not been solved yet, 19.5 is still the best value an integer solution could possibly assume. Let us denote this as the “best bound”.
Next, let us solve P2.
sol_p2 = linprog(c, A_ub=A_ub_p2, b_ub=b_ub_p2)And now we obtained our first integer solution [2.0, 2.0] with an objective value of 18.0. Visually, the problem so far looks like this:
Notice that, when branching on some decision variable, one or both of the subproblems generated might be infeasible. In such situation, the node is no longer of interest and one must choose another node to be explored. The same might occur if the partial solution of a node has a worse objective than the previous best integer solution found.
As now both subproblems derived from the root node have been explored, we know that the optimal integer solution is not worse than 18.0 (our current best integer solution) or better than 18.33 (the best bound). Let us then create subproblems from the node corresponding to P1 to see if we could possibly overcome 18.0. As in the solution of P1, only x₂ assumed fractional values, let us branch on x₂.
- P3: P1 with new constraint x_2 ≤ 3.
- P4: P1 with new constraint x_2 ≥ 4.
In Python,
# P3:
A_ub_p3 = np.vstack((A_ub_p1, np.atleast_2d([0.0, 1.0])))
b_ub_p3 = np.append(b_ub_p1, np.floor(sol_p1.x[1]))
# P4:
A_ub_p4 = np.vstack((A_ub_p1, np.atleast_2d([0.0, -1.0])))
b_ub_p4 = np.append(b_ub_p1, -np.ceil(sol_p1.x[1]))Both these problems lead to integer solutions, [1.0, 3.0] and [0.0, 4.0], although they are both worse than the one obtained from P2. Therefore, we have proved its optimality as the tree was fully explored and no better integer solution was found (nor can possibly be found). And we can visually represent it by the following figures.
Further reading
There are several types of problems that can be solved using mathematical modeling, of which some are formulations more usual. Those interested in learning more about how to describe some usual rules using mathematical expressions can refer to Winston & Goldberg (2004). For a deeper understanding of integer programming algorithms, probably Wolsey (2020) is the best reference.
It is likely that, to describe real-world problems, one should use some algebraic modeling language (AML). Pyomo (Bynum et al., 2021) is an interesting Python alternative to do so, as it is open-source and compatible with several solvers. You can see a quick Pyomo example in my previous article about the multi-dimensional knapsack problem.
Other problems can not be described by linear expressions, such as engineering design, financial modeling, and curve fitting. They lead to very distinct solution techniques compared to linear optimization problems. The most usual alternative to solve them is to use some of the classic convex nonlinear programming algorithms, which I present in this other article.
Yet sometimes it is not sufficient, as one might be in pursuit of multiple objectives, or the problem violates some of convex algorithm assumptions. Then, one should resort to other solution techniques, of which metaheuristics can be much helpful. Pymoo (Blank & Deb, 2020) is an amazing alternative to do so, as several of these algorithms are available there.
Recently, I published the Python library pymoode, an extension of pymoo, which I developed throughout the styrene reactor optimization project (Leite et al., 2022). One can find a complete tutorial in this other article.
I intend to fully incorporate pymoode algorithms into pymoo soon, but for now, they are available as a separate package.
Conclusions
Throughout this article, fundamental concepts of linear integer programming were presented with a general overview of the Branch & Bound algorithm. A simple problem was used to explain the algorithm step-by-step with graphical illustrations. The complete Python code used is available to the reader in an example notebook.
References
Blank, J. & Deb, K., 2020. pymoo: Multi-Objective Optimization in Python. IEEE Access, Volume 8, pp. 89497–89509.
Bynum, M. L. et al., 2021. Pyomo-optimization modeling in python. Springer.
Leite, B., Costa, A., Costa Junior, E., 2022. Multi-objective optimization of adiabatic styrene reactors using Generalized Differential Evolution 3 (GDE3). Chem. Eng. Sci. 118196. doi: 10.1016/j.ces.2022.118196.
Nocedal, J. & Wright, S. J., 2006. Numerical Optimization. 2nd ed. New York: Springer New York.
Winston, W. L. & Goldberg, J. B., 2004. Operations research: applications and algorithms. 4th ed. Belmont, CA: Thomson Brooks/Cole Belmont.
Wolsey, L. A., 2020. Integer Programming. 2nd ed. John Wiley & Sons.
