Solving IPs: It's all about bounding

• Assume we're maximizing (we could also minimize)
• Dual bound: a value we can't do better than (optimal <= this)
• Primal solution: any feasible solution (optimal >= this)
• If dual bound = primal bound: we're optimal

Generating Primal Bounds

• Any feasible solution will do.
• One common method: greedy algorithm
• Greedy for knapsack problem:
  • order the list by descending value
  • put in as many items as will fit
• Sometimes finding a feasible solution is as hard finding the optimal solution...

Generating Dual Bounds

• Typically this means relaxing the problem.
• The LP relaxation of an IP assumes the variables can be fractional. We can do no better than this.
• Important: If the optimal solution of an LP relaxation has integer variables, then that is the optimal solution of the IP!

Cutting Planes

Add constraints that make the LP have (more) integral extreme points.

If all extreme points are integral, then Simplex will return an optimal integer point.

Branch & Bound

If the LP relaxation returns fractional variables, divide it into two LPs.

Discard the fractional area.

Now we have two sub-LPs, each with a dual bound. We can use this to guide our search.

Preprocessing

Sometimes we can make logical implications about integer variables just by looking at the model. Consider:

x1 + x2 >= 1
x1 + x2 + x3 <= 1
xi in {0, 1} for all i

We can set x₃ to zero and remove it from our problem entirely. (Why?)

OK, but how do I really do it?

All these techniques (and many more) exist in commercial and open source optimization solvers.

The ZIB Optimization Suite is a set of non-commercial tools for solving MIPs. It contains:

• ZIMPL: algebraic modeling language like AMPL
• SoPlex: linear programming solver
• SCIP: mixed integer programming solver

ZIMPL Example: SpamCo

The ZIMPL version of our model is a lot like the algebraic one:

var x1 >= 0; # pots of type 1
var x2 >= 0; # pots of type 2

maximize profit:  3*x1 + 4*x2;
subto byproducts: 2*x1 + 3*x2 <= 10;
subto curing:     3*x1 +  x2  <=  6;

Algebraic Modeling Languages

• The quickest way to get started
• Produce models that closely resemble mathematical descriptions of problemsi
• Separate data from model
• Aren't full programming languages. It is difficult to:
  • Change the model depending on solver output
  • Combine data from disparate sources
    (CSV, RDBMS, arbitrarily structured text...)
  • Base parts of the model on logic or user input

Solver APIs

• Give all the existing power of a solver
• Allow you to build your own problem-specific optimization components
• Are usually in C...

SpamCo in Python

from zibopt import scip
solver = scip.solver()

x1 = solver.variable() # pots of type 1
x2 = solver.variable() # pots of type 2

solver += 2*x1 + 3*x2 <= 10 # animal byproducts
solver += 3*x1 +   x2 <=  6 # curing time

solution = solver.maximize(objective=3*x1 + 4*x2)
if solution:
    print 'profit:', solution.objective
    print 'pots of type 1:', solution[x1]
    print 'pots of type 2:', solution[x2]
else:
    print 'invalid problem'

A Larger Example: Sudoku

We can model Sudoku as a Binary Integer Program (BIP) with 9*9*9 = 729 variables.

# 0 indicates a cell value is not given
problem = [
    [0, 0, 0,   6, 9, 2,   0, 4, 0],
    [7, 0, 0,   0, 0, 0,   8, 9, 0],
    [0, 0, 0,   0, 0, 0,   0, 0, 6],

    [0, 0, 9,   0, 1, 7,   0, 0, 3],
    [0, 0, 7,   0, 8, 0,   5, 0, 0],
    [8, 0, 0,   4, 6, 0,   1, 0, 0],

    [5, 0, 0,   0, 0, 0,   0, 0, 0],
    [0, 8, 6,   0, 0, 0,   0, 0, 1],
    [0, 3, 0,   7, 2, 8,   0, 0, 0]
]

Sudoku: Data & Indexes

Indexes for rows, columns, 3x3 groups, and values cells can take on:

# Indexes for creating variables and constraints
rows = range(len(problem))
cols = rows
vals = range(1, 10)

9x9x9 structure for storing our binary variables. .

x = dict((i, dict((j, {}) for j in cols)) for i in rows)

Sudoku: Creating the Variables

Create a variable for each cell and value. If row i, column j requires value k, then x[i][j][k] = 1.

from itertools import product
from zibopt import scip

# Initialize one binary variable per cell value
solver = scip.solver()
for i, j, k in product(rows, cols, vals):
    x[i][j][k] = solver.variable(
        vartype = scip.BINARY,
        lower   = 1 if problem[i][j] == k else 0
    )

Sudoku: Cell Values, Columns, & Rows

Each cell takes exactly one value:

for i, j in product(rows, cols):
    solver += sum(x[i][j][k] for k in vals) == 1

Each value occurs in each row once:

for i, k in product(rows, vals):
    solver += sum(x[i][j][k] for j in cols) == 1

Each value occurs in each column once:

for j, k in product(cols, vals):
    solver += sum(x[i][j][k] for i in rows) == 1

Sudoku: 3x3 Groups

Each 3x3 group contains each value once:

groups = (0, 3, 6)

for m, n, k in product(groups, groups, vals):
    solver += sum(
        x[i][j][k] for i, j in product(
            range(m,m+3), range(n,n+3)
        )
    ) == 1

Sudoku: Finding the Solution

Find the solution and fill in our problem array with appropriate values.

solution = solver.maximize()
if solution:
    for i, j in product(rows, cols):
        for k in vals:
            if solution[x[i][j][k]]:
                problem[i][j] = k
                break

    for row in problem:
        print row

else:
    print 'invalid'

Sudoku: Output

[3, 1, 8, 6, 9, 2, 7, 4, 5]
[7, 6, 4, 3, 5, 1, 8, 9, 2]
[2, 9, 5, 8, 7, 4, 3, 1, 6]
[6, 2, 9, 5, 1, 7, 4, 8, 3]
[1, 4, 7, 2, 8, 3, 5, 6, 9]
[8, 5, 3, 4, 6, 9, 1, 2, 7]
[5, 7, 2, 1, 4, 6, 9, 3, 8]
[4, 8, 6, 9, 3, 5, 2, 7, 1]
[9, 3, 1, 7, 2, 8, 6, 5, 4]

Conclusions

• Just because a problem is intractable doesn't mean you can't solve it (or at least know bounds for the solution you have)
• ZIB has an excellent suite of open source tools for just this purpose
• Optimization is for Python programmers too!

References

• ZIB Optimization Suite
  http://zibopt.zib.de/
• ZIMPL User Guide
  http://zimpl.zib.de/download/zimpl.pdf
• Constraint Integer Programming
  http://opus.kobv.de/tuberlin/volltexte/2007/1611/
• python-zibopt
  http://code.google.com/p/python-zibopt/

Questions