changes from Christophe-Marie Duquesne supporting PYGLPK and YAPOSIB

This commit is contained in:
Stuart Mitchell
2013-01-07 22:01:15 +13:00
3 changed files with 111 additions and 71 deletions

View File

@@ -1,2 +1,3 @@
Roy, J.S Roy, J.S
Mitchell, Stuart A Mitchell, Stuart A
Christophe-Marie Duquesne (PYGLPK and YAPOSIB bindings)

View File

@@ -956,7 +956,7 @@ else:
""" """
Solve a well formulated lp problem Solve a well formulated lp problem
creates a gurobi model, variables and constraints and attaches creates a cplex model, variables and constraints and attaches
them to the lp model which it then solves them to the lp model which it then solves
""" """
self.buildSolverModel(lp) self.buildSolverModel(lp)
@@ -1900,6 +1900,8 @@ class PYGLPK(LpSolver):
""" """
The glpk LP/MIP solver (via its python interface) The glpk LP/MIP solver (via its python interface)
Copyright Christophe-Marie Duquesne 2012
The glpk variables are available (after a solve) in var.solverVar The glpk variables are available (after a solve) in var.solverVar
The glpk constraints are available in constraint.solverConstraint The glpk constraints are available in constraint.solverConstraint
The Model is in prob.solverModel The Model is in prob.solverModel
@@ -1907,7 +1909,7 @@ class PYGLPK(LpSolver):
try: try:
#import the model into the global scope #import the model into the global scope
global glpk global glpk
import glpk import glpk.glpkpi as glpk
except ImportError: except ImportError:
def available(self): def available(self):
"""True if the solver is available""" """True if the solver is available"""
@@ -1927,42 +1929,42 @@ class PYGLPK(LpSolver):
@param mip: if False the solver will solve a MIP as an LP @param mip: if False the solver will solve a MIP as an LP
@param msg: displays information from the solver to stdout @param msg: displays information from the solver to stdout
@param timeLimit: not handled by glpk @param timeLimit: not handled
@param epgap: sets the integer bound gap @param epgap: not handled
@param solverParams: not handled @param solverParams: not handled
""" """
LpSolver.__init__(self, mip, msg) LpSolver.__init__(self, mip, msg)
self.timeLimit = timeLimit # time limits are not handled
self.epgap = epgap
if not self.msg: if not self.msg:
glpk.env.term_on = False glpk.glp_term_out(glpk.GLP_OFF)
def findSolutionValues(self, lp): def findSolutionValues(self, lp):
model = lp.solverModel prob = lp.solverModel
solutionStatus = model.status if self.mip and self.hasMIPConstraints(lp.solverModel):
glpkLpStatus = {"opt": LpStatusOptimal, solutionStatus = glpk.glp_mip_status(prob)
"undef": LpStatusUndefined, else:
"feas": LpStatusNotSolved, solutionStatus = glpk.glp_get_status(prob)
"infeas": LpStatusInfeasible, glpkLpStatus = {glpk.GLP_OPT: LpStatusOptimal,
"nofeas": LpStatusInfeasible, glpk.GLP_UNDEF: LpStatusUndefined,
"unbnd": LpStatusUnbounded glpk.GLP_FEAS: LpStatusNotSolved,
glpk.GLP_INFEAS: LpStatusInfeasible,
glpk.GLP_NOFEAS: LpStatusInfeasible,
glpk.GLP_UNBND: LpStatusUnbounded
} }
#populate pulp solution values #populate pulp solution values
for var in lp.variables(): for var in lp.variables():
var.varValue = var.solverVar.primal if self.mip and self.hasMIPConstraints(lp.solverModel):
try: var.varValue = glpk.glp_mip_col_val(prob, var.glpk_index)
var.dj = var.solverVar.dual else:
except RuntimeError: var.varValue = glpk.glp_get_col_prim(prob, var.glpk_index)
var.dj = None var.dj = glpk.glp_get_col_dual(prob, var.glpk_index)
#put pi and slack variables against the constraints #put pi and slack variables against the constraints
for constr in lp.constraints.values(): for constr in lp.constraints.values():
try: if self.mip and self.hasMIPConstraints(lp.solverModel):
constr.pi = constr.solverConstraint.dual row_val = glpk.glp_mip_row_val(prob, constr.glpk_index)
except RuntimeError: else:
constr.pi = None row_val = glpk.glp_get_row_prim(prob, constr.glpk_index)
constr.slack = constr.solverConstraint.primal constr.slack = -constr.constant - row_val
if self.msg: constr.pi = glpk.glp_get_row_dual(prob, constr.glpk_index)
print "glpk status=", solutionStatus
lp.resolveOK = True lp.resolveOK = True
for var in lp.variables(): for var in lp.variables():
var.isModified = False var.isModified = False
@@ -1974,17 +1976,20 @@ class PYGLPK(LpSolver):
"""True if the solver is available""" """True if the solver is available"""
return True return True
def hasMIPConstraints(self, solverModel):
return (glpk.glp_get_num_int(solverModel) > 0 or
glpk.glp_get_num_bin(solverModel) > 0)
def callSolver(self, lp, callback = None): def callSolver(self, lp, callback = None):
"""Solves the problem with glpk """Solves the problem with glpk
""" """
self.solveTime = -clock() self.solveTime = -clock()
lp.solverModel.simplex() glpk.glp_adv_basis(lp.solverModel, 0)
if self.mip: glpk.glp_simplex(lp.solverModel, None)
if (lp.solverModel.status != "infeas" if self.mip and self.hasMIPConstraints(lp.solverModel):
and lp.solverModel.status != "nofeas" status = glpk.glp_get_status(lp.solverModel)
and lp.solverModel.status != "unbnd" if status in (glpk.GLP_OPT, glpk.GLP_UNDEF, glpk.GLP_FEAS):
): glpk.glp_intopt(lp.solverModel, None)
lp.solverModel.integer()
self.solveTime += clock() self.solveTime += clock()
def buildSolverModel(self, lp): def buildSolverModel(self, lp):
@@ -1992,45 +1997,69 @@ class PYGLPK(LpSolver):
Takes the pulp lp model and translates it into a glpk model Takes the pulp lp model and translates it into a glpk model
""" """
log.debug("create the glpk model") log.debug("create the glpk model")
lp.solverModel = glpk.LPX() prob = glpk.glp_create_prob()
lp.solverModel.name = lp.name glpk.glp_set_prob_name(prob, lp.name)
log.debug("set the sense of the problem") log.debug("set the sense of the problem")
if lp.sense == LpMaximize: if lp.sense == LpMaximize:
lp.solverModel.obj.maximize = True glpk.glp_set_obj_dir(prob, glpk.GLP_MAX)
log.debug("add the Constraints to the problem") log.debug("add the constraints to the problem")
lp.solverModel.rows.add(len(lp.constraints.keys())) glpk.glp_add_rows(prob, len(lp.constraints.keys()))
i = 0 for i, v in enumerate(lp.constraints.items(), start=1):
for name, constraint in lp.constraints.items(): name, constraint = v
row = lp.solverModel.rows[i] glpk.glp_set_row_name(prob, i, name)
row.name = name
if constraint.sense == LpConstraintLE: if constraint.sense == LpConstraintLE:
row.bounds = None,-constraint.constant glpk.glp_set_row_bnds(prob, i, glpk.GLP_UP,
0.0, -constraint.constant)
elif constraint.sense == LpConstraintGE: elif constraint.sense == LpConstraintGE:
row.bounds = -constraint.constant, None glpk.glp_set_row_bnds(prob, i, glpk.GLP_LO,
-constraint.constant, 0.0)
elif constraint.sense == LpConstraintEQ: elif constraint.sense == LpConstraintEQ:
row.bounds = -constraint.constant,-constraint.constant glpk.glp_set_row_bnds(prob, i, glpk.GLP_FX,
-constraint.constant, -constraint.constant)
else: else:
raise PulpSolverError, 'Detected an invalid constraint type' raise PulpSolverError, 'Detected an invalid constraint type'
i += 1 constraint.glpk_index = i
constraint.solverConstraint = row
log.debug("add the variables to the problem") log.debug("add the variables to the problem")
lp.solverModel.cols.add(len(lp.variables())) glpk.glp_add_cols(prob, len(lp.variables()))
j = 0 for j, var in enumerate(lp.variables(), start=1):
for var in lp.variables(): glpk.glp_set_col_name(prob, j, var.name)
col = lp.solverModel.cols[j] lb = 0.0
col.name = var.name ub = 0.0
col.bounds = var.lowBound,var.upBound t = glpk.GLP_FR
if not var.lowBound is None:
lb = var.lowBound
t = glpk.GLP_LO
if not var.upBound is None:
ub = var.upBound
t = glpk.GLP_UP
if not var.upBound is None and not var.lowBound is None:
if ub == lb:
t = glpk.GLP_FX
else:
t = glpk.GLP_DB
glpk.glp_set_col_bnds(prob, j, t, lb, ub)
if var.cat == LpInteger: if var.cat == LpInteger:
col.kind = int glpk.glp_set_col_kind(prob, j, glpk.GLP_IV)
var.solverVar = col assert glpk.glp_get_col_kind(prob, j) == glpk.GLP_IV
j += 1 var.glpk_index = j
log.debug("set the objective function") log.debug("set the objective function")
lp.solverModel.obj[:] = [lp.objective.get(var, 0.0) for var in for var in lp.variables():
lp.variables()] value = lp.objective.get(var)
if value:
glpk.glp_set_obj_coef(prob, var.glpk_index, value)
log.debug("set the problem matrix") log.debug("set the problem matrix")
for name,constraint in lp.constraints.items(): for constraint in lp.constraints.values():
constraint.solverConstraint.matrix =[(var.solverVar.index, l = len(constraint.items())
value ) for var, value in constraint.items()] ind = glpk.intArray(l + 1)
val = glpk.doubleArray(l + 1)
for j, v in enumerate(constraint.items(), start=1):
var, value = v
ind[j] = var.glpk_index
val[j] = value
glpk.glp_set_mat_row(prob, constraint.glpk_index, l, ind,
val)
lp.solverModel = prob
#glpk.glp_write_lp(prob, None, "glpk.lp")
def actualSolve(self, lp, callback = None): def actualSolve(self, lp, callback = None):
""" """
@@ -2060,14 +2089,17 @@ class PYGLPK(LpSolver):
""" """
log.debug("Resolve the Model using glpk") log.debug("Resolve the Model using glpk")
for constraint in lp.constraints.values(): for constraint in lp.constraints.values():
row = constraint.solverConstraint i = constraint.glpk_index
if constraint.modified: if constraint.modified:
if constraint.sense == LpConstraintLE: if constraint.sense == LpConstraintLE:
row.bounds = None,-constraint.constant glpk.glp_set_row_bnds(prob, i, glpk.GLP_UP,
0.0, -constraint.constant)
elif constraint.sense == LpConstraintGE: elif constraint.sense == LpConstraintGE:
row.bounds = -constraint.constant, None glpk.glp_set_row_bnds(prob, i, glpk.GLP_LO,
-constraint.constant, 0.0)
elif constraint.sense == LpConstraintEQ: elif constraint.sense == LpConstraintEQ:
row.bounds = -constraint.constant,-constraint.constant glpk.glp_set_row_bnds(prob, i, glpk.GLP_FX,
-constraint.constant, -constraint.constant)
else: else:
raise PulpSolverError, 'Detected an invalid constraint type' raise PulpSolverError, 'Detected an invalid constraint type'
self.callSolver(lp, callback = callback) self.callSolver(lp, callback = callback)
@@ -2084,6 +2116,8 @@ class YAPOSIB(LpSolver):
""" """
COIN OSI (via its python interface) COIN OSI (via its python interface)
Copyright Christophe-Marie Duquesne 2012
The yaposib variables are available (after a solve) in var.solverVar The yaposib variables are available (after a solve) in var.solverVar
The yaposib constraints are available in constraint.solverConstraint The yaposib constraints are available in constraint.solverConstraint
The Model is in prob.solverModel The Model is in prob.solverModel
@@ -2105,7 +2139,7 @@ class YAPOSIB(LpSolver):
msg = True, msg = True,
timeLimit = None, timeLimit = None,
epgap = None, epgap = None,
solverName = "Clp", solverName = None,
**solverParams): **solverParams):
""" """
Initializes the yaposib solver. Initializes the yaposib solver.
@@ -2119,7 +2153,10 @@ class YAPOSIB(LpSolver):
@param solverParams: not supported @param solverParams: not supported
""" """
LpSolver.__init__(self, mip, msg) LpSolver.__init__(self, mip, msg)
self.solverName = solverName if solverName:
self.solverName = solverName
else:
self.solverName = yaposib.available_solvers()[0]
def findSolutionValues(self, lp): def findSolutionValues(self, lp):
model = lp.solverModel model = lp.solverModel
@@ -2137,7 +2174,7 @@ class YAPOSIB(LpSolver):
#put pi and slack variables against the constraints #put pi and slack variables against the constraints
for constr in lp.constraints.values(): for constr in lp.constraints.values():
constr.pi = constr.solverConstraint.dual constr.pi = constr.solverConstraint.dual
constr.slack = constr.solverConstraint.activity constr.slack = -constr.constant - constr.solverConstraint.activity
if self.msg: if self.msg:
print "yaposib status=", solutionStatus print "yaposib status=", solutionStatus
lp.resolveOK = True lp.resolveOK = True

View File

@@ -357,7 +357,8 @@ def pulpTest075(solver):
x = LpVariable("x", 0, 4, LpContinuous, obj + b) x = LpVariable("x", 0, 4, LpContinuous, obj + b)
y = LpVariable("y", -1, 1, LpContinuous, 4*obj - c) y = LpVariable("y", -1, 1, LpContinuous, 4*obj - c)
z = LpVariable("z", 0, None, LpContinuous, 9*obj + b + c) z = LpVariable("z", 0, None, LpContinuous, 9*obj + b + c)
if solver.__class__ in [CPLEX_DLL, CPLEX_CMD, COINMP_DLL]: if solver.__class__ in [CPLEX_DLL, CPLEX_CMD, COINMP_DLL, YAPOSIB,
PYGLPK]:
print "\t Testing column based modelling with empty constraints" print "\t Testing column based modelling with empty constraints"
pulpTestCheck(prob, solver, [LpStatusOptimal], {x:4, y:-1, z:6}) pulpTestCheck(prob, solver, [LpStatusOptimal], {x:4, y:-1, z:6})
@@ -378,7 +379,8 @@ def pulpTest080(solver):
prob += c2,"c2" prob += c2,"c2"
prob += c3,"c3" prob += c3,"c3"
if solver.__class__ in [CPLEX_DLL, CPLEX_CMD, COINMP_DLL, PULP_CBC_CMD]: if solver.__class__ in [CPLEX_DLL, CPLEX_CMD, COINMP_DLL,
PULP_CBC_CMD, YAPOSIB, PYGLPK]:
print "\t Testing dual variables and slacks reporting" print "\t Testing dual variables and slacks reporting"
pulpTestCheck(prob, solver, [LpStatusOptimal], pulpTestCheck(prob, solver, [LpStatusOptimal],
sol = {x:4, y:-1, z:6}, sol = {x:4, y:-1, z:6},