Optimization modeling using Solver Foundation and Sho
Two days ago, my boss's boss's boss's boss's boss's boss Soma announced the launch of several Technical Computing projects on Dev Labs. One of these projects, Sho, is an Iron Python-based environment for doing math, stats, optimization, and the like. Sho’s optimization libraries are powered by Solver Foundation. The Team of Sho has written an “optimizer” package that wraps the Simplex, Interior Point, and Compact Quasi-Newton solvers in a way that is convenient for Sho programmers. If you have used other scientific computing environments then these matrix-based APIs are pretty easy to learn. Even if you don't care about optimization I encourage you to try it out!
Download Sho from this site, and then start the Sho Console. Now try running these Sho commands to see how you can use the LP object to solve the “Petrochem” example I have beaten to death in previous posts:
help(LP)
lp = LP()
lp.A = DoubleArray.From([[0.3, 0.4], [0.4, 0.2], [0.2, 0.3]])
lp.c = DoubleArray.From([20.0, 15.0])
lp.blo = DoubleArray.From([2000, 1500, 500])
lp.bhi = DoubleArray.From([100000, 100000, 100000])
lp.xlo = DoubleArray.From([0, 0])
lp.xhi = DoubleArray.From([9000, 6000])
lp.Solve()
lp.x
lp.OptVal
But wait, there’s more! The Express version Solver Foundation is included with Sho, and since Solver Foundation is a .Net library it is very easy to import and use within Sho.
ShoLoadAssembly("..\packages\Optimizer\Microsoft.Solver.Foundation.dll")
from Microsoft.SolverFoundation.Services import *
from Microsoft.SolverFoundation.Common import *
c = SolverContext.GetContext()
Let’s take this a step further. I find using some of the Solver Foundation Services classes a bit unwieldy in Iron Python. It “works”, but I have a hard time getting anything done in real life. So this morning I wrote a simple “model” class in Python to wrap Solver Foundation. The concept is that the model class wraps a Solver Foundation Model object and corresponding SolverContext. I provide simple wrappers for decisions, constraints, and goals. These wrappers combine the steps of creating and adding entities to a model, and handle type conversion issues. As I have written about in the past, Solver Foundation Services relies heavily on implicit type conversion and operator overloading to provide a productive programming environment for .Net programmers. Not all of these features carry over to Iron Python, so the module presented in this post is meant to hide those limitations. The complete code for my model module is at the end of the post.
The point is that with this module it is a bit easier to build models in Sho using Solver Foundation. Let’s use it to write a function that solves the N-Queens problem. We want to find a placement of queens on an N x N chessboard so that no queen can capture any other in a single move. Because of the rules of chess, we know that we will have to place one and only one queen on each rank (row). Therefore we have one decision per rank: on which file (column) should we place a queen? We have three additional constraints:
- No two queens can be placed on the same file.
- No two queens can be placed on the same \ diagonal.
- No two queens can be placed on the same / diagonal.
Writing a Solver Foundation model for this problem is very easy using the model class. Note that you will need to paste my model code first before trying this sample. I hope I didn’t goof it up:
def queens(count):
m = model()
q = [m.int(0, count - 1) for i in range(0, count)]
m.alldifferent(q)
m.alldifferent([q[i] + const(i) for i in range(0, count)])
m.alldifferent([q[i] - const(i) for i in range(0, count)])
m.solve()
b = zeros(count, count)
for i in range(0, count): b[i, q[i].GetDouble()] = 1
return b
The function creates the model, solves it, and then returns a matrix with the placement of queens. To get the results we use the Decision.GetDouble() method.
Now we can use the built-in visualization functions in Sho to see the output:
imageview(queens(50))
Here is the code for my simple model wrapper. If you have difficulties understanding this code, try consulting:
- the Python language reference,
- the Iron Python documentation,
- the “Book of Sho” included with the Installation of Sho.
from sho import *
from System import Array
ShoLoadAssembly("C:\Program Files (x86)\Sho 2.0 for .NET 4\packages\Optimizer\Microsoft.Solver.Foundation.dll")
from Microsoft.SolverFoundation.Services import *
from Microsoft.SolverFoundation.Common import *
def const(i): return Term.op_Implicit(i)
def term(i):
if isinstance(i, Term):
return i
else:
return const(i)
class model:
"""Represents an optimization model consisting of decision, constraint,
and goal objects. Use the real or int functions to add new decisions to the model.
Use lessequal, greaterequal, equal, alldifferent to add constraints. Use minimize
or maximize to add goals. Use solve to solve the model. Consult the Solver Foundation
Services documentation on MSDN for more.
"""
def __init__(self):
self.c = SolverContext.GetContext()
self.c.ClearModel()
self.m = self.c.CreateModel()
self.i = 0
def _nextname(self, prefix):
name = prefix + ("%s" % self.i)
self.i = self.i + 1
return name
"""Adds a real decision variable."""
def real(self, lower=Rational.NegativeInfinity, upper=Rational.PositiveInfinity):
d = Decision(Domain.RealRange(lower, upper), self._nextname("d"))
self.m.AddDecision(d)
return d
"""Adds an integer decision variable."""
def int(self, lower=Rational.NegativeInfinity, upper=Rational.PositiveInfinity):
d = Decision(Domain.IntegerRange(lower, upper), self._nextname("d"))
self.m.AddDecision(d)
return d
"""Adds a constraint. The input must be a Solver Foundation Term. Use the Model
object to create Terms, e.g. Model.Sum. """
def constraint(self, expr):
return self.m.AddConstraint(self._nextname("c"), expr)
"""Adds an 'all different' constraint from a list of terms."""
def alldifferent(self, exprs):
return self.constraint(Model.AllDifferent(Array[Term](exprs)))
"""Adds an equality constraint."""
def equal(self, expr1, expr2):
return self.constraint(Model.Equal(term(expr1), term(expr2)))
"""Adds an inequality constraint."""
def notequal(self, expr1, expr2):
return self.constraint(Model.Not(Model.Equal(term(expr1), term(expr2))))
"""Adds an inequality constraint."""
def greaterequal(self, expr1, expr2):
return self.constraint(Model.GreaterEqual(term(expr1), term(expr2)))
"""Adds an inequality constraint."""
def lessequal(self, expr1, expr2):
return self.constraint(Model.LessEqual(term(expr1), term(expr2)))
"""Adds a minimization goal."""
def minimize(self, expr):
return self.m.AddGoal(self._nextname("g"), GoalKind.Minimize, expr)
"""Adds a maximization goal."""
def maximize(self, expr):
return self.m.AddGoal(self._nextname("g"), GoalKind.Maximize, expr)
"""Solves the model, returning a solution object."""
def solve(self, verbose = False):
sln = self.c.Solve()
if verbose: print(sln.GetReport().ToString())
return sln
"""Return the OML representation of the model."""
def oml(self):
sw = System.IO.StringWriter()
self.c.SaveModel(FileFormat.OML, sw)
return sw.ToString()