Commit 36289259 authored by David Bommes's avatar David Bommes

added GUROBISolver

added VariableType
added is_linear function to NConstraints bas class

git-svn-id: http://www.openflipper.org/svnrepo/CoMISo/trunk@101 1355f012-dd97-4b2f-ae87-10fa9f823a57
parent 15572ed0
......@@ -57,6 +57,7 @@ public:
virtual void eval_gradient ( const double* _x, SVectorNC& _g ) { _g.resize(n_); _g.coeffRef(idx_) = 1.0; }
virtual void eval_hessian ( const double* _x, SMatrixNC& _h ) { _h.clear(); _h.resize(n_,n_); }
virtual bool is_linear() { return true;}
// set/get values
unsigned int& idx() {return idx_;}
......
//=============================================================================
//
// CLASS GUROBISolver - IMPLEMENTATION
//
//=============================================================================
//== INCLUDES =================================================================
//== COMPILE-TIME PACKAGE REQUIREMENTS ========================================
#include <CoMISo/Config/config.hh>
#if COMISO_GUROBI_AVAILABLE
//=============================================================================
#include "GUROBISolver.hh"
//== NAMESPACES ===============================================================
namespace COMISO {
//== IMPLEMENTATION ==========================================================
//// Constructor
//GUROBISolver::
//GUROBISolver()
//{
//
//}
//-----------------------------------------------------------------------------
// ********** SOLVE **************** //
bool
GUROBISolver::
solve(NProblemInterface* _problem,
std::vector<NConstraintInterface*>& _constraints,
std::vector<PairUiV>& _discrete_constraints,
const double _time_limit)
{
try
{
//----------------------------------------------
// 0. set up environment
//----------------------------------------------
GRBEnv env = GRBEnv();
GRBModel model = GRBModel(env);
model.getEnv().set(GRB_DoubleParam_TimeLimit, _time_limit);
//----------------------------------------------
// 1. allocate variables
//----------------------------------------------
// determine variable types: 0->real, 1->integer, 2->bool
std::vector<char> vtypes(_problem->n_unknowns(),0);
for(unsigned int i=0; i<_discrete_constraints.size(); ++i)
switch(_discrete_constraints[i].second)
{
case Integer: vtypes[_discrete_constraints[i].first] = 1; break;
case Binary : vtypes[_discrete_constraints[i].first] = 2; break;
default : break;
}
// GUROBI variables
std::vector<GRBVar> vars;
// first all
for( int i=0; i<_problem->n_unknowns(); ++i)
switch(vtypes[i])
{
case 0 : vars.push_back( model.addVar(-GRB_INFINITY, GRB_INFINITY, 0.0, GRB_CONTINUOUS) ); break;
case 1 : vars.push_back( model.addVar(-GRB_INFINITY, GRB_INFINITY, 0.0, GRB_INTEGER ) ); break;
case 2 : vars.push_back( model.addVar(-GRB_INFINITY, GRB_INFINITY, 0.0, GRB_BINARY ) ); break;
}
// Integrate new variables
model.update();
//----------------------------------------------
// 2. setup constraints
//----------------------------------------------
// get zero vector
std::vector<double> x(_problem->n_unknowns(), 0.0);
for(unsigned int i=0; i<_constraints.size(); ++i)
{
if(!_constraints[i]->is_linear())
std::cerr << "Warning: GUROBISolver received a problem with non-linear constraints!!!" << std::endl;
GRBLinExpr lin_expr;
NConstraintInterface::SVectorNC gc;
_constraints[i]->eval_gradient(P(x), gc);
NConstraintInterface::SVectorNC::InnerIterator v_it(gc);
for(; v_it; ++v_it)
// lin_expr += v_it.value()*vars[v_it.index()];
lin_expr = lin_expr + vars[v_it.index()]*v_it.value();
double b = _constraints[i]->eval_constraint(P(x));
switch(_constraints[i]->constraint_type())
{
case NConstraintInterface::NC_EQUAL : model.addConstr(lin_expr + b == 0); break;
case NConstraintInterface::NC_LESS_EQUAL : model.addConstr(lin_expr + b <= 0); break;
case NConstraintInterface::NC_GREATER_EQUAL : model.addConstr(lin_expr + b >= 0); break;
}
}
model.update();
//----------------------------------------------
// 3. setup energy
//----------------------------------------------
if(!_problem->constant_hessian())
std::cerr << "Warning: GUROBISolver received a problem with non-constant hessian!!!" << std::endl;
GRBQuadExpr objective;
// add quadratic part
NProblemInterface::SMatrixNP H;
_problem->eval_hessian(P(x), H);
for( int i=0; i<H.outerSize(); ++i)
for (NProblemInterface::SMatrixNP::InnerIterator it(H,i); it; ++it)
objective += 0.5*it.value()*vars[it.row()]*vars[it.col()];
// add linear part
std::vector<double> g(_problem->n_unknowns());
_problem->eval_gradient(P(x), P(g));
for(unsigned int i=0; i<g.size(); ++i)
objective += g[i]*vars[i];
// add constant part
objective += _problem->eval_f(P(x));
model.set(GRB_IntAttr_ModelSense, 1);
model.setObjective(objective);
//----------------------------------------------
// 4. solve problem
//----------------------------------------------
model.optimize();
//----------------------------------------------
// 5. store result
//----------------------------------------------
for(unsigned int i=0; i<vars.size(); ++i)
x[i] = vars[i].get(GRB_DoubleAttr_X);
_problem->store_result(P(x));
std::cout << "GUROBI Objective: " << model.get(GRB_DoubleAttr_ObjVal) << std::endl;
return true;
}
catch(GRBException e)
{
std::cout << "Error code = " << e.getErrorCode() << std::endl;
std::cout << e.getMessage() << std::endl;
return false;
}
catch(...)
{
std::cout << "Exception during optimization" << std::endl;
return false;
}
return false;
}
//=============================================================================
} // namespace COMISO
//=============================================================================
#endif // COMISO_GUROBI_AVAILABLE
//=============================================================================
//=============================================================================
//
// CLASS GUROBISolver
//
//=============================================================================
#ifndef COMISO_GUROBISOLVER_HH
#define COMISO_GUROBISOLVER_HH
//== COMPILE-TIME PACKAGE REQUIREMENTS ========================================
#include <CoMISo/Config/config.hh>
#if COMISO_GUROBI_AVAILABLE
//== INCLUDES =================================================================
#include <CoMISo/Config/CoMISoDefines.hh>
#include <vector>
#include "NProblemInterface.hh"
#include "NConstraintInterface.hh"
#include "VariableType.hh"
#include <gurobi_c++.h>
//== FORWARDDECLARATIONS ======================================================
//== NAMESPACES ===============================================================
namespace COMISO {
//== CLASS DEFINITION =========================================================
/** \class NewtonSolver GUROBISolver.hh
Brief Description.
A more elaborate description follows.
*/
class COMISODLLEXPORT GUROBISolver
{
public:
typedef std::pair<unsigned int, VariableType> PairUiV;
/// Default constructor -> set up IpOptApplication
GUROBISolver() {}
/// Destructor
~GUROBISolver() {}
// ********** SOLVE **************** //
bool solve(NProblemInterface* _problem, // problem instance
std::vector<NConstraintInterface*>& _constraints, // linear constraints
std::vector<PairUiV>& _discrete_constraints, // discrete constraints
const double _time_limit = 60 ); // time limit in seconds
protected:
double* P(std::vector<double>& _v)
{
if( !_v.empty())
return ((double*)&_v[0]);
else
return 0;
}
private:
};
//=============================================================================
} // namespace COMISO
//=============================================================================
#endif // COMISO_GUROBI_AVAILABLE
//=============================================================================
#endif // ACG_GUROBISOLVER_HH defined
//=============================================================================
......@@ -65,6 +65,8 @@ public:
virtual void eval_hessian ( const double* _x, SMatrixNC& _h );
virtual bool is_linear() { return true;}
// inherited from base
// virtual ConstraintType constraint_type ( ) { return type_; }
......
......@@ -67,6 +67,8 @@ public:
return false;
}
virtual bool is_linear() { return false;}
virtual double gradient_update_factor( const double* _x, double _eps = 1e-6 )
{
double val = eval_constraint(_x);
......
//=============================================================================
//
// ENUM VariableType
//
//=============================================================================
#ifndef COMISO_VARIABLETYPE_HH
#define COMISO_VARIABLETYPE_HH
//== INCLUDES =================================================================
#include <CoMISo/Config/CoMISoDefines.hh>
//== FORWARDDECLARATIONS ======================================================
//== NAMESPACES ===============================================================
namespace COMISO {
//== CLASS DEFINITION =========================================================
enum VariableType { Real, Integer, Binary};
//=============================================================================
} // namespace COMISO
//=============================================================================
#endif // COMISO_VARIABLETYPE_HH defined
//=============================================================================
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment