Modeling OR Problems in SAGE

Sage is an open source project providing a comprehensive environment for mathematical computing. It's a big package with big ambitions "Creating a viable free open source alternative to Magma, Maple, Mathematica and Matlab."  See SIAM Review (Vol 41, No. 4, 2009) for a  recent review.

Sage incorporates Python as a scripting language from which one can generate objects MixedIntegerLinearProgram objects.  Both GLPK and CBC from the COIN-OR project can be installed as optional Sage packages and used as solvers.  Coupled with the native features of Python, this approaches provides an interesting alternative to developing models using AMPL/GMPL.

# Dictionary of task durations indexed by (Job, Machine) tuples
dur = {
    ('Paper 1','Blue')   : 45,
    ('Paper 1','Yellow') : 10,
    ('Paper 2','Blue')   : 20,
    ('Paper 2','Green')  : 10,
    ('Paper 2','Yellow') : 34,
    ('Paper 3','Blue')   : 12,
    ('Paper 3','Green')  : 17,
    ('Paper 3','Yellow') : 28 }

# Task sequencing indicated by pairs of (Job, Machine) tuples
order = [
    (('Paper 1','Blue'),  ('Paper 1','Yellow')),
    (('Paper 2','Green'), ('Paper 2','Blue')),
    (('Paper 2','Blue'),  ('Paper 2','Yellow')),
    (('Paper 3','Yellow'),('Paper 3','Blue')),
    (('Paper 3','Blue'),  ('Paper 3','Green')) ]

# TASKS are a list of (job,machine) tuples found from the keys for dur{}
TASKS = dur.keys()

# JOBS and MACHS are unduplicated lists of jobs and machines 
JOBS = sorted(list(set(zip(*TASKS)[0])))
MACHS = sorted(list(set(zip(*TASKS)[1])))

# Create MILP
JobShop = MixedIntegerLinearProgram(maximization=False)

# The objective is to minimize Makespan
MakeSpan = JobShop.new_variable()
JobShop.set_objective(MakeSpan[0])

# MakeSpan is an upper bound on all tasks
x = JobShop.new_variable()
[JobShop.add_constraint(x[t] + dur[t] - MakeSpan[0], max = 0) for t in TASKS]

# Tasks must be complete in the specified order 
[JobShop.add_constraint(x[s] + dur[s] - x[t], max=0) for s,t in order]

# Disjunctive constraints to avoid machine conflicts. This uses a
# 'Big M' technique and a binary variable where y[s][t] = 0 implies
# task s must finish before task t can start and y[s][t] = 1 implies
# task t finishes before task s starts.

CONFLICTS = [(s,t) for s in TASKS for t in TASKS if s<t and s[1]==t[1]]

y = JobShop.new_variable(dim=2,vtype=JobShop.__BINARY)

BigM = sum(dur[t] for t in TASKS)

[JobShop.add_constraint(BigM*y[s][t] + x[t] - x[s] - dur[s], min=0)
    for s,t in CONFLICTS]
[JobShop.add_constraint(BigM*(1-y[s][t]) + x[s] - x[t] - dur[t], min=0)
    for s,t in CONFLICTS]
JobShop.solve('Coin')