Navigation

This document explains the use of linear programming (LP) – and of
mixed integer linear programming (MILP) – in Sage by illustrating it
with several problems it can solve. Most of the examples given are
motivated by graph-theoretic concerns, and should be understandable
without any specific knowledge of this field. As a tool in
Combinatorics, using linear programming amounts to understanding how
to reformulate an optimization (or existence) problem through linear
constraints.

Here we present the usual definition of what a linear program is: it
is defined by a matrix \(A: \mathbb{R}^m \mapsto \mathbb{R}^n\), along
with two vectors \(b,c \in \mathbb{R}^n\). Solving a linear program is
searching for a vector \(x\) maximizing an objective function and
satisfying a set of constraints, i.e.

\[c^t x = \max_{x' \text{ such that } Ax' \leq b} c^t x'\]

where the ordering \(u \leq u'\) between two vectors means that the
entries of \(u'\) are pairwise greater than the entries of \(u\). We also
write:

Equivalently, we can also say that solving a linear program amounts to
maximizing a linear function defined over a polytope (preimage or
\(A^{-1} (\leq b)\)). These definitions, however, do not tell us how to
use linear programming in combinatorics. In the following, we will
show how to solve optimization problems like the Knapsack problem, the
Maximum Matching problem, and a Flow problem.

There are bad news coming along with this definition of linear
programming: an LP can be solved in polynomial time. This is indeed
bad news, because this would mean that unless we define LP of
exponential size, we cannot expect LP to solve NP-complete problems,
which would be a disappointment. On a brighter side, it becomes
NP-complete to solve a linear program if we are allowed to specify
constraints of a different kind: requiring that some variables be
integers instead of real values. Such an LP is actually called a “mixed
integer linear program” (some variables can be integers, some other
reals). Hence, we can expect to find in the MILP framework a wide
range of expressivity.

The MILP class in Sage represents a MILP! It is also used to
solve regular LP. It has a very small number of methods, meant to
define our set of constraints and variables, then to read the solution
found by the solvers once computed. It is also possible to export a
MILP defined with Sage to a .lp or .mps file, understood by most
solvers.

In the previous example, we obtained variables through v['x'], v['y']
and v['z']. This being said, larger LP/MILP will require us to associate an
LP variable to many Sage objects, which can be integers, strings, or even the
vertices and edges of a graph. For example:

sage: x=p.new_variable(real=True,nonnegative=True)

With this new object x we can now write constraints using
x[1],...,x[15].

sage: p.add_constraint(x[1]+x[12]-x[14]>=8)

Notice that we did not need to define the “length” of x. Actually, x
would accept any immutable object as a key, as a dictionary would. We can now
write

sage: p.add_constraint(x["I am a valid key"]+....: x[("a",pi)]<=3)

And because any immutable object can be used as a key, doubly indexed variables
\(x^{1,1}, ..., x^{1,15}, x^{2,1}, ..., x^{15,15}\) can be referenced by
x[1,1],...,x[1,15],x[2,1],...,x[15,15]

Types : If you want a variable to assume only integer or binary values, use
the integer=True or binary=True arguments of the new_variable
method. Alternatively, call the set_integer and set_binary methods.

Bounds : If you want your variables to only take nonnegative values, you can
say so when calling new_variable with the argument nonnegative=True. If
you want to set a different upper/lower bound on a variable, add a constraint or
use the set_min, set_max methods.

The Knapsack problem is the following: given a collection of items
having both a weight and a usefulness, we would like to fill a bag
whose capacity is constrained while maximizing the usefulness of the
items contained in the bag (we will consider the sum of the items’
usefulness). For the purpose of this tutorial, we set the restriction
that the bag can only carry a certain total weight.

To achieve this, we have to associate to each object \(o\) of our
collection \(C\) a binary variable taken[o], set to 1 when the
object is in the bag, and to 0 otherwise. We are trying to solve the
following MILP

Given a graph \(G\), a matching is a set of pairwise disjoint edges. The
empty set is a trivial matching. So we focus our attention on maximum
matchings: we want to find in a graph a matching whose cardinality is
maximal. Computing the maximum matching in a graph is a polynomial
problem, which is a famous result of Edmonds. Edmonds’ algorithm is
based on local improvements and the proof that a given matching is
maximum if it cannot be improved. This algorithm is not the hardest to
implement among those graph theory can offer, though this problem can
be modeled with a very simple MILP.

To do it, we need – as previously – to associate a binary variable
to each one of our objects: the edges of our graph (a value of 1
meaning that the corresponding edge is included in the maximum
matching). Our constraint on the edges taken being that they are
disjoint, it is enough to require that, \(x\) and \(y\) being two edges
and \(m_x, m_y\) their associated variables, the inequality \(m_x + m_y
\leq 1\) is satisfied, as we are sure that the two of them cannot both
belong to the matching. Hence, we are able to write the MILP we
want. However, the number of inequalities can be easily decreased by
noticing that two edges cannot be taken simultaneously inside a
matching if and only if they have a common endpoint \(v\). We can then
require instead that at most one edge incident to \(v\) be taken inside
the matching, which is a linear constraint. We will be solving:

Yet another fundamental algorithm in graph theory: maximum flow! It
consists, given a directed graph and two vertices \(s, t\), in sending a
maximum flow from \(s\) to \(t\) using the edges of \(G\), each of them
having a maximal capacity.

The definition of this problem is almost its LP formulation. We are
looking for real values associated to each edge, which would
represent the intensity of flow going through them, under two types of
constraints:

The amount of flow arriving on a vertex (different from \(s\) or \(t\))
is equal to the amount of flow leaving it.

The amount of flow going through an edge is bounded by the capacity
of this edge.

This being said, we have to maximize the amount of flow leaving
\(s\): all of it will end up in \(t\), as the other vertices are sending
just as much as they receive. We can model the flow problem with the
following LP

Sage solves linear programs by calling specific libraries. The
following libraries are currently supported:

CBC: A solver from
COIN-OR,
provided under the Eclipse Public License (EPL), which is an open source
license but incompatible with GPL. CBC and the Sage CBC backend can be
installed using the shell command:

$ sage -i -c sage_numerical_backends_coin

CPLEX:
Proprietary, but available for free for researchers and students through
IBM’s Academic Initiative. Since trac ticket #27790, only versions 12.8 and
above are supported.

Install CPLEX according to the instructions on the
website, which includes obtaining a license key.

Then find the installation directory of your ILOG CPLEX Studio installation, which contains subdirectories cplex, doc, opl, etc.
Set the environment variable CPLEX_HOME to this directory; for example using the following shell command (on macOS):

$ export CPLEX_HOME=/Applications/CPLEX_Studio1210

or (on Linux):

$ export CPLEX_HOME=/opt/ibm/ILOG/CPLEX_Studio1210

Now verify that the CPLEX binary that you will find in the subdirectory
cplex/bin/ARCH-OS starts correctly, for example:

Licensed under the GPLv3. This solver is always installed, as the default one, in Sage.

Gurobi:
Proprietary, but available for free for researchers and students via Gurobi’s
Academic Program.

Install Gurobi according to the instructions on the website,
which includes obtaining a license key. The installation should make the
interactive Gurobi shell gurobi.sh available in your PATH.
Verify this by typing the shell command gurobi.sh: