1.1 --- a/lemon/lp_glpk.cc Tue Dec 02 21:40:33 2008 +0100
1.2 +++ b/lemon/lp_glpk.cc Tue Dec 02 22:48:28 2008 +0100
1.3 @@ -17,628 +17,936 @@
1.4 */
1.5
1.6 ///\file
1.7 -///\brief Implementation of the LEMON-GLPK lp solver interface.
1.8 +///\brief Implementation of the LEMON GLPK LP and MIP solver interface.
1.9
1.10 #include <lemon/lp_glpk.h>
1.11 -//#include <iostream>
1.12 +#include <glpk.h>
1.13
1.14 -extern "C" {
1.15 -#include <glpk.h>
1.16 -}
1.17 -
1.18 -#if GLP_MAJOR_VERSION > 4 || (GLP_MAJOR_VERSION == 4 && GLP_MINOR_VERSION > 15)
1.19 -#define LEMON_glp(func) (glp_##func)
1.20 -#define LEMON_lpx(func) (lpx_##func)
1.21 -
1.22 -#define LEMON_GLP(def) (GLP_##def)
1.23 -#define LEMON_LPX(def) (LPX_##def)
1.24 -
1.25 -#else
1.26 -
1.27 -#define LEMON_glp(func) (lpx_##func)
1.28 -#define LEMON_lpx(func) (lpx_##func)
1.29 -
1.30 -#define LEMON_GLP(def) (LPX_##def)
1.31 -#define LEMON_LPX(def) (LPX_##def)
1.32 -
1.33 -#endif
1.34 +#include <lemon/assert.h>
1.35
1.36 namespace lemon {
1.37
1.38 - LpGlpk::LpGlpk() : Parent() {
1.39 - solved = false;
1.40 - rows = _lp_bits::LpId(1);
1.41 - cols = _lp_bits::LpId(1);
1.42 - lp = LEMON_glp(create_prob)();
1.43 - LEMON_glp(create_index)(lp);
1.44 - messageLevel(0);
1.45 + // GlpkBase members
1.46 +
1.47 + GlpkBase::GlpkBase() : LpBase() {
1.48 + lp = glp_create_prob();
1.49 + glp_create_index(lp);
1.50 }
1.51
1.52 - LpGlpk::LpGlpk(const LpGlpk &glp) : Parent() {
1.53 - solved = false;
1.54 - rows = _lp_bits::LpId(1);
1.55 - cols = _lp_bits::LpId(1);
1.56 - lp = LEMON_glp(create_prob)();
1.57 - LEMON_glp(create_index)(lp);
1.58 - messageLevel(0);
1.59 - //Coefficient matrix, row bounds
1.60 - LEMON_glp(add_rows)(lp, LEMON_glp(get_num_rows)(glp.lp));
1.61 - LEMON_glp(add_cols)(lp, LEMON_glp(get_num_cols)(glp.lp));
1.62 - int len;
1.63 - std::vector<int> ind(1+LEMON_glp(get_num_cols)(glp.lp));
1.64 - std::vector<Value> val(1+LEMON_glp(get_num_cols)(glp.lp));
1.65 - for (int i=1;i<=LEMON_glp(get_num_rows)(glp.lp);++i)
1.66 - {
1.67 - len=LEMON_glp(get_mat_row)(glp.lp,i,&*ind.begin(),&*val.begin());
1.68 - LEMON_glp(set_mat_row)(lp, i,len,&*ind.begin(),&*val.begin());
1.69 - LEMON_glp(set_row_bnds)(lp,i,
1.70 - LEMON_glp(get_row_type)(glp.lp,i),
1.71 - LEMON_glp(get_row_lb)(glp.lp,i),
1.72 - LEMON_glp(get_row_ub)(glp.lp,i));
1.73 - }
1.74 -
1.75 - //Objective function, coloumn bounds
1.76 - LEMON_glp(set_obj_dir)(lp, LEMON_glp(get_obj_dir)(glp.lp));
1.77 - //Objectif function's constant term treated separately
1.78 - LEMON_glp(set_obj_coef)(lp,0,LEMON_glp(get_obj_coef)(glp.lp,0));
1.79 - for (int i=1;i<=LEMON_glp(get_num_cols)(glp.lp);++i)
1.80 - {
1.81 - LEMON_glp(set_obj_coef)(lp,i,
1.82 - LEMON_glp(get_obj_coef)(glp.lp,i));
1.83 - LEMON_glp(set_col_bnds)(lp,i,
1.84 - LEMON_glp(get_col_type)(glp.lp,i),
1.85 - LEMON_glp(get_col_lb)(glp.lp,i),
1.86 - LEMON_glp(get_col_ub)(glp.lp,i));
1.87 - }
1.88 - rows = glp.rows;
1.89 - cols = glp.cols;
1.90 + GlpkBase::GlpkBase(const GlpkBase &other) : LpBase() {
1.91 + lp = glp_create_prob();
1.92 + glp_copy_prob(lp, other.lp, GLP_ON);
1.93 + glp_create_index(lp);
1.94 + rows = other.rows;
1.95 + cols = other.cols;
1.96 }
1.97
1.98 - LpGlpk::~LpGlpk() {
1.99 - LEMON_glp(delete_prob)(lp);
1.100 + GlpkBase::~GlpkBase() {
1.101 + glp_delete_prob(lp);
1.102 }
1.103
1.104 - int LpGlpk::_addCol() {
1.105 - int i=LEMON_glp(add_cols)(lp, 1);
1.106 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), 0.0, 0.0);
1.107 - solved = false;
1.108 + int GlpkBase::_addCol() {
1.109 + int i = glp_add_cols(lp, 1);
1.110 + glp_set_col_bnds(lp, i, GLP_FR, 0.0, 0.0);
1.111 return i;
1.112 }
1.113
1.114 - ///\e
1.115 -
1.116 -
1.117 - LpSolverBase* LpGlpk::_newLp()
1.118 - {
1.119 - LpGlpk* newlp = new LpGlpk;
1.120 - return newlp;
1.121 - }
1.122 -
1.123 - ///\e
1.124 -
1.125 - LpSolverBase* LpGlpk::_copyLp()
1.126 - {
1.127 - LpGlpk *newlp = new LpGlpk(*this);
1.128 - return newlp;
1.129 - }
1.130 -
1.131 - int LpGlpk::_addRow() {
1.132 - int i=LEMON_glp(add_rows)(lp, 1);
1.133 - solved = false;
1.134 + int GlpkBase::_addRow() {
1.135 + int i = glp_add_rows(lp, 1);
1.136 + glp_set_row_bnds(lp, i, GLP_FR, 0.0, 0.0);
1.137 return i;
1.138 }
1.139
1.140 -
1.141 - void LpGlpk::_eraseCol(int i) {
1.142 + void GlpkBase::_eraseCol(int i) {
1.143 int ca[2];
1.144 - ca[1]=i;
1.145 - LEMON_glp(del_cols)(lp, 1, ca);
1.146 - solved = false;
1.147 + ca[1] = i;
1.148 + glp_del_cols(lp, 1, ca);
1.149 }
1.150
1.151 - void LpGlpk::_eraseRow(int i) {
1.152 + void GlpkBase::_eraseRow(int i) {
1.153 int ra[2];
1.154 - ra[1]=i;
1.155 - LEMON_glp(del_rows)(lp, 1, ra);
1.156 - solved = false;
1.157 + ra[1] = i;
1.158 + glp_del_rows(lp, 1, ra);
1.159 }
1.160
1.161 - void LpGlpk::_getColName(int c, std::string & name) const
1.162 - {
1.163 -
1.164 - const char *n = LEMON_glp(get_col_name)(lp,c);
1.165 - name = n?n:"";
1.166 + void GlpkBase::_eraseColId(int i) {
1.167 + cols.eraseIndex(i);
1.168 + cols.shiftIndices(i);
1.169 }
1.170
1.171 + void GlpkBase::_eraseRowId(int i) {
1.172 + rows.eraseIndex(i);
1.173 + rows.shiftIndices(i);
1.174 + }
1.175
1.176 - void LpGlpk::_setColName(int c, const std::string & name)
1.177 - {
1.178 - LEMON_glp(set_col_name)(lp,c,const_cast<char*>(name.c_str()));
1.179 + void GlpkBase::_getColName(int c, std::string& name) const {
1.180 + const char *str = glp_get_col_name(lp, c);
1.181 + if (str) name = str;
1.182 + else name.clear();
1.183 + }
1.184 +
1.185 + void GlpkBase::_setColName(int c, const std::string & name) {
1.186 + glp_set_col_name(lp, c, const_cast<char*>(name.c_str()));
1.187
1.188 }
1.189
1.190 - int LpGlpk::_colByName(const std::string& name) const
1.191 - {
1.192 - int k = LEMON_glp(find_col)(lp, const_cast<char*>(name.c_str()));
1.193 + int GlpkBase::_colByName(const std::string& name) const {
1.194 + int k = glp_find_col(lp, const_cast<char*>(name.c_str()));
1.195 return k > 0 ? k : -1;
1.196 }
1.197
1.198 + void GlpkBase::_getRowName(int r, std::string& name) const {
1.199 + const char *str = glp_get_row_name(lp, r);
1.200 + if (str) name = str;
1.201 + else name.clear();
1.202 + }
1.203
1.204 - void LpGlpk::_setRowCoeffs(int i, ConstRowIterator b, ConstRowIterator e)
1.205 - {
1.206 - std::vector<int> indices;
1.207 + void GlpkBase::_setRowName(int r, const std::string & name) {
1.208 + glp_set_row_name(lp, r, const_cast<char*>(name.c_str()));
1.209 +
1.210 + }
1.211 +
1.212 + int GlpkBase::_rowByName(const std::string& name) const {
1.213 + int k = glp_find_row(lp, const_cast<char*>(name.c_str()));
1.214 + return k > 0 ? k : -1;
1.215 + }
1.216 +
1.217 + void GlpkBase::_setRowCoeffs(int i, ExprIterator b, ExprIterator e) {
1.218 + std::vector<int> indexes;
1.219 std::vector<Value> values;
1.220
1.221 - indices.push_back(0);
1.222 + indexes.push_back(0);
1.223 values.push_back(0);
1.224
1.225 - for(ConstRowIterator it=b; it!=e; ++it) {
1.226 - indices.push_back(it->first);
1.227 + for(ExprIterator it = b; it != e; ++it) {
1.228 + indexes.push_back(it->first);
1.229 values.push_back(it->second);
1.230 }
1.231
1.232 - LEMON_glp(set_mat_row)(lp, i, values.size() - 1,
1.233 - &indices[0], &values[0]);
1.234 -
1.235 - solved = false;
1.236 + glp_set_mat_row(lp, i, values.size() - 1,
1.237 + &indexes.front(), &values.front());
1.238 }
1.239
1.240 - void LpGlpk::_getRowCoeffs(int ix, RowIterator b) const
1.241 - {
1.242 - int length = LEMON_glp(get_mat_row)(lp, ix, 0, 0);
1.243 + void GlpkBase::_getRowCoeffs(int ix, InsertIterator b) const {
1.244 + int length = glp_get_mat_row(lp, ix, 0, 0);
1.245
1.246 - std::vector<int> indices(length + 1);
1.247 + std::vector<int> indexes(length + 1);
1.248 std::vector<Value> values(length + 1);
1.249
1.250 - LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);
1.251 + glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
1.252
1.253 for (int i = 1; i <= length; ++i) {
1.254 - *b = std::make_pair(indices[i], values[i]);
1.255 + *b = std::make_pair(indexes[i], values[i]);
1.256 ++b;
1.257 }
1.258 }
1.259
1.260 - void LpGlpk::_setColCoeffs(int ix, ConstColIterator b, ConstColIterator e) {
1.261 + void GlpkBase::_setColCoeffs(int ix, ExprIterator b,
1.262 + ExprIterator e) {
1.263
1.264 - std::vector<int> indices;
1.265 + std::vector<int> indexes;
1.266 std::vector<Value> values;
1.267
1.268 - indices.push_back(0);
1.269 + indexes.push_back(0);
1.270 values.push_back(0);
1.271
1.272 - for(ConstColIterator it=b; it!=e; ++it) {
1.273 - indices.push_back(it->first);
1.274 + for(ExprIterator it = b; it != e; ++it) {
1.275 + indexes.push_back(it->first);
1.276 values.push_back(it->second);
1.277 }
1.278
1.279 - LEMON_glp(set_mat_col)(lp, ix, values.size() - 1,
1.280 - &indices[0], &values[0]);
1.281 -
1.282 - solved = false;
1.283 + glp_set_mat_col(lp, ix, values.size() - 1,
1.284 + &indexes.front(), &values.front());
1.285 }
1.286
1.287 - void LpGlpk::_getColCoeffs(int ix, ColIterator b) const
1.288 - {
1.289 - int length = LEMON_glp(get_mat_col)(lp, ix, 0, 0);
1.290 + void GlpkBase::_getColCoeffs(int ix, InsertIterator b) const {
1.291 + int length = glp_get_mat_col(lp, ix, 0, 0);
1.292
1.293 - std::vector<int> indices(length + 1);
1.294 + std::vector<int> indexes(length + 1);
1.295 std::vector<Value> values(length + 1);
1.296
1.297 - LEMON_glp(get_mat_col)(lp, ix, &indices[0], &values[0]);
1.298 + glp_get_mat_col(lp, ix, &indexes.front(), &values.front());
1.299
1.300 - for (int i = 1; i <= length; ++i) {
1.301 - *b = std::make_pair(indices[i], values[i]);
1.302 + for (int i = 1; i <= length; ++i) {
1.303 + *b = std::make_pair(indexes[i], values[i]);
1.304 ++b;
1.305 }
1.306 }
1.307
1.308 - void LpGlpk::_setCoeff(int ix, int jx, Value value)
1.309 - {
1.310 + void GlpkBase::_setCoeff(int ix, int jx, Value value) {
1.311
1.312 - if (LEMON_glp(get_num_cols)(lp) < LEMON_glp(get_num_rows)(lp)) {
1.313 + if (glp_get_num_cols(lp) < glp_get_num_rows(lp)) {
1.314
1.315 - int length=LEMON_glp(get_mat_row)(lp, ix, 0, 0);
1.316 + int length = glp_get_mat_row(lp, ix, 0, 0);
1.317
1.318 - std::vector<int> indices(length + 2);
1.319 + std::vector<int> indexes(length + 2);
1.320 std::vector<Value> values(length + 2);
1.321
1.322 - LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);
1.323 + glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
1.324
1.325 //The following code does not suppose that the elements of the
1.326 - //array indices are sorted
1.327 - bool found=false;
1.328 - for (int i = 1; i <= length; ++i) {
1.329 - if (indices[i]==jx){
1.330 - found=true;
1.331 - values[i]=value;
1.332 + //array indexes are sorted
1.333 + bool found = false;
1.334 + for (int i = 1; i <= length; ++i) {
1.335 + if (indexes[i] == jx) {
1.336 + found = true;
1.337 + values[i] = value;
1.338 break;
1.339 }
1.340 }
1.341 - if (!found){
1.342 + if (!found) {
1.343 ++length;
1.344 - indices[length]=jx;
1.345 - values[length]=value;
1.346 + indexes[length] = jx;
1.347 + values[length] = value;
1.348 }
1.349
1.350 - LEMON_glp(set_mat_row)(lp, ix, length, &indices[0], &values[0]);
1.351 + glp_set_mat_row(lp, ix, length, &indexes.front(), &values.front());
1.352
1.353 } else {
1.354
1.355 - int length=LEMON_glp(get_mat_col)(lp, jx, 0, 0);
1.356 + int length = glp_get_mat_col(lp, jx, 0, 0);
1.357
1.358 - std::vector<int> indices(length + 2);
1.359 + std::vector<int> indexes(length + 2);
1.360 std::vector<Value> values(length + 2);
1.361
1.362 - LEMON_glp(get_mat_col)(lp, jx, &indices[0], &values[0]);
1.363 + glp_get_mat_col(lp, jx, &indexes.front(), &values.front());
1.364
1.365 //The following code does not suppose that the elements of the
1.366 - //array indices are sorted
1.367 - bool found=false;
1.368 + //array indexes are sorted
1.369 + bool found = false;
1.370 for (int i = 1; i <= length; ++i) {
1.371 - if (indices[i]==ix){
1.372 - found=true;
1.373 - values[i]=value;
1.374 + if (indexes[i] == ix) {
1.375 + found = true;
1.376 + values[i] = value;
1.377 break;
1.378 }
1.379 }
1.380 - if (!found){
1.381 + if (!found) {
1.382 ++length;
1.383 - indices[length]=ix;
1.384 - values[length]=value;
1.385 + indexes[length] = ix;
1.386 + values[length] = value;
1.387 }
1.388
1.389 - LEMON_glp(set_mat_col)(lp, jx, length, &indices[0], &values[0]);
1.390 + glp_set_mat_col(lp, jx, length, &indexes.front(), &values.front());
1.391 }
1.392
1.393 - solved = false;
1.394 }
1.395
1.396 - LpGlpk::Value LpGlpk::_getCoeff(int ix, int jx) const
1.397 - {
1.398 + GlpkBase::Value GlpkBase::_getCoeff(int ix, int jx) const {
1.399
1.400 - int length=LEMON_glp(get_mat_row)(lp, ix, 0, 0);
1.401 + int length = glp_get_mat_row(lp, ix, 0, 0);
1.402
1.403 - std::vector<int> indices(length + 1);
1.404 + std::vector<int> indexes(length + 1);
1.405 std::vector<Value> values(length + 1);
1.406
1.407 - LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);
1.408 + glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
1.409
1.410 - //The following code does not suppose that the elements of the
1.411 - //array indices are sorted
1.412 - for (int i = 1; i <= length; ++i) {
1.413 - if (indices[i]==jx){
1.414 + for (int i = 1; i <= length; ++i) {
1.415 + if (indexes[i] == jx) {
1.416 return values[i];
1.417 }
1.418 }
1.419 +
1.420 return 0;
1.421 + }
1.422 +
1.423 + void GlpkBase::_setColLowerBound(int i, Value lo) {
1.424 + LEMON_ASSERT(lo != INF, "Invalid bound");
1.425 +
1.426 + int b = glp_get_col_type(lp, i);
1.427 + double up = glp_get_col_ub(lp, i);
1.428 + if (lo == -INF) {
1.429 + switch (b) {
1.430 + case GLP_FR:
1.431 + case GLP_LO:
1.432 + glp_set_col_bnds(lp, i, GLP_FR, lo, up);
1.433 + break;
1.434 + case GLP_UP:
1.435 + break;
1.436 + case GLP_DB:
1.437 + case GLP_FX:
1.438 + glp_set_col_bnds(lp, i, GLP_UP, lo, up);
1.439 + break;
1.440 + default:
1.441 + break;
1.442 + }
1.443 + } else {
1.444 + switch (b) {
1.445 + case GLP_FR:
1.446 + case GLP_LO:
1.447 + glp_set_col_bnds(lp, i, GLP_LO, lo, up);
1.448 + break;
1.449 + case GLP_UP:
1.450 + case GLP_DB:
1.451 + case GLP_FX:
1.452 + if (lo == up)
1.453 + glp_set_col_bnds(lp, i, GLP_FX, lo, up);
1.454 + else
1.455 + glp_set_col_bnds(lp, i, GLP_DB, lo, up);
1.456 + break;
1.457 + default:
1.458 + break;
1.459 + }
1.460 + }
1.461 + }
1.462 +
1.463 + GlpkBase::Value GlpkBase::_getColLowerBound(int i) const {
1.464 + int b = glp_get_col_type(lp, i);
1.465 + switch (b) {
1.466 + case GLP_LO:
1.467 + case GLP_DB:
1.468 + case GLP_FX:
1.469 + return glp_get_col_lb(lp, i);
1.470 + default:
1.471 + return -INF;
1.472 + }
1.473 + }
1.474 +
1.475 + void GlpkBase::_setColUpperBound(int i, Value up) {
1.476 + LEMON_ASSERT(up != -INF, "Invalid bound");
1.477 +
1.478 + int b = glp_get_col_type(lp, i);
1.479 + double lo = glp_get_col_lb(lp, i);
1.480 + if (up == INF) {
1.481 + switch (b) {
1.482 + case GLP_FR:
1.483 + case GLP_LO:
1.484 + break;
1.485 + case GLP_UP:
1.486 + glp_set_col_bnds(lp, i, GLP_FR, lo, up);
1.487 + break;
1.488 + case GLP_DB:
1.489 + case GLP_FX:
1.490 + glp_set_col_bnds(lp, i, GLP_LO, lo, up);
1.491 + break;
1.492 + default:
1.493 + break;
1.494 + }
1.495 + } else {
1.496 + switch (b) {
1.497 + case GLP_FR:
1.498 + glp_set_col_bnds(lp, i, GLP_UP, lo, up);
1.499 + break;
1.500 + case GLP_UP:
1.501 + glp_set_col_bnds(lp, i, GLP_UP, lo, up);
1.502 + break;
1.503 + case GLP_LO:
1.504 + case GLP_DB:
1.505 + case GLP_FX:
1.506 + if (lo == up)
1.507 + glp_set_col_bnds(lp, i, GLP_FX, lo, up);
1.508 + else
1.509 + glp_set_col_bnds(lp, i, GLP_DB, lo, up);
1.510 + break;
1.511 + default:
1.512 + break;
1.513 + }
1.514 + }
1.515
1.516 }
1.517
1.518 -
1.519 - void LpGlpk::_setColLowerBound(int i, Value lo)
1.520 - {
1.521 - if (lo==INF) {
1.522 - //FIXME error
1.523 - }
1.524 - int b=LEMON_glp(get_col_type)(lp, i);
1.525 - double up=LEMON_glp(get_col_ub)(lp, i);
1.526 - if (lo==-INF) {
1.527 + GlpkBase::Value GlpkBase::_getColUpperBound(int i) const {
1.528 + int b = glp_get_col_type(lp, i);
1.529 switch (b) {
1.530 - case LEMON_GLP(FR):
1.531 - case LEMON_GLP(LO):
1.532 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), lo, up);
1.533 - break;
1.534 - case LEMON_GLP(UP):
1.535 - break;
1.536 - case LEMON_GLP(DB):
1.537 - case LEMON_GLP(FX):
1.538 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);
1.539 - break;
1.540 - default: ;
1.541 - //FIXME error
1.542 - }
1.543 - } else {
1.544 - switch (b) {
1.545 - case LEMON_GLP(FR):
1.546 - case LEMON_GLP(LO):
1.547 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(LO), lo, up);
1.548 - break;
1.549 - case LEMON_GLP(UP):
1.550 - case LEMON_GLP(DB):
1.551 - case LEMON_GLP(FX):
1.552 - if (lo==up)
1.553 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FX), lo, up);
1.554 - else
1.555 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(DB), lo, up);
1.556 - break;
1.557 - default: ;
1.558 - //FIXME error
1.559 - }
1.560 - }
1.561 -
1.562 - solved = false;
1.563 - }
1.564 -
1.565 - LpGlpk::Value LpGlpk::_getColLowerBound(int i) const
1.566 - {
1.567 - int b=LEMON_glp(get_col_type)(lp, i);
1.568 - switch (b) {
1.569 - case LEMON_GLP(LO):
1.570 - case LEMON_GLP(DB):
1.571 - case LEMON_GLP(FX):
1.572 - return LEMON_glp(get_col_lb)(lp, i);
1.573 - default: ;
1.574 - return -INF;
1.575 - }
1.576 - }
1.577 -
1.578 - void LpGlpk::_setColUpperBound(int i, Value up)
1.579 - {
1.580 - if (up==-INF) {
1.581 - //FIXME error
1.582 - }
1.583 - int b=LEMON_glp(get_col_type)(lp, i);
1.584 - double lo=LEMON_glp(get_col_lb)(lp, i);
1.585 - if (up==INF) {
1.586 - switch (b) {
1.587 - case LEMON_GLP(FR):
1.588 - case LEMON_GLP(LO):
1.589 - break;
1.590 - case LEMON_GLP(UP):
1.591 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), lo, up);
1.592 - break;
1.593 - case LEMON_GLP(DB):
1.594 - case LEMON_GLP(FX):
1.595 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(LO), lo, up);
1.596 - break;
1.597 - default: ;
1.598 - //FIXME error
1.599 - }
1.600 - } else {
1.601 - switch (b) {
1.602 - case LEMON_GLP(FR):
1.603 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);
1.604 - break;
1.605 - case LEMON_GLP(UP):
1.606 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);
1.607 - break;
1.608 - case LEMON_GLP(LO):
1.609 - case LEMON_GLP(DB):
1.610 - case LEMON_GLP(FX):
1.611 - if (lo==up)
1.612 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FX), lo, up);
1.613 - else
1.614 - LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(DB), lo, up);
1.615 - break;
1.616 - default: ;
1.617 - //FIXME error
1.618 - }
1.619 - }
1.620 -
1.621 - solved = false;
1.622 - }
1.623 -
1.624 - LpGlpk::Value LpGlpk::_getColUpperBound(int i) const
1.625 - {
1.626 - int b=LEMON_glp(get_col_type)(lp, i);
1.627 - switch (b) {
1.628 - case LEMON_GLP(UP):
1.629 - case LEMON_GLP(DB):
1.630 - case LEMON_GLP(FX):
1.631 - return LEMON_glp(get_col_ub)(lp, i);
1.632 - default: ;
1.633 + case GLP_UP:
1.634 + case GLP_DB:
1.635 + case GLP_FX:
1.636 + return glp_get_col_ub(lp, i);
1.637 + default:
1.638 return INF;
1.639 }
1.640 }
1.641
1.642 - void LpGlpk::_setRowBounds(int i, Value lb, Value ub)
1.643 - {
1.644 - //Bad parameter
1.645 - if (lb==INF || ub==-INF) {
1.646 - //FIXME error
1.647 - }
1.648 + void GlpkBase::_setRowLowerBound(int i, Value lo) {
1.649 + LEMON_ASSERT(lo != INF, "Invalid bound");
1.650
1.651 - if (lb == -INF){
1.652 - if (ub == INF){
1.653 - LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(FR), lb, ub);
1.654 + int b = glp_get_row_type(lp, i);
1.655 + double up = glp_get_row_ub(lp, i);
1.656 + if (lo == -INF) {
1.657 + switch (b) {
1.658 + case GLP_FR:
1.659 + case GLP_LO:
1.660 + glp_set_row_bnds(lp, i, GLP_FR, lo, up);
1.661 + break;
1.662 + case GLP_UP:
1.663 + break;
1.664 + case GLP_DB:
1.665 + case GLP_FX:
1.666 + glp_set_row_bnds(lp, i, GLP_UP, lo, up);
1.667 + break;
1.668 + default:
1.669 + break;
1.670 }
1.671 - else{
1.672 - LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(UP), lb, ub);
1.673 + } else {
1.674 + switch (b) {
1.675 + case GLP_FR:
1.676 + case GLP_LO:
1.677 + glp_set_row_bnds(lp, i, GLP_LO, lo, up);
1.678 + break;
1.679 + case GLP_UP:
1.680 + case GLP_DB:
1.681 + case GLP_FX:
1.682 + if (lo == up)
1.683 + glp_set_row_bnds(lp, i, GLP_FX, lo, up);
1.684 + else
1.685 + glp_set_row_bnds(lp, i, GLP_DB, lo, up);
1.686 + break;
1.687 + default:
1.688 + break;
1.689 }
1.690 }
1.691 - else{
1.692 - if (ub==INF){
1.693 - LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(LO), lb, ub);
1.694 -
1.695 - }
1.696 - else{
1.697 - if (lb == ub){
1.698 - LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(FX), lb, ub);
1.699 - }
1.700 - else{
1.701 - LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(DB), lb, ub);
1.702 - }
1.703 - }
1.704 - }
1.705 -
1.706 - solved = false;
1.707 - }
1.708 -
1.709 - void LpGlpk::_getRowBounds(int i, Value &lb, Value &ub) const
1.710 - {
1.711 -
1.712 - int b=LEMON_glp(get_row_type)(lp, i);
1.713 - switch (b) {
1.714 - case LEMON_GLP(FR):
1.715 - case LEMON_GLP(UP):
1.716 - lb = -INF;
1.717 - break;
1.718 - default:
1.719 - lb=LEMON_glp(get_row_lb)(lp, i);
1.720 - }
1.721 -
1.722 - switch (b) {
1.723 - case LEMON_GLP(FR):
1.724 - case LEMON_GLP(LO):
1.725 - ub = INF;
1.726 - break;
1.727 - default:
1.728 - ub=LEMON_glp(get_row_ub)(lp, i);
1.729 - }
1.730
1.731 }
1.732
1.733 - void LpGlpk::_setObjCoeff(int i, Value obj_coef)
1.734 - {
1.735 - //i=0 means the constant term (shift)
1.736 - LEMON_glp(set_obj_coef)(lp, i, obj_coef);
1.737 -
1.738 - solved = false;
1.739 - }
1.740 -
1.741 - LpGlpk::Value LpGlpk::_getObjCoeff(int i) const {
1.742 - //i=0 means the constant term (shift)
1.743 - return LEMON_glp(get_obj_coef)(lp, i);
1.744 - }
1.745 -
1.746 - void LpGlpk::_clearObj()
1.747 - {
1.748 - for (int i=0;i<=LEMON_glp(get_num_cols)(lp);++i){
1.749 - LEMON_glp(set_obj_coef)(lp, i, 0);
1.750 - }
1.751 -
1.752 - solved = false;
1.753 - }
1.754 -
1.755 - LpGlpk::SolveExitStatus LpGlpk::_solve()
1.756 - {
1.757 - // A way to check the problem to be solved
1.758 - //LEMON_glp(write_cpxlp(lp,"naittvan.cpx");
1.759 -
1.760 - LEMON_lpx(std_basis)(lp);
1.761 - int i = LEMON_lpx(simplex)(lp);
1.762 -
1.763 - switch (i) {
1.764 - case LEMON_LPX(E_OK):
1.765 - solved = true;
1.766 - return SOLVED;
1.767 + GlpkBase::Value GlpkBase::_getRowLowerBound(int i) const {
1.768 + int b = glp_get_row_type(lp, i);
1.769 + switch (b) {
1.770 + case GLP_LO:
1.771 + case GLP_DB:
1.772 + case GLP_FX:
1.773 + return glp_get_row_lb(lp, i);
1.774 default:
1.775 - return UNSOLVED;
1.776 + return -INF;
1.777 }
1.778 }
1.779
1.780 - LpGlpk::Value LpGlpk::_getPrimal(int i) const
1.781 - {
1.782 - return LEMON_glp(get_col_prim)(lp,i);
1.783 - }
1.784 + void GlpkBase::_setRowUpperBound(int i, Value up) {
1.785 + LEMON_ASSERT(up != -INF, "Invalid bound");
1.786
1.787 - LpGlpk::Value LpGlpk::_getDual(int i) const
1.788 - {
1.789 - return LEMON_glp(get_row_dual)(lp,i);
1.790 - }
1.791 -
1.792 - LpGlpk::Value LpGlpk::_getPrimalValue() const
1.793 - {
1.794 - return LEMON_glp(get_obj_val)(lp);
1.795 - }
1.796 - bool LpGlpk::_isBasicCol(int i) const
1.797 - {
1.798 - return (LEMON_glp(get_col_stat)(lp, i)==LEMON_GLP(BS));
1.799 - }
1.800 -
1.801 -
1.802 - LpGlpk::SolutionStatus LpGlpk::_getPrimalStatus() const
1.803 - {
1.804 - if (!solved) return UNDEFINED;
1.805 - int stat= LEMON_lpx(get_status)(lp);
1.806 - switch (stat) {
1.807 - case LEMON_LPX(UNDEF)://Undefined (no solve has been run yet)
1.808 - return UNDEFINED;
1.809 - case LEMON_LPX(NOFEAS)://There is no feasible solution (primal, I guess)
1.810 - case LEMON_LPX(INFEAS)://Infeasible
1.811 - return INFEASIBLE;
1.812 - case LEMON_LPX(UNBND)://Unbounded
1.813 - return INFINITE;
1.814 - case LEMON_LPX(FEAS)://Feasible
1.815 - return FEASIBLE;
1.816 - case LEMON_LPX(OPT)://Feasible
1.817 - return OPTIMAL;
1.818 - default:
1.819 - return UNDEFINED; //to avoid gcc warning
1.820 - //FIXME error
1.821 + int b = glp_get_row_type(lp, i);
1.822 + double lo = glp_get_row_lb(lp, i);
1.823 + if (up == INF) {
1.824 + switch (b) {
1.825 + case GLP_FR:
1.826 + case GLP_LO:
1.827 + break;
1.828 + case GLP_UP:
1.829 + glp_set_row_bnds(lp, i, GLP_FR, lo, up);
1.830 + break;
1.831 + case GLP_DB:
1.832 + case GLP_FX:
1.833 + glp_set_row_bnds(lp, i, GLP_LO, lo, up);
1.834 + break;
1.835 + default:
1.836 + break;
1.837 + }
1.838 + } else {
1.839 + switch (b) {
1.840 + case GLP_FR:
1.841 + glp_set_row_bnds(lp, i, GLP_UP, lo, up);
1.842 + break;
1.843 + case GLP_UP:
1.844 + glp_set_row_bnds(lp, i, GLP_UP, lo, up);
1.845 + break;
1.846 + case GLP_LO:
1.847 + case GLP_DB:
1.848 + case GLP_FX:
1.849 + if (lo == up)
1.850 + glp_set_row_bnds(lp, i, GLP_FX, lo, up);
1.851 + else
1.852 + glp_set_row_bnds(lp, i, GLP_DB, lo, up);
1.853 + break;
1.854 + default:
1.855 + break;
1.856 + }
1.857 }
1.858 }
1.859
1.860 - LpGlpk::SolutionStatus LpGlpk::_getDualStatus() const
1.861 - {
1.862 - if (!solved) return UNDEFINED;
1.863 - switch (LEMON_lpx(get_dual_stat)(lp)) {
1.864 - case LEMON_LPX(D_UNDEF)://Undefined (no solve has been run yet)
1.865 - return UNDEFINED;
1.866 - case LEMON_LPX(D_NOFEAS)://There is no dual feasible solution
1.867 -// case LEMON_LPX(D_INFEAS://Infeasible
1.868 - return INFEASIBLE;
1.869 - case LEMON_LPX(D_FEAS)://Feasible
1.870 - switch (LEMON_lpx(get_status)(lp)) {
1.871 - case LEMON_LPX(NOFEAS):
1.872 - return INFINITE;
1.873 - case LEMON_LPX(OPT):
1.874 - return OPTIMAL;
1.875 - default:
1.876 - return FEASIBLE;
1.877 - }
1.878 + GlpkBase::Value GlpkBase::_getRowUpperBound(int i) const {
1.879 + int b = glp_get_row_type(lp, i);
1.880 + switch (b) {
1.881 + case GLP_UP:
1.882 + case GLP_DB:
1.883 + case GLP_FX:
1.884 + return glp_get_row_ub(lp, i);
1.885 default:
1.886 - return UNDEFINED; //to avoid gcc warning
1.887 - //FIXME error
1.888 + return INF;
1.889 }
1.890 }
1.891
1.892 - LpGlpk::ProblemTypes LpGlpk::_getProblemType() const
1.893 - {
1.894 - if (!solved) return UNKNOWN;
1.895 - //int stat= LEMON_glp(get_status(lp);
1.896 - int statp= LEMON_lpx(get_prim_stat)(lp);
1.897 - int statd= LEMON_lpx(get_dual_stat)(lp);
1.898 - if (statp==LEMON_LPX(P_FEAS) && statd==LEMON_LPX(D_FEAS))
1.899 - return PRIMAL_DUAL_FEASIBLE;
1.900 - if (statp==LEMON_LPX(P_FEAS) && statd==LEMON_LPX(D_NOFEAS))
1.901 - return PRIMAL_FEASIBLE_DUAL_INFEASIBLE;
1.902 - if (statp==LEMON_LPX(P_NOFEAS) && statd==LEMON_LPX(D_FEAS))
1.903 - return PRIMAL_INFEASIBLE_DUAL_FEASIBLE;
1.904 - if (statp==LEMON_LPX(P_NOFEAS) && statd==LEMON_LPX(D_NOFEAS))
1.905 - return PRIMAL_DUAL_INFEASIBLE;
1.906 - //In all other cases
1.907 - return UNKNOWN;
1.908 + void GlpkBase::_setObjCoeffs(ExprIterator b, ExprIterator e) {
1.909 + for (int i = 1; i <= glp_get_num_cols(lp); ++i) {
1.910 + glp_set_obj_coef(lp, i, 0.0);
1.911 + }
1.912 + for (ExprIterator it = b; it != e; ++it) {
1.913 + glp_set_obj_coef(lp, it->first, it->second);
1.914 + }
1.915 }
1.916
1.917 - void LpGlpk::_setMax()
1.918 - {
1.919 - solved = false;
1.920 - LEMON_glp(set_obj_dir)(lp, LEMON_GLP(MAX));
1.921 + void GlpkBase::_getObjCoeffs(InsertIterator b) const {
1.922 + for (int i = 1; i <= glp_get_num_cols(lp); ++i) {
1.923 + Value val = glp_get_obj_coef(lp, i);
1.924 + if (val != 0.0) {
1.925 + *b = std::make_pair(i, val);
1.926 + ++b;
1.927 + }
1.928 + }
1.929 }
1.930
1.931 - void LpGlpk::_setMin()
1.932 - {
1.933 - solved = false;
1.934 - LEMON_glp(set_obj_dir)(lp, LEMON_GLP(MIN));
1.935 + void GlpkBase::_setObjCoeff(int i, Value obj_coef) {
1.936 + //i = 0 means the constant term (shift)
1.937 + glp_set_obj_coef(lp, i, obj_coef);
1.938 }
1.939
1.940 - bool LpGlpk::_isMax() const
1.941 - {
1.942 - return (LEMON_glp(get_obj_dir)(lp)==LEMON_GLP(MAX));
1.943 + GlpkBase::Value GlpkBase::_getObjCoeff(int i) const {
1.944 + //i = 0 means the constant term (shift)
1.945 + return glp_get_obj_coef(lp, i);
1.946 }
1.947
1.948 -
1.949 -
1.950 - void LpGlpk::messageLevel(int m)
1.951 - {
1.952 - LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_MSGLEV), m);
1.953 + void GlpkBase::_setSense(GlpkBase::Sense sense) {
1.954 + switch (sense) {
1.955 + case MIN:
1.956 + glp_set_obj_dir(lp, GLP_MIN);
1.957 + break;
1.958 + case MAX:
1.959 + glp_set_obj_dir(lp, GLP_MAX);
1.960 + break;
1.961 + }
1.962 }
1.963
1.964 - void LpGlpk::presolver(bool b)
1.965 - {
1.966 - LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_PRESOL), b);
1.967 + GlpkBase::Sense GlpkBase::_getSense() const {
1.968 + switch(glp_get_obj_dir(lp)) {
1.969 + case GLP_MIN:
1.970 + return MIN;
1.971 + case GLP_MAX:
1.972 + return MAX;
1.973 + default:
1.974 + LEMON_ASSERT(false, "Wrong sense");
1.975 + return GlpkBase::Sense();
1.976 + }
1.977 }
1.978
1.979 + void GlpkBase::_clear() {
1.980 + glp_erase_prob(lp);
1.981 + rows.clear();
1.982 + cols.clear();
1.983 + }
1.984 +
1.985 + // LpGlpk members
1.986 +
1.987 + LpGlpk::LpGlpk()
1.988 + : LpBase(), GlpkBase(), LpSolver() {
1.989 + messageLevel(MESSAGE_NO_OUTPUT);
1.990 + }
1.991 +
1.992 + LpGlpk::LpGlpk(const LpGlpk& other)
1.993 + : LpBase(other), GlpkBase(other), LpSolver(other) {
1.994 + messageLevel(MESSAGE_NO_OUTPUT);
1.995 + }
1.996 +
1.997 + LpGlpk* LpGlpk::_newSolver() const { return new LpGlpk; }
1.998 + LpGlpk* LpGlpk::_cloneSolver() const { return new LpGlpk(*this); }
1.999 +
1.1000 + const char* LpGlpk::_solverName() const { return "LpGlpk"; }
1.1001 +
1.1002 + void LpGlpk::_clear_temporals() {
1.1003 + _primal_ray.clear();
1.1004 + _dual_ray.clear();
1.1005 + }
1.1006 +
1.1007 + LpGlpk::SolveExitStatus LpGlpk::_solve() {
1.1008 + return solvePrimal();
1.1009 + }
1.1010 +
1.1011 + LpGlpk::SolveExitStatus LpGlpk::solvePrimal() {
1.1012 + _clear_temporals();
1.1013 +
1.1014 + glp_smcp smcp;
1.1015 + glp_init_smcp(&smcp);
1.1016 +
1.1017 + switch (_message_level) {
1.1018 + case MESSAGE_NO_OUTPUT:
1.1019 + smcp.msg_lev = GLP_MSG_OFF;
1.1020 + break;
1.1021 + case MESSAGE_ERROR_MESSAGE:
1.1022 + smcp.msg_lev = GLP_MSG_ERR;
1.1023 + break;
1.1024 + case MESSAGE_NORMAL_OUTPUT:
1.1025 + smcp.msg_lev = GLP_MSG_ON;
1.1026 + break;
1.1027 + case MESSAGE_FULL_OUTPUT:
1.1028 + smcp.msg_lev = GLP_MSG_ALL;
1.1029 + break;
1.1030 + }
1.1031 +
1.1032 + if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
1.1033 + return SOLVED;
1.1034 + }
1.1035 +
1.1036 + LpGlpk::SolveExitStatus LpGlpk::solveDual() {
1.1037 + _clear_temporals();
1.1038 +
1.1039 + glp_smcp smcp;
1.1040 + glp_init_smcp(&smcp);
1.1041 +
1.1042 + switch (_message_level) {
1.1043 + case MESSAGE_NO_OUTPUT:
1.1044 + smcp.msg_lev = GLP_MSG_OFF;
1.1045 + break;
1.1046 + case MESSAGE_ERROR_MESSAGE:
1.1047 + smcp.msg_lev = GLP_MSG_ERR;
1.1048 + break;
1.1049 + case MESSAGE_NORMAL_OUTPUT:
1.1050 + smcp.msg_lev = GLP_MSG_ON;
1.1051 + break;
1.1052 + case MESSAGE_FULL_OUTPUT:
1.1053 + smcp.msg_lev = GLP_MSG_ALL;
1.1054 + break;
1.1055 + }
1.1056 + smcp.meth = GLP_DUAL;
1.1057 +
1.1058 + if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
1.1059 + return SOLVED;
1.1060 + }
1.1061 +
1.1062 + LpGlpk::Value LpGlpk::_getPrimal(int i) const {
1.1063 + return glp_get_col_prim(lp, i);
1.1064 + }
1.1065 +
1.1066 + LpGlpk::Value LpGlpk::_getDual(int i) const {
1.1067 + return glp_get_row_dual(lp, i);
1.1068 + }
1.1069 +
1.1070 + LpGlpk::Value LpGlpk::_getPrimalValue() const {
1.1071 + return glp_get_obj_val(lp);
1.1072 + }
1.1073 +
1.1074 + LpGlpk::VarStatus LpGlpk::_getColStatus(int i) const {
1.1075 + switch (glp_get_col_stat(lp, i)) {
1.1076 + case GLP_BS:
1.1077 + return BASIC;
1.1078 + case GLP_UP:
1.1079 + return UPPER;
1.1080 + case GLP_LO:
1.1081 + return LOWER;
1.1082 + case GLP_NF:
1.1083 + return FREE;
1.1084 + case GLP_NS:
1.1085 + return FIXED;
1.1086 + default:
1.1087 + LEMON_ASSERT(false, "Wrong column status");
1.1088 + return LpGlpk::VarStatus();
1.1089 + }
1.1090 + }
1.1091 +
1.1092 + LpGlpk::VarStatus LpGlpk::_getRowStatus(int i) const {
1.1093 + switch (glp_get_row_stat(lp, i)) {
1.1094 + case GLP_BS:
1.1095 + return BASIC;
1.1096 + case GLP_UP:
1.1097 + return UPPER;
1.1098 + case GLP_LO:
1.1099 + return LOWER;
1.1100 + case GLP_NF:
1.1101 + return FREE;
1.1102 + case GLP_NS:
1.1103 + return FIXED;
1.1104 + default:
1.1105 + LEMON_ASSERT(false, "Wrong row status");
1.1106 + return LpGlpk::VarStatus();
1.1107 + }
1.1108 + }
1.1109 +
1.1110 + LpGlpk::Value LpGlpk::_getPrimalRay(int i) const {
1.1111 + if (_primal_ray.empty()) {
1.1112 + int row_num = glp_get_num_rows(lp);
1.1113 + int col_num = glp_get_num_cols(lp);
1.1114 +
1.1115 + _primal_ray.resize(col_num + 1, 0.0);
1.1116 +
1.1117 + int index = glp_get_unbnd_ray(lp);
1.1118 + if (index != 0) {
1.1119 + // The primal ray is found in primal simplex second phase
1.1120 + LEMON_ASSERT((index <= row_num ? glp_get_row_stat(lp, index) :
1.1121 + glp_get_col_stat(lp, index - row_num)) != GLP_BS,
1.1122 + "Wrong primal ray");
1.1123 +
1.1124 + bool negate = glp_get_obj_dir(lp) == GLP_MAX;
1.1125 +
1.1126 + if (index > row_num) {
1.1127 + _primal_ray[index - row_num] = 1.0;
1.1128 + if (glp_get_col_dual(lp, index - row_num) > 0) {
1.1129 + negate = !negate;
1.1130 + }
1.1131 + } else {
1.1132 + if (glp_get_row_dual(lp, index) > 0) {
1.1133 + negate = !negate;
1.1134 + }
1.1135 + }
1.1136 +
1.1137 + std::vector<int> ray_indexes(row_num + 1);
1.1138 + std::vector<Value> ray_values(row_num + 1);
1.1139 + int ray_length = glp_eval_tab_col(lp, index, &ray_indexes.front(),
1.1140 + &ray_values.front());
1.1141 +
1.1142 + for (int i = 1; i <= ray_length; ++i) {
1.1143 + if (ray_indexes[i] > row_num) {
1.1144 + _primal_ray[ray_indexes[i] - row_num] = ray_values[i];
1.1145 + }
1.1146 + }
1.1147 +
1.1148 + if (negate) {
1.1149 + for (int i = 1; i <= col_num; ++i) {
1.1150 + _primal_ray[i] = - _primal_ray[i];
1.1151 + }
1.1152 + }
1.1153 + } else {
1.1154 + for (int i = 1; i <= col_num; ++i) {
1.1155 + _primal_ray[i] = glp_get_col_prim(lp, i);
1.1156 + }
1.1157 + }
1.1158 + }
1.1159 + return _primal_ray[i];
1.1160 + }
1.1161 +
1.1162 + LpGlpk::Value LpGlpk::_getDualRay(int i) const {
1.1163 + if (_dual_ray.empty()) {
1.1164 + int row_num = glp_get_num_rows(lp);
1.1165 +
1.1166 + _dual_ray.resize(row_num + 1, 0.0);
1.1167 +
1.1168 + int index = glp_get_unbnd_ray(lp);
1.1169 + if (index != 0) {
1.1170 + // The dual ray is found in dual simplex second phase
1.1171 + LEMON_ASSERT((index <= row_num ? glp_get_row_stat(lp, index) :
1.1172 + glp_get_col_stat(lp, index - row_num)) == GLP_BS,
1.1173 +
1.1174 + "Wrong dual ray");
1.1175 +
1.1176 + int idx;
1.1177 + bool negate = false;
1.1178 +
1.1179 + if (index > row_num) {
1.1180 + idx = glp_get_col_bind(lp, index - row_num);
1.1181 + if (glp_get_col_prim(lp, index - row_num) >
1.1182 + glp_get_col_ub(lp, index - row_num)) {
1.1183 + negate = true;
1.1184 + }
1.1185 + } else {
1.1186 + idx = glp_get_row_bind(lp, index);
1.1187 + if (glp_get_row_prim(lp, index) > glp_get_row_ub(lp, index)) {
1.1188 + negate = true;
1.1189 + }
1.1190 + }
1.1191 +
1.1192 + _dual_ray[idx] = negate ? - 1.0 : 1.0;
1.1193 +
1.1194 + glp_btran(lp, &_dual_ray.front());
1.1195 + } else {
1.1196 + double eps = 1e-7;
1.1197 + // The dual ray is found in primal simplex first phase
1.1198 + // We assume that the glpk minimizes the slack to get feasible solution
1.1199 + for (int i = 1; i <= row_num; ++i) {
1.1200 + int index = glp_get_bhead(lp, i);
1.1201 + if (index <= row_num) {
1.1202 + double res = glp_get_row_prim(lp, index);
1.1203 + if (res > glp_get_row_ub(lp, index) + eps) {
1.1204 + _dual_ray[i] = -1;
1.1205 + } else if (res < glp_get_row_lb(lp, index) - eps) {
1.1206 + _dual_ray[i] = 1;
1.1207 + } else {
1.1208 + _dual_ray[i] = 0;
1.1209 + }
1.1210 + _dual_ray[i] *= glp_get_rii(lp, index);
1.1211 + } else {
1.1212 + double res = glp_get_col_prim(lp, index - row_num);
1.1213 + if (res > glp_get_col_ub(lp, index - row_num) + eps) {
1.1214 + _dual_ray[i] = -1;
1.1215 + } else if (res < glp_get_col_lb(lp, index - row_num) - eps) {
1.1216 + _dual_ray[i] = 1;
1.1217 + } else {
1.1218 + _dual_ray[i] = 0;
1.1219 + }
1.1220 + _dual_ray[i] /= glp_get_sjj(lp, index - row_num);
1.1221 + }
1.1222 + }
1.1223 +
1.1224 + glp_btran(lp, &_dual_ray.front());
1.1225 +
1.1226 + for (int i = 1; i <= row_num; ++i) {
1.1227 + _dual_ray[i] /= glp_get_rii(lp, i);
1.1228 + }
1.1229 + }
1.1230 + }
1.1231 + return _dual_ray[i];
1.1232 + }
1.1233 +
1.1234 + LpGlpk::ProblemType LpGlpk::_getPrimalType() const {
1.1235 + if (glp_get_status(lp) == GLP_OPT)
1.1236 + return OPTIMAL;
1.1237 + switch (glp_get_prim_stat(lp)) {
1.1238 + case GLP_UNDEF:
1.1239 + return UNDEFINED;
1.1240 + case GLP_FEAS:
1.1241 + case GLP_INFEAS:
1.1242 + if (glp_get_dual_stat(lp) == GLP_NOFEAS) {
1.1243 + return UNBOUNDED;
1.1244 + } else {
1.1245 + return UNDEFINED;
1.1246 + }
1.1247 + case GLP_NOFEAS:
1.1248 + return INFEASIBLE;
1.1249 + default:
1.1250 + LEMON_ASSERT(false, "Wrong primal type");
1.1251 + return LpGlpk::ProblemType();
1.1252 + }
1.1253 + }
1.1254 +
1.1255 + LpGlpk::ProblemType LpGlpk::_getDualType() const {
1.1256 + if (glp_get_status(lp) == GLP_OPT)
1.1257 + return OPTIMAL;
1.1258 + switch (glp_get_dual_stat(lp)) {
1.1259 + case GLP_UNDEF:
1.1260 + return UNDEFINED;
1.1261 + case GLP_FEAS:
1.1262 + case GLP_INFEAS:
1.1263 + if (glp_get_prim_stat(lp) == GLP_NOFEAS) {
1.1264 + return UNBOUNDED;
1.1265 + } else {
1.1266 + return UNDEFINED;
1.1267 + }
1.1268 + case GLP_NOFEAS:
1.1269 + return INFEASIBLE;
1.1270 + default:
1.1271 + LEMON_ASSERT(false, "Wrong primal type");
1.1272 + return LpGlpk::ProblemType();
1.1273 + }
1.1274 + }
1.1275 +
1.1276 + void LpGlpk::presolver(bool b) {
1.1277 + lpx_set_int_parm(lp, LPX_K_PRESOL, b ? 1 : 0);
1.1278 + }
1.1279 +
1.1280 + void LpGlpk::messageLevel(MessageLevel m) {
1.1281 + _message_level = m;
1.1282 + }
1.1283 +
1.1284 + // MipGlpk members
1.1285 +
1.1286 + MipGlpk::MipGlpk()
1.1287 + : LpBase(), GlpkBase(), MipSolver() {
1.1288 + messageLevel(MESSAGE_NO_OUTPUT);
1.1289 + }
1.1290 +
1.1291 + MipGlpk::MipGlpk(const MipGlpk& other)
1.1292 + : LpBase(), GlpkBase(other), MipSolver() {
1.1293 + messageLevel(MESSAGE_NO_OUTPUT);
1.1294 + }
1.1295 +
1.1296 + void MipGlpk::_setColType(int i, MipGlpk::ColTypes col_type) {
1.1297 + switch (col_type) {
1.1298 + case INTEGER:
1.1299 + glp_set_col_kind(lp, i, GLP_IV);
1.1300 + break;
1.1301 + case REAL:
1.1302 + glp_set_col_kind(lp, i, GLP_CV);
1.1303 + break;
1.1304 + }
1.1305 + }
1.1306 +
1.1307 + MipGlpk::ColTypes MipGlpk::_getColType(int i) const {
1.1308 + switch (glp_get_col_kind(lp, i)) {
1.1309 + case GLP_IV:
1.1310 + case GLP_BV:
1.1311 + return INTEGER;
1.1312 + default:
1.1313 + return REAL;
1.1314 + }
1.1315 +
1.1316 + }
1.1317 +
1.1318 + MipGlpk::SolveExitStatus MipGlpk::_solve() {
1.1319 + glp_smcp smcp;
1.1320 + glp_init_smcp(&smcp);
1.1321 +
1.1322 + switch (_message_level) {
1.1323 + case MESSAGE_NO_OUTPUT:
1.1324 + smcp.msg_lev = GLP_MSG_OFF;
1.1325 + break;
1.1326 + case MESSAGE_ERROR_MESSAGE:
1.1327 + smcp.msg_lev = GLP_MSG_ERR;
1.1328 + break;
1.1329 + case MESSAGE_NORMAL_OUTPUT:
1.1330 + smcp.msg_lev = GLP_MSG_ON;
1.1331 + break;
1.1332 + case MESSAGE_FULL_OUTPUT:
1.1333 + smcp.msg_lev = GLP_MSG_ALL;
1.1334 + break;
1.1335 + }
1.1336 + smcp.meth = GLP_DUAL;
1.1337 +
1.1338 + if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
1.1339 + if (glp_get_status(lp) != GLP_OPT) return SOLVED;
1.1340 +
1.1341 + glp_iocp iocp;
1.1342 + glp_init_iocp(&iocp);
1.1343 +
1.1344 + switch (_message_level) {
1.1345 + case MESSAGE_NO_OUTPUT:
1.1346 + iocp.msg_lev = GLP_MSG_OFF;
1.1347 + break;
1.1348 + case MESSAGE_ERROR_MESSAGE:
1.1349 + iocp.msg_lev = GLP_MSG_ERR;
1.1350 + break;
1.1351 + case MESSAGE_NORMAL_OUTPUT:
1.1352 + iocp.msg_lev = GLP_MSG_ON;
1.1353 + break;
1.1354 + case MESSAGE_FULL_OUTPUT:
1.1355 + iocp.msg_lev = GLP_MSG_ALL;
1.1356 + break;
1.1357 + }
1.1358 +
1.1359 + if (glp_intopt(lp, &iocp) != 0) return UNSOLVED;
1.1360 + return SOLVED;
1.1361 + }
1.1362 +
1.1363 +
1.1364 + MipGlpk::ProblemType MipGlpk::_getType() const {
1.1365 + switch (glp_get_status(lp)) {
1.1366 + case GLP_OPT:
1.1367 + switch (glp_mip_status(lp)) {
1.1368 + case GLP_UNDEF:
1.1369 + return UNDEFINED;
1.1370 + case GLP_NOFEAS:
1.1371 + return INFEASIBLE;
1.1372 + case GLP_FEAS:
1.1373 + return FEASIBLE;
1.1374 + case GLP_OPT:
1.1375 + return OPTIMAL;
1.1376 + default:
1.1377 + LEMON_ASSERT(false, "Wrong problem type.");
1.1378 + return MipGlpk::ProblemType();
1.1379 + }
1.1380 + case GLP_NOFEAS:
1.1381 + return INFEASIBLE;
1.1382 + case GLP_INFEAS:
1.1383 + case GLP_FEAS:
1.1384 + if (glp_get_dual_stat(lp) == GLP_NOFEAS) {
1.1385 + return UNBOUNDED;
1.1386 + } else {
1.1387 + return UNDEFINED;
1.1388 + }
1.1389 + default:
1.1390 + LEMON_ASSERT(false, "Wrong problem type.");
1.1391 + return MipGlpk::ProblemType();
1.1392 + }
1.1393 + }
1.1394 +
1.1395 + MipGlpk::Value MipGlpk::_getSol(int i) const {
1.1396 + return glp_mip_col_val(lp, i);
1.1397 + }
1.1398 +
1.1399 + MipGlpk::Value MipGlpk::_getSolValue() const {
1.1400 + return glp_mip_obj_val(lp);
1.1401 + }
1.1402 +
1.1403 + MipGlpk* MipGlpk::_newSolver() const { return new MipGlpk; }
1.1404 + MipGlpk* MipGlpk::_cloneSolver() const {return new MipGlpk(*this); }
1.1405 +
1.1406 + const char* MipGlpk::_solverName() const { return "MipGlpk"; }
1.1407 +
1.1408 + void MipGlpk::messageLevel(MessageLevel m) {
1.1409 + _message_level = m;
1.1410 + }
1.1411
1.1412 } //END OF NAMESPACE LEMON