scipy.optimize.milp#
- scipy.optimize.milp(c, *, integrality=None, bounds=None, constraints=None, options=None)[source]#
Mixed-integer linear programming
Solves problems of the following form:
\[\begin{split}\min_x \ & c^T x \\ \mbox{such that} \ & b_l \leq A x \leq b_u,\\ & l \leq x \leq u, \\ & x_i \in \mathbb{Z}, i \in X_i\end{split}\]where \(x\) is a vector of decision variables; \(c\), \(b_l\), \(b_u\), \(l\), and \(u\) are vectors; \(A\) is a matrix, and \(X_i\) is the set of indices of decision variables that must be integral. (In this context, a variable that can assume only integer values is said to be “integral”; it has an “integrality” constraint.)
Alternatively, that’s:
minimize:
c @ x
such that:
b_l <= A @ x <= b_u l <= x <= u Specified elements of x must be integers
By default,
l = 0
andu = np.inf
unless specified withbounds
.- Parameters
- c1D array_like
The coefficients of the linear objective function to be minimized. c is converted to a double precision array before the problem is solved.
- integrality1D array_like, optional
Indicates the type of integrality constraint on each decision variable.
0
: Continuous variable; no integrality constraint.1
: Integer variable; decision variable must be an integer within bounds.2
: Semi-continuous variable; decision variable must be within bounds or take value0
.3
: Semi-integer variable; decision variable must be an integer within bounds or take value0
.By default, all variables are continuous. integrality is converted to an array of integers before the problem is solved.
- boundsscipy.optimize.Bounds, optional
Bounds on the decision variables. Lower and upper bounds are converted to double precision arrays before the problem is solved. The
keep_feasible
parameter of theBounds
object is ignored. If not specified, all decision variables are constrained to be non-negative.- constraintssequence of scipy.optimize.LinearConstraint, optional
Linear constraints of the optimization problem. Arguments may be one of the following:
A single
LinearConstraint
objectA single tuple that can be converted to a
LinearConstraint
object asLinearConstraint(*constraints)
A sequence composed entirely of objects of type 1. and 2.
Before the problem is solved, all values are converted to double precision, and the matrices of constraint coefficients are converted to instances of
scipy.sparse.csc_array
. Thekeep_feasible
parameter ofLinearConstraint
objects is ignored.- optionsdict, optional
A dictionary of solver options. The following keys are recognized.
- dispbool (default:
False
) Set to
True
if indicators of optimization status are to be printed to the console during optimization.- presolvebool (default:
True
) Presolve attempts to identify trivial infeasibilities, identify trivial unboundedness, and simplify the problem before sending it to the main solver.
- time_limitfloat, optional
The maximum number of seconds allotted to solve the problem. Default is no time limit.
- dispbool (default:
- Returns
- resOptimizeResult
An instance of
scipy.optimize.OptimizeResult
. The object is guaranteed to have the following attributes.- statusint
An integer representing the exit status of the algorithm.
0
: Optimal solution found.1
: Iteration or time limit reached.2
: Problem is infeasible.3
: Problem is unbounded.4
: Other; see message for details.- successbool
True
when an optimal solution is found andFalse
otherwise.- messagestr
A string descriptor of the exit status of the algorithm.
The following attributes will also be present, but the values may be
None
, depending on the solution status.- xndarray
The values of the decision variables that minimize the objective function while satisfying the constraints.
- funfloat
The optimal value of the objective function
c @ x
.- mip_node_countint
The number of subproblems or “nodes” solved by the MILP solver.
- mip_dual_boundfloat
The MILP solver’s final estimate of the lower bound on the optimal solution.
- mip_gapfloat
The difference between the final objective function value and the final dual bound.
Notes
milp
is a wrapper of the HiGHS linear optimization software [1]. The algorithm is deterministic, and it typically finds the global optimum of moderately challenging mixed-integer linear programs (when it exists).References
- 1
Huangfu, Q., Galabova, I., Feldmeier, M., and Hall, J. A. J. “HiGHS - high performance software for linear optimization.” Accessed 12/25/2021 at https://www.maths.ed.ac.uk/hall/HiGHS/#guide
- 2
Huangfu, Q. and Hall, J. A. J. “Parallelizing the dual revised simplex method.” Mathematical Programming Computation, 10 (1), 119-142, 2018. DOI: 10.1007/s12532-017-0130-5
Examples
Consider the problem at https://en.wikipedia.org/wiki/Integer_programming#Example, which is expressed as a maximization problem of two variables. Since
milp
requires that the problem be expressed as a minimization problem, the objective function coefficients on the decision variables are:>>> c = -np.array([0, 1])
Note the negative sign: we maximize the original objective function by minimizing the negative of the objective function.
We collect the coefficients of the constraints into arrays like:
>>> A = np.array([[-1, 1], [3, 2], [2, 3]]) >>> b_u = np.array([1, 12, 12]) >>> b_l = np.full_like(b_u, -np.inf)
Because there is no lower limit on these constraints, we have defined a variable
b_l
full of values representing negative infinity. This may be unfamiliar to users ofscipy.optimize.linprog
, which only accepts “less than” (or “upper bound”) inequality constraints of the formA_ub @ x <= b_u
. By accepting bothb_l
andb_u
of constraintsb_l <= A_ub @ x <= b_u
,milp
makes it easy to specify “greater than” inequality constraints, “less than” inequality constraints, and equality constraints concisely.These arrays are collected into a single
LinearConstraint
object like:>>> from scipy.optimize import LinearConstraint >>> constraints = LinearConstraint(A, b_l, b_u)
The non-negativity bounds on the decision variables are enforced by default, so we do not need to provide an argument for bounds.
Finally, the problem states that both decision variables must be integers:
>>> integrality = np.ones_like(c)
We solve the problem like:
>>> from scipy.optimize import milp >>> res = milp(c=c, constraints=constraints, integrality=integrality) >>> res.x [1.0, 2.0]
Note that had we solved the relaxed problem (without integrality constraints):
>>> res = milp(c=c, constraints=constraints) # OR: >>> # from scipy.optimize import linprog; res = linprog(c, A, b_u) >>> res.x [1.8, 2.8]
we would not have obtained the correct solution by rounding to the nearest integers.