Mixed Integer Programming in IronPython QuickStart Sample

Illustrates how to solve mixed integer programming by solving Sudoku puzzles using the linear programming solver in IronPython.

This sample is also available in: C#, Visual Basic, F#.

Overview

This QuickStart sample demonstrates how to use mixed integer programming capabilities to solve complex constraint satisfaction problems, using a Sudoku puzzle as an engaging example.

The sample shows how to formulate the Sudoku puzzle as a mixed integer programming problem by:

Creating binary variables to represent possible digit placements
Expressing Sudoku rules as linear constraints
Setting up constraints for each position, row, column, and 3x3 block
Using bound constraints to fix known values
Solving the system to find a valid solution

The code specifically tackles the puzzle known as “the world’s hardest Sudoku” to demonstrate the power of the approach. It showcases key features of the LinearProgram class including:

Adding binary variables
Creating linear constraints
Setting variable bounds
Solving mixed integer programs

This example illustrates how mathematical programming can be used to solve constraint satisfaction problems that might be difficult to solve using traditional algorithmic approaches.

The code

import numerics

from System import Array

# The linear programming classes reside in their own namespace.
from Extreme.Mathematics.Optimization import *
# Vectors and matrices are in the Extreme.Mathematics namespace
from Extreme.Mathematics import *

# Illustrates solving mixed integer programming problems
# using the classes in the Extreme.Mathematics.Optimization
# namespace of Extreme Numerics.NET.

# In this QuickStart sample, we'll use the Mixed Integer
# programming capabilities to solve Sudoku puzzles.
# The rules of Sudoku will be4 expressed in terms of
# linear constraints on binary variables.

# First, create an empty linear program.
lp = LinearProgram()

# Create an array of binary variables that indicate whether
# the cell at a specific row and column contain a specific digit.
# - The first index corresponds to the row.
# - The second index corresponds to the column.
# - The third index corresponds to the digit.
variables = Array.CreateInstance(LinearProgramVariable, 9, 9, 9)

# Create a binary variable for each digit in each row and column.
# The AddBinaryVariable method creates a variable that can have values of 0 or 1.
for row in range(9):
    for column in range(9):
        for digit in range(9):
            variables[row, column, digit] = lp.AddBinaryVariable("x{0}{1}{2}".format(row, column, digit), 0.0)

# To add integer variables, you can use the AddIntegerVariable method.
# To add real variables, you can use the AddVariable method.

# Now add constraints that represent the rules of Sudoku.

# There are 4 rules in Sudoku. They are all of the kind
# where only one of a certain set of combinations 
# of (row, column, digit) can occur at the same time.
# We can express this by stating that the sum of the corresponding
# binary variables must be one.

# AddConstraints is a helper function defined below.
# For each combination of the first two arguments, 
# it builds a constraint by iterating over the third argument.
coefficients = Vector([ 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 ])
def AddConstraints(lp, variable):
    for i in range(9):
        for j in range(9):
            variables = Array.CreateInstance(LinearProgramVariable, 9)
            for k in range(9):
                variables[k] = variable(i, j, k)
            lp.AddLinearConstraint(variables, coefficients, ConstraintType.Equal, 1.0)

# Rule 1: each posiion contains exactly one digit
AddConstraints(lp, lambda row, column, digit: variables[row, column, digit])
# Rule 2: each digit appears once in each row
AddConstraints(lp, lambda row, digit, column: variables[row, column, digit])
# Rule 3: each digit appears once in each column
AddConstraints(lp, lambda column, digit, row: variables[row, column, digit])
# Rule 4: each digit appears exactly once in each block
AddConstraints(lp, lambda block, digit, index: \
    variables[3 * (block % 3) + (index % 3), 3 * (block / 3) + (index / 3), digit])

# We represent the board with a 9x9 sparse matrix.
# The nonzero entries correspond to the numbers
# already on the board.

# Let's see if we can solve "the world's hardest Sudoku" puzzle:
# http:#www.mirror.co.uk/fun-games/sudoku/2010/08/19/world-s-hardest-sudoku-can-you-solve-dr-arto-inkala-s-puzzle-115875-22496946/
rows = Array[int]([ 0,0,1,1,2,2,2,3,3,3,4,4,4,5,5,5,6,6,6,7,7,8,8 ])
columns = Array[int]([ 2,3,0,7,1,4,6,0,5,6,1,4,8,2,3,7,1,3,8,2,7,5,6 ])
digits = Array[float]([ 5,3,8,2,7,1,5,4,5,3,1,7,6,3,2,8,6,5,9,4,3,9,7 ])
board = Matrix.CreateSparse(9, 9, rows, columns, digits)

# Now fix the variables for the for the digits that are already on the board.
# We do this by setting the lower bound equal to the upper bound:
for triplet in board.NonzeroComponents:
    variables[triplet.Row, triplet.Column, triplet.Value - 1].LowerBound = 1.0

# Solve the linear program.
solution = lp.Solve()

# Scan the variables and print the digit if the value is 1.
for row in range(9):
    for column in range(9):
        theDigit = 0
        for digit in range(9):
            if variables[row, column, digit].Value == 1.0:
                theDigit = digit + 1
                break
        print theDigit.ToString() if theDigit > 0 else ".",
    print