add API to query activity of a row

This commit is contained in:
lperron@google.com
2011-08-11 18:44:23 +00:00
parent d0926c5188
commit 3794f40fba
9 changed files with 133 additions and 32 deletions

View File

@@ -61,14 +61,17 @@ public class LinearSolverExample {
System.out.println("Optimal value = " + solver.objectiveValue());
System.out.println("x1 = " + x1.solutionValue() +
" + reduced cost = " + x1.reducedCost());
" reduced cost = " + x1.reducedCost());
System.out.println("x2 = " + x2.solutionValue() +
" + reduced cost = " + x2.reducedCost());
" reduced cost = " + x2.reducedCost());
System.out.println("x3 = " + x3.solutionValue() +
" + reduced cost = " + x3.reducedCost());
System.out.println("c0 dual value = " + c0.dualValue());
System.out.println("c1 dual value = " + c1.dualValue());
System.out.println("c2 dual value = " + c2.dualValue());
" reduced cost = " + x3.reducedCost());
System.out.println("c0 dual value = " + c0.dualValue() +
" activity = " + c0.activity());
System.out.println("c1 dual value = " + c1.dualValue() +
" activity = " + c1.activity());
System.out.println("c2 dual value = " + c2.dualValue() +
" activity = " + c2.activity());
}
public static void main(String[] args) throws Exception {

View File

@@ -51,9 +51,12 @@ void BuildLinearProgrammingMaxExample(MPSolver::OptimizationProblemType type) {
<< ", reduced cost = " << x2->reduced_cost();
LOG(INFO) << "x3 = " << x3->solution_value()
<< ", reduced cost = " << x3->reduced_cost();
LOG(INFO) << "c0 dual value = " << c0->dual_value();
LOG(INFO) << "c1 dual value = " << c1->dual_value();
LOG(INFO) << "c2 dual value = " << c2->dual_value();
LOG(INFO) << "c0 dual value = " << c0->dual_value()
<< " activity = " << c0->activity();
LOG(INFO) << "c1 dual value = " << c1->dual_value()
<< " activity = " << c1->activity();
LOG(INFO) << "c2 dual value = " << c2->dual_value()
<< " activity = " << c2->activity();
}
void RunAllExamples() {

View File

@@ -142,6 +142,14 @@ class CBCInterface : public MPSolverInterface {
virtual string SolverVersion() const {
return PACKAGE_STRING;
}
// TODO(user): Maybe we should expose the CbcModel build from osi_
// instead, but a new CbcModel is built every time Solve is called,
// so it is not possible right now.
virtual void* underlying_solver() {
return reinterpret_cast<void*>(&osi_);
}
private:
// Reset best objective bound to +/- infinity depending on the
// optimization direction.
@@ -400,7 +408,6 @@ MPSolver::ResultStatus CBCInterface::Solve(const MPSolverParameters& param) {
objective_value_ = model.getObjValue();
VLOG(1) << "objective=" << objective_value_;
const double* const values = model.bestSolution();
if (values != NULL) {
// if optimal or feasible solution is found.
for (int i = 0; i < solver_->variables_.size(); ++i) {
@@ -414,6 +421,18 @@ MPSolver::ResultStatus CBCInterface::Solve(const MPSolverParameters& param) {
VLOG(1) << "No feasible solution found.";
}
const double* const row_activities = model.getRowActivity();
if (row_activities != NULL) {
for (int i = 0; i < solver_->constraints_.size(); ++i) {
MPConstraint* const ct = solver_->constraints_[i];
const int constraint_index = ct->index();
const double row_activity = row_activities[constraint_index];
ct->set_activity(row_activity);
VLOG(4) << "row " << ct->index()
<< ": activity = " << row_activity;
}
}
// Check the status: optimal, infeasible, etc.
int tmp_status = model.status();

View File

@@ -114,6 +114,10 @@ class CLPInterface : public MPSolverInterface {
return PACKAGE_STRING;
}
virtual void* underlying_solver() {
return reinterpret_cast<void*>(clp_.get());
}
private:
// Create dummy variable to be able to create empty constraints.
void CreateDummyVariableForEmptyConstraints();
@@ -453,7 +457,7 @@ MPSolver::ResultStatus CLPInterface::Solve(const MPSolverParameters& param) {
const double* const reduced_costs = clp_->getReducedCost();
for (int i = 0; i < solver_->variables_.size(); ++i) {
MPVariable* const var = solver_->variables_[i];
int var_index = var->index();
const int var_index = var->index();
double val = values[var_index];
var->set_solution_value(val);
VLOG(3) << var->name() << ": value = " << val;
@@ -462,12 +466,17 @@ MPSolver::ResultStatus CLPInterface::Solve(const MPSolverParameters& param) {
VLOG(4) << var->name() << ": reduced cost = " << reduced_cost;
}
const double* const dual_values = clp_->getRowPrice();
const double* const row_activities = clp_->getRowActivity();
for (int i = 0; i < solver_->constraints_.size(); ++i) {
MPConstraint* const ct = solver_->constraints_[i];
double dual_value = 0.0;
dual_value = dual_values[i];
const int constraint_index = ct->index();
const double row_activity = row_activities[constraint_index];
ct->set_activity(row_activity);
const double dual_value = dual_values[constraint_index];
ct->set_dual_value(dual_value);
VLOG(4) << "row " << ct->index() << ": dual value = " << dual_value;
VLOG(4) << "row " << ct->index()
<< ": activity = " << row_activity
<< " dual value = " << dual_value;
}
// Check the status: optimal, infeasible, etc.

View File

@@ -160,6 +160,10 @@ class GLPKInterface : public MPSolverInterface {
return StringPrintf("GLPK %s", glp_version());
}
virtual void* underlying_solver() {
return reinterpret_cast<void*>(lp_);
}
private:
// Configure the solver's parameters.
void ConfigureGLPKParameters(const MPSolverParameters& param);
@@ -573,12 +577,21 @@ MPSolver::ResultStatus GLPKInterface::Solve(const MPSolverParameters& param) {
VLOG(4) << var->name() << ": reduced cost = " << reduced_cost;
}
}
if (!mip_) {
for (int i = 0; i < solver_->constraints_.size(); ++i) {
MPConstraint* const ct = solver_->constraints_[i];
double dual_value = glp_get_row_dual(lp_, ct->index());
for (int i = 0; i < solver_->constraints_.size(); ++i) {
MPConstraint* const ct = solver_->constraints_[i];
if (mip_) {
const double row_activity = glp_mip_row_val(lp_, ct->index());
ct->set_activity(row_activity);
VLOG(4) << "row " << ct->index()
<< ": activity = " << row_activity;
} else {
const double row_activity = glp_get_row_prim(lp_, ct->index());
ct->set_activity(row_activity);
const double dual_value = glp_get_row_dual(lp_, ct->index());
ct->set_dual_value(dual_value);
VLOG(4) << "row " << ct->index() << ": dual value = " << dual_value;
VLOG(4) << "row " << ct->index()
<< ": activity = " << row_activity
<< ": dual value = " << dual_value;
}
}

View File

@@ -138,6 +138,12 @@ double MPConstraint::dual_value() const {
return dual_value_;
}
double MPConstraint::activity() const {
interface_->CheckSolutionIsSynchronized();
interface_->CheckSolutionExists();
return activity_;
}
bool MPConstraint::ContainsNewVariables() {
const int last_variable_index = interface_->last_variable_index();
for (ConstIter<hash_map<MPVariable*, double> > it(coefficients_);
@@ -280,6 +286,12 @@ string MPSolver::SolverVersion() const {
return interface_->SolverVersion();
}
// ---- Underlying solver ----
void* MPSolver::underlying_solver() {
return interface_->underlying_solver();
}
// ----- Solver -----
#if defined(USE_CLP) || defined(USE_CBC)

View File

@@ -165,6 +165,9 @@ class MPConstraint {
void SetUB(double ub) { SetBounds(lb_, ub); }
void SetBounds(double lb, double ub);
// Returns the constraint's activity in the current solution:
// sum over all terms of (coefficient * variable value)
double activity() const;
// Only available for continuous problems.
double dual_value() const;
@@ -183,10 +186,11 @@ class MPConstraint {
const string& name,
MPSolverInterface* const interface)
: lb_(lb), ub_(ub), name_(name), index_(-1), dual_value_(0.0),
interface_(interface) {}
activity_(0.0), interface_(interface) {}
void set_index(int index) { index_ = index; }
void set_dual_value(double dual_value) { dual_value_ = dual_value; }
void set_activity(double activity) { activity_ = activity; }
private:
// Returns true if the constraint contains variables that have not
// been extracted yet.
@@ -203,6 +207,7 @@ class MPConstraint {
const string name_;
int index_;
double dual_value_;
double activity_;
MPSolverInterface* const interface_;
DISALLOW_COPY_AND_ASSIGN(MPConstraint);
};
@@ -447,6 +452,19 @@ class MPSolver {
// return a string describing the engine used.
string SolverVersion() const;
// Returns the underlying solver so that the user can use
// solver-specific features or features that are not exposed in the
// simple API of MPSolver. This method is for advanced users, use at
// your own risk! In particular, if you modify the model or the
// solution by accessing the underlying solver directly, then the
// underlying solver will be out of sync with the information kept
// in the wrapper (MPSolver, MPVariable, MPConstraint,
// MPObjective). You need to cast the void* returned back to its
// original type that depends on the interface (CBC:
// OsiClpSolverInterface*, CLP: ClpSimplex*, GLPK: glp_prob*, SCIP:
// SCIP*).
void* underlying_solver();
friend class GLPKInterface;
friend class CLPInterface;
friend class CBCInterface;
@@ -707,6 +725,9 @@ class MPSolverInterface {
// Returns a string describing the solver.
virtual string SolverVersion() const = 0;
// Returns the underlying solver.
virtual void* underlying_solver() = 0;
friend class MPSolver;
protected:
MPSolver* const solver_;

View File

@@ -114,6 +114,10 @@ class SCIPInterface : public MPSolverInterface {
SCIPtechVersion(), SCIPlpiGetSolverName());
}
virtual void* underlying_solver() {
return reinterpret_cast<void*>(scip_);
}
private:
// Set all parameters in the underlying solver.
virtual void SetParameters(const MPSolverParameters& param);
@@ -484,6 +488,17 @@ MPSolver::ResultStatus SCIPInterface::Solve(const MPSolverParameters& param) {
var->set_solution_value(val);
VLOG(3) << var->name() << "=" << val;
}
for (int i = 0; i < solver_->constraints_.size(); ++i) {
MPConstraint* const ct = solver_->constraints_[i];
const int constraint_index = ct->index();
const double row_activity =
SCIPgetActivityLinear(scip_,
scip_constraints_[constraint_index],
solution);
ct->set_activity(row_activity);
VLOG(4) << "row " << ct->index()
<< ": activity = " << row_activity;
}
} else {
VLOG(1) << "No feasible solution found.";
}

View File

@@ -33,7 +33,7 @@ class PyWrapLPExamples(object):
x1 = solver.NumVar(0.0, infinity, 'x1')
x2 = solver.NumVar(0.0, infinity, 'x2')
x3 = solver.NumVar(0.0, infinity, 'x3')
print solver.NumVariables()
print 'number of variables = ', solver.NumVariables()
solver.AddObjectiveTerm(x1, 10)
solver.AddObjectiveTerm(x2, 6)
@@ -52,18 +52,21 @@ class PyWrapLPExamples(object):
c2.AddTerm(x1, 2)
c2.AddTerm(x2, 2)
c2.AddTerm(x3, 6)
print solver.NumConstraints()
print 'number of constraints = ', solver.NumConstraints()
# The problem has an optimal solution.
solver.Solve()
print solver.objective_value()
print 'optimal objective value = ', solver.objective_value()
print 'x1 = ', x1.solution_value(), ', reduced_cost = ', x1.reduced_cost()
print 'x2 = ', x2.solution_value(), ', reduced_cost = ', x2.reduced_cost()
print 'x3 = ', x3.solution_value(), ', reduced_cost = ', x3.reduced_cost()
print 'c0 dual value = ', c0.dual_value()
print 'c0 activity = ', c0.activity()
print 'c1 dual value = ', c1.dual_value()
print 'c1 activity = ', c1.activity()
print 'c2 dual value = ', c2.dual_value()
print 'c2 activity = ', c2.activity()
def RunAllFirstLinearExample(self):
self.RunFirstLinearExample(pywraplp.Solver.GLPK_LINEAR_PROGRAMMING)
@@ -76,7 +79,7 @@ class PyWrapLPExamples(object):
x1 = solver.NumVar(0.0, infinity, 'x1')
x2 = solver.NumVar(0.0, infinity, 'x2')
x3 = solver.NumVar(0.0, infinity, 'x3')
print solver.NumVariables()
print 'number of variables = ', solver.NumVariables()
solver.Maximize(10 * x1 + 6 * x2 + 4 * x3)
@@ -84,18 +87,21 @@ class PyWrapLPExamples(object):
c1 = solver.Add(10 * x1 + 4 * x2 + 5 * x3 <= 600)
c2 = solver.Add(2 * x1 + 2 * x2 + 6 * x3 <= 300)
print solver.NumConstraints()
print 'number of constraints = ', solver.NumConstraints()
# The problem has an optimal solution.
solver.Solve()
print solver.objective_value()
print 'optimal objective value = ', solver.objective_value()
print 'x1 = ', x1.solution_value(), ', reduced_cost = ', x1.reduced_cost()
print 'x2 = ', x2.solution_value(), ', reduced_cost = ', x2.reduced_cost()
print 'x3 = ', x3.solution_value(), ', reduced_cost = ', x3.reduced_cost()
print 'c0 dual value = ', c0.dual_value()
print 'c0 activity = ', c0.activity()
print 'c1 dual value = ', c1.dual_value()
print 'c1 activity = ', c1.activity()
print 'c2 dual value = ', c2.dual_value()
print 'c2 activity = ', c2.activity()
def RunAllFirstLinearExampleNewAPI(self):
self.RunFirstLinearExampleNewAPI(pywraplp.Solver.GLPK_LINEAR_PROGRAMMING)
@@ -111,22 +117,22 @@ class PyWrapLPExamples(object):
# Check the solution
solver.Solve()
print x1.solution_value()
print x2.solution_value()
print 'x1 = ', x1.solution_value()
print 'x2 = ', x2.solution_value()
solver.Maximize(x2)
# Check the solution
solver.Solve()
print x1.solution_value()
print x2.solution_value()
print 'x1 = ', x1.solution_value()
print 'x2 = ', x2.solution_value()
solver.Minimize(-x1)
# Check the solution
solver.Solve()
print x1.solution_value()
print x2.solution_value()
print 'x1 = ', x1.solution_value()
print 'x2 = ', x2.solution_value()
def RunAllSuccessiveObjectives(self):
self.RunSuccessiveObjectives(pywraplp.Solver.GLPK_LINEAR_PROGRAMMING)