1.1 --- a/lemon/lp_glpk.cc Mon Jan 12 12:25:55 2009 +0000
1.2 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000
1.3 @@ -1,952 +0,0 @@
1.4 -/* -*- mode: C++; indent-tabs-mode: nil; -*-
1.5 - *
1.6 - * This file is a part of LEMON, a generic C++ optimization library.
1.7 - *
1.8 - * Copyright (C) 2003-2008
1.9 - * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
1.10 - * (Egervary Research Group on Combinatorial Optimization, EGRES).
1.11 - *
1.12 - * Permission to use, modify and distribute this software is granted
1.13 - * provided that this copyright notice appears in all copies. For
1.14 - * precise terms see the accompanying LICENSE file.
1.15 - *
1.16 - * This software is provided "AS IS" with no warranty of any kind,
1.17 - * express or implied, and with no claim as to its suitability for any
1.18 - * purpose.
1.19 - *
1.20 - */
1.21 -
1.22 -///\file
1.23 -///\brief Implementation of the LEMON GLPK LP and MIP solver interface.
1.24 -
1.25 -#include <lemon/lp_glpk.h>
1.26 -#include <glpk.h>
1.27 -
1.28 -#include <lemon/assert.h>
1.29 -
1.30 -namespace lemon {
1.31 -
1.32 - // GlpkBase members
1.33 -
1.34 - GlpkBase::GlpkBase() : LpBase() {
1.35 - lp = glp_create_prob();
1.36 - glp_create_index(lp);
1.37 - }
1.38 -
1.39 - GlpkBase::GlpkBase(const GlpkBase &other) : LpBase() {
1.40 - lp = glp_create_prob();
1.41 - glp_copy_prob(lp, other.lp, GLP_ON);
1.42 - glp_create_index(lp);
1.43 - rows = other.rows;
1.44 - cols = other.cols;
1.45 - }
1.46 -
1.47 - GlpkBase::~GlpkBase() {
1.48 - glp_delete_prob(lp);
1.49 - }
1.50 -
1.51 - int GlpkBase::_addCol() {
1.52 - int i = glp_add_cols(lp, 1);
1.53 - glp_set_col_bnds(lp, i, GLP_FR, 0.0, 0.0);
1.54 - return i;
1.55 - }
1.56 -
1.57 - int GlpkBase::_addRow() {
1.58 - int i = glp_add_rows(lp, 1);
1.59 - glp_set_row_bnds(lp, i, GLP_FR, 0.0, 0.0);
1.60 - return i;
1.61 - }
1.62 -
1.63 - void GlpkBase::_eraseCol(int i) {
1.64 - int ca[2];
1.65 - ca[1] = i;
1.66 - glp_del_cols(lp, 1, ca);
1.67 - }
1.68 -
1.69 - void GlpkBase::_eraseRow(int i) {
1.70 - int ra[2];
1.71 - ra[1] = i;
1.72 - glp_del_rows(lp, 1, ra);
1.73 - }
1.74 -
1.75 - void GlpkBase::_eraseColId(int i) {
1.76 - cols.eraseIndex(i);
1.77 - cols.shiftIndices(i);
1.78 - }
1.79 -
1.80 - void GlpkBase::_eraseRowId(int i) {
1.81 - rows.eraseIndex(i);
1.82 - rows.shiftIndices(i);
1.83 - }
1.84 -
1.85 - void GlpkBase::_getColName(int c, std::string& name) const {
1.86 - const char *str = glp_get_col_name(lp, c);
1.87 - if (str) name = str;
1.88 - else name.clear();
1.89 - }
1.90 -
1.91 - void GlpkBase::_setColName(int c, const std::string & name) {
1.92 - glp_set_col_name(lp, c, const_cast<char*>(name.c_str()));
1.93 -
1.94 - }
1.95 -
1.96 - int GlpkBase::_colByName(const std::string& name) const {
1.97 - int k = glp_find_col(lp, const_cast<char*>(name.c_str()));
1.98 - return k > 0 ? k : -1;
1.99 - }
1.100 -
1.101 - void GlpkBase::_getRowName(int r, std::string& name) const {
1.102 - const char *str = glp_get_row_name(lp, r);
1.103 - if (str) name = str;
1.104 - else name.clear();
1.105 - }
1.106 -
1.107 - void GlpkBase::_setRowName(int r, const std::string & name) {
1.108 - glp_set_row_name(lp, r, const_cast<char*>(name.c_str()));
1.109 -
1.110 - }
1.111 -
1.112 - int GlpkBase::_rowByName(const std::string& name) const {
1.113 - int k = glp_find_row(lp, const_cast<char*>(name.c_str()));
1.114 - return k > 0 ? k : -1;
1.115 - }
1.116 -
1.117 - void GlpkBase::_setRowCoeffs(int i, ExprIterator b, ExprIterator e) {
1.118 - std::vector<int> indexes;
1.119 - std::vector<Value> values;
1.120 -
1.121 - indexes.push_back(0);
1.122 - values.push_back(0);
1.123 -
1.124 - for(ExprIterator it = b; it != e; ++it) {
1.125 - indexes.push_back(it->first);
1.126 - values.push_back(it->second);
1.127 - }
1.128 -
1.129 - glp_set_mat_row(lp, i, values.size() - 1,
1.130 - &indexes.front(), &values.front());
1.131 - }
1.132 -
1.133 - void GlpkBase::_getRowCoeffs(int ix, InsertIterator b) const {
1.134 - int length = glp_get_mat_row(lp, ix, 0, 0);
1.135 -
1.136 - std::vector<int> indexes(length + 1);
1.137 - std::vector<Value> values(length + 1);
1.138 -
1.139 - glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
1.140 -
1.141 - for (int i = 1; i <= length; ++i) {
1.142 - *b = std::make_pair(indexes[i], values[i]);
1.143 - ++b;
1.144 - }
1.145 - }
1.146 -
1.147 - void GlpkBase::_setColCoeffs(int ix, ExprIterator b,
1.148 - ExprIterator e) {
1.149 -
1.150 - std::vector<int> indexes;
1.151 - std::vector<Value> values;
1.152 -
1.153 - indexes.push_back(0);
1.154 - values.push_back(0);
1.155 -
1.156 - for(ExprIterator it = b; it != e; ++it) {
1.157 - indexes.push_back(it->first);
1.158 - values.push_back(it->second);
1.159 - }
1.160 -
1.161 - glp_set_mat_col(lp, ix, values.size() - 1,
1.162 - &indexes.front(), &values.front());
1.163 - }
1.164 -
1.165 - void GlpkBase::_getColCoeffs(int ix, InsertIterator b) const {
1.166 - int length = glp_get_mat_col(lp, ix, 0, 0);
1.167 -
1.168 - std::vector<int> indexes(length + 1);
1.169 - std::vector<Value> values(length + 1);
1.170 -
1.171 - glp_get_mat_col(lp, ix, &indexes.front(), &values.front());
1.172 -
1.173 - for (int i = 1; i <= length; ++i) {
1.174 - *b = std::make_pair(indexes[i], values[i]);
1.175 - ++b;
1.176 - }
1.177 - }
1.178 -
1.179 - void GlpkBase::_setCoeff(int ix, int jx, Value value) {
1.180 -
1.181 - if (glp_get_num_cols(lp) < glp_get_num_rows(lp)) {
1.182 -
1.183 - int length = glp_get_mat_row(lp, ix, 0, 0);
1.184 -
1.185 - std::vector<int> indexes(length + 2);
1.186 - std::vector<Value> values(length + 2);
1.187 -
1.188 - glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
1.189 -
1.190 - //The following code does not suppose that the elements of the
1.191 - //array indexes are sorted
1.192 - bool found = false;
1.193 - for (int i = 1; i <= length; ++i) {
1.194 - if (indexes[i] == jx) {
1.195 - found = true;
1.196 - values[i] = value;
1.197 - break;
1.198 - }
1.199 - }
1.200 - if (!found) {
1.201 - ++length;
1.202 - indexes[length] = jx;
1.203 - values[length] = value;
1.204 - }
1.205 -
1.206 - glp_set_mat_row(lp, ix, length, &indexes.front(), &values.front());
1.207 -
1.208 - } else {
1.209 -
1.210 - int length = glp_get_mat_col(lp, jx, 0, 0);
1.211 -
1.212 - std::vector<int> indexes(length + 2);
1.213 - std::vector<Value> values(length + 2);
1.214 -
1.215 - glp_get_mat_col(lp, jx, &indexes.front(), &values.front());
1.216 -
1.217 - //The following code does not suppose that the elements of the
1.218 - //array indexes are sorted
1.219 - bool found = false;
1.220 - for (int i = 1; i <= length; ++i) {
1.221 - if (indexes[i] == ix) {
1.222 - found = true;
1.223 - values[i] = value;
1.224 - break;
1.225 - }
1.226 - }
1.227 - if (!found) {
1.228 - ++length;
1.229 - indexes[length] = ix;
1.230 - values[length] = value;
1.231 - }
1.232 -
1.233 - glp_set_mat_col(lp, jx, length, &indexes.front(), &values.front());
1.234 - }
1.235 -
1.236 - }
1.237 -
1.238 - GlpkBase::Value GlpkBase::_getCoeff(int ix, int jx) const {
1.239 -
1.240 - int length = glp_get_mat_row(lp, ix, 0, 0);
1.241 -
1.242 - std::vector<int> indexes(length + 1);
1.243 - std::vector<Value> values(length + 1);
1.244 -
1.245 - glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
1.246 -
1.247 - for (int i = 1; i <= length; ++i) {
1.248 - if (indexes[i] == jx) {
1.249 - return values[i];
1.250 - }
1.251 - }
1.252 -
1.253 - return 0;
1.254 - }
1.255 -
1.256 - void GlpkBase::_setColLowerBound(int i, Value lo) {
1.257 - LEMON_ASSERT(lo != INF, "Invalid bound");
1.258 -
1.259 - int b = glp_get_col_type(lp, i);
1.260 - double up = glp_get_col_ub(lp, i);
1.261 - if (lo == -INF) {
1.262 - switch (b) {
1.263 - case GLP_FR:
1.264 - case GLP_LO:
1.265 - glp_set_col_bnds(lp, i, GLP_FR, lo, up);
1.266 - break;
1.267 - case GLP_UP:
1.268 - break;
1.269 - case GLP_DB:
1.270 - case GLP_FX:
1.271 - glp_set_col_bnds(lp, i, GLP_UP, lo, up);
1.272 - break;
1.273 - default:
1.274 - break;
1.275 - }
1.276 - } else {
1.277 - switch (b) {
1.278 - case GLP_FR:
1.279 - case GLP_LO:
1.280 - glp_set_col_bnds(lp, i, GLP_LO, lo, up);
1.281 - break;
1.282 - case GLP_UP:
1.283 - case GLP_DB:
1.284 - case GLP_FX:
1.285 - if (lo == up)
1.286 - glp_set_col_bnds(lp, i, GLP_FX, lo, up);
1.287 - else
1.288 - glp_set_col_bnds(lp, i, GLP_DB, lo, up);
1.289 - break;
1.290 - default:
1.291 - break;
1.292 - }
1.293 - }
1.294 - }
1.295 -
1.296 - GlpkBase::Value GlpkBase::_getColLowerBound(int i) const {
1.297 - int b = glp_get_col_type(lp, i);
1.298 - switch (b) {
1.299 - case GLP_LO:
1.300 - case GLP_DB:
1.301 - case GLP_FX:
1.302 - return glp_get_col_lb(lp, i);
1.303 - default:
1.304 - return -INF;
1.305 - }
1.306 - }
1.307 -
1.308 - void GlpkBase::_setColUpperBound(int i, Value up) {
1.309 - LEMON_ASSERT(up != -INF, "Invalid bound");
1.310 -
1.311 - int b = glp_get_col_type(lp, i);
1.312 - double lo = glp_get_col_lb(lp, i);
1.313 - if (up == INF) {
1.314 - switch (b) {
1.315 - case GLP_FR:
1.316 - case GLP_LO:
1.317 - break;
1.318 - case GLP_UP:
1.319 - glp_set_col_bnds(lp, i, GLP_FR, lo, up);
1.320 - break;
1.321 - case GLP_DB:
1.322 - case GLP_FX:
1.323 - glp_set_col_bnds(lp, i, GLP_LO, lo, up);
1.324 - break;
1.325 - default:
1.326 - break;
1.327 - }
1.328 - } else {
1.329 - switch (b) {
1.330 - case GLP_FR:
1.331 - glp_set_col_bnds(lp, i, GLP_UP, lo, up);
1.332 - break;
1.333 - case GLP_UP:
1.334 - glp_set_col_bnds(lp, i, GLP_UP, lo, up);
1.335 - break;
1.336 - case GLP_LO:
1.337 - case GLP_DB:
1.338 - case GLP_FX:
1.339 - if (lo == up)
1.340 - glp_set_col_bnds(lp, i, GLP_FX, lo, up);
1.341 - else
1.342 - glp_set_col_bnds(lp, i, GLP_DB, lo, up);
1.343 - break;
1.344 - default:
1.345 - break;
1.346 - }
1.347 - }
1.348 -
1.349 - }
1.350 -
1.351 - GlpkBase::Value GlpkBase::_getColUpperBound(int i) const {
1.352 - int b = glp_get_col_type(lp, i);
1.353 - switch (b) {
1.354 - case GLP_UP:
1.355 - case GLP_DB:
1.356 - case GLP_FX:
1.357 - return glp_get_col_ub(lp, i);
1.358 - default:
1.359 - return INF;
1.360 - }
1.361 - }
1.362 -
1.363 - void GlpkBase::_setRowLowerBound(int i, Value lo) {
1.364 - LEMON_ASSERT(lo != INF, "Invalid bound");
1.365 -
1.366 - int b = glp_get_row_type(lp, i);
1.367 - double up = glp_get_row_ub(lp, i);
1.368 - if (lo == -INF) {
1.369 - switch (b) {
1.370 - case GLP_FR:
1.371 - case GLP_LO:
1.372 - glp_set_row_bnds(lp, i, GLP_FR, lo, up);
1.373 - break;
1.374 - case GLP_UP:
1.375 - break;
1.376 - case GLP_DB:
1.377 - case GLP_FX:
1.378 - glp_set_row_bnds(lp, i, GLP_UP, lo, up);
1.379 - break;
1.380 - default:
1.381 - break;
1.382 - }
1.383 - } else {
1.384 - switch (b) {
1.385 - case GLP_FR:
1.386 - case GLP_LO:
1.387 - glp_set_row_bnds(lp, i, GLP_LO, lo, up);
1.388 - break;
1.389 - case GLP_UP:
1.390 - case GLP_DB:
1.391 - case GLP_FX:
1.392 - if (lo == up)
1.393 - glp_set_row_bnds(lp, i, GLP_FX, lo, up);
1.394 - else
1.395 - glp_set_row_bnds(lp, i, GLP_DB, lo, up);
1.396 - break;
1.397 - default:
1.398 - break;
1.399 - }
1.400 - }
1.401 -
1.402 - }
1.403 -
1.404 - GlpkBase::Value GlpkBase::_getRowLowerBound(int i) const {
1.405 - int b = glp_get_row_type(lp, i);
1.406 - switch (b) {
1.407 - case GLP_LO:
1.408 - case GLP_DB:
1.409 - case GLP_FX:
1.410 - return glp_get_row_lb(lp, i);
1.411 - default:
1.412 - return -INF;
1.413 - }
1.414 - }
1.415 -
1.416 - void GlpkBase::_setRowUpperBound(int i, Value up) {
1.417 - LEMON_ASSERT(up != -INF, "Invalid bound");
1.418 -
1.419 - int b = glp_get_row_type(lp, i);
1.420 - double lo = glp_get_row_lb(lp, i);
1.421 - if (up == INF) {
1.422 - switch (b) {
1.423 - case GLP_FR:
1.424 - case GLP_LO:
1.425 - break;
1.426 - case GLP_UP:
1.427 - glp_set_row_bnds(lp, i, GLP_FR, lo, up);
1.428 - break;
1.429 - case GLP_DB:
1.430 - case GLP_FX:
1.431 - glp_set_row_bnds(lp, i, GLP_LO, lo, up);
1.432 - break;
1.433 - default:
1.434 - break;
1.435 - }
1.436 - } else {
1.437 - switch (b) {
1.438 - case GLP_FR:
1.439 - glp_set_row_bnds(lp, i, GLP_UP, lo, up);
1.440 - break;
1.441 - case GLP_UP:
1.442 - glp_set_row_bnds(lp, i, GLP_UP, lo, up);
1.443 - break;
1.444 - case GLP_LO:
1.445 - case GLP_DB:
1.446 - case GLP_FX:
1.447 - if (lo == up)
1.448 - glp_set_row_bnds(lp, i, GLP_FX, lo, up);
1.449 - else
1.450 - glp_set_row_bnds(lp, i, GLP_DB, lo, up);
1.451 - break;
1.452 - default:
1.453 - break;
1.454 - }
1.455 - }
1.456 - }
1.457 -
1.458 - GlpkBase::Value GlpkBase::_getRowUpperBound(int i) const {
1.459 - int b = glp_get_row_type(lp, i);
1.460 - switch (b) {
1.461 - case GLP_UP:
1.462 - case GLP_DB:
1.463 - case GLP_FX:
1.464 - return glp_get_row_ub(lp, i);
1.465 - default:
1.466 - return INF;
1.467 - }
1.468 - }
1.469 -
1.470 - void GlpkBase::_setObjCoeffs(ExprIterator b, ExprIterator e) {
1.471 - for (int i = 1; i <= glp_get_num_cols(lp); ++i) {
1.472 - glp_set_obj_coef(lp, i, 0.0);
1.473 - }
1.474 - for (ExprIterator it = b; it != e; ++it) {
1.475 - glp_set_obj_coef(lp, it->first, it->second);
1.476 - }
1.477 - }
1.478 -
1.479 - void GlpkBase::_getObjCoeffs(InsertIterator b) const {
1.480 - for (int i = 1; i <= glp_get_num_cols(lp); ++i) {
1.481 - Value val = glp_get_obj_coef(lp, i);
1.482 - if (val != 0.0) {
1.483 - *b = std::make_pair(i, val);
1.484 - ++b;
1.485 - }
1.486 - }
1.487 - }
1.488 -
1.489 - void GlpkBase::_setObjCoeff(int i, Value obj_coef) {
1.490 - //i = 0 means the constant term (shift)
1.491 - glp_set_obj_coef(lp, i, obj_coef);
1.492 - }
1.493 -
1.494 - GlpkBase::Value GlpkBase::_getObjCoeff(int i) const {
1.495 - //i = 0 means the constant term (shift)
1.496 - return glp_get_obj_coef(lp, i);
1.497 - }
1.498 -
1.499 - void GlpkBase::_setSense(GlpkBase::Sense sense) {
1.500 - switch (sense) {
1.501 - case MIN:
1.502 - glp_set_obj_dir(lp, GLP_MIN);
1.503 - break;
1.504 - case MAX:
1.505 - glp_set_obj_dir(lp, GLP_MAX);
1.506 - break;
1.507 - }
1.508 - }
1.509 -
1.510 - GlpkBase::Sense GlpkBase::_getSense() const {
1.511 - switch(glp_get_obj_dir(lp)) {
1.512 - case GLP_MIN:
1.513 - return MIN;
1.514 - case GLP_MAX:
1.515 - return MAX;
1.516 - default:
1.517 - LEMON_ASSERT(false, "Wrong sense");
1.518 - return GlpkBase::Sense();
1.519 - }
1.520 - }
1.521 -
1.522 - void GlpkBase::_clear() {
1.523 - glp_erase_prob(lp);
1.524 - rows.clear();
1.525 - cols.clear();
1.526 - }
1.527 -
1.528 - // LpGlpk members
1.529 -
1.530 - LpGlpk::LpGlpk()
1.531 - : LpBase(), GlpkBase(), LpSolver() {
1.532 - messageLevel(MESSAGE_NO_OUTPUT);
1.533 - }
1.534 -
1.535 - LpGlpk::LpGlpk(const LpGlpk& other)
1.536 - : LpBase(other), GlpkBase(other), LpSolver(other) {
1.537 - messageLevel(MESSAGE_NO_OUTPUT);
1.538 - }
1.539 -
1.540 - LpGlpk* LpGlpk::_newSolver() const { return new LpGlpk; }
1.541 - LpGlpk* LpGlpk::_cloneSolver() const { return new LpGlpk(*this); }
1.542 -
1.543 - const char* LpGlpk::_solverName() const { return "LpGlpk"; }
1.544 -
1.545 - void LpGlpk::_clear_temporals() {
1.546 - _primal_ray.clear();
1.547 - _dual_ray.clear();
1.548 - }
1.549 -
1.550 - LpGlpk::SolveExitStatus LpGlpk::_solve() {
1.551 - return solvePrimal();
1.552 - }
1.553 -
1.554 - LpGlpk::SolveExitStatus LpGlpk::solvePrimal() {
1.555 - _clear_temporals();
1.556 -
1.557 - glp_smcp smcp;
1.558 - glp_init_smcp(&smcp);
1.559 -
1.560 - switch (_message_level) {
1.561 - case MESSAGE_NO_OUTPUT:
1.562 - smcp.msg_lev = GLP_MSG_OFF;
1.563 - break;
1.564 - case MESSAGE_ERROR_MESSAGE:
1.565 - smcp.msg_lev = GLP_MSG_ERR;
1.566 - break;
1.567 - case MESSAGE_NORMAL_OUTPUT:
1.568 - smcp.msg_lev = GLP_MSG_ON;
1.569 - break;
1.570 - case MESSAGE_FULL_OUTPUT:
1.571 - smcp.msg_lev = GLP_MSG_ALL;
1.572 - break;
1.573 - }
1.574 -
1.575 - if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
1.576 - return SOLVED;
1.577 - }
1.578 -
1.579 - LpGlpk::SolveExitStatus LpGlpk::solveDual() {
1.580 - _clear_temporals();
1.581 -
1.582 - glp_smcp smcp;
1.583 - glp_init_smcp(&smcp);
1.584 -
1.585 - switch (_message_level) {
1.586 - case MESSAGE_NO_OUTPUT:
1.587 - smcp.msg_lev = GLP_MSG_OFF;
1.588 - break;
1.589 - case MESSAGE_ERROR_MESSAGE:
1.590 - smcp.msg_lev = GLP_MSG_ERR;
1.591 - break;
1.592 - case MESSAGE_NORMAL_OUTPUT:
1.593 - smcp.msg_lev = GLP_MSG_ON;
1.594 - break;
1.595 - case MESSAGE_FULL_OUTPUT:
1.596 - smcp.msg_lev = GLP_MSG_ALL;
1.597 - break;
1.598 - }
1.599 - smcp.meth = GLP_DUAL;
1.600 -
1.601 - if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
1.602 - return SOLVED;
1.603 - }
1.604 -
1.605 - LpGlpk::Value LpGlpk::_getPrimal(int i) const {
1.606 - return glp_get_col_prim(lp, i);
1.607 - }
1.608 -
1.609 - LpGlpk::Value LpGlpk::_getDual(int i) const {
1.610 - return glp_get_row_dual(lp, i);
1.611 - }
1.612 -
1.613 - LpGlpk::Value LpGlpk::_getPrimalValue() const {
1.614 - return glp_get_obj_val(lp);
1.615 - }
1.616 -
1.617 - LpGlpk::VarStatus LpGlpk::_getColStatus(int i) const {
1.618 - switch (glp_get_col_stat(lp, i)) {
1.619 - case GLP_BS:
1.620 - return BASIC;
1.621 - case GLP_UP:
1.622 - return UPPER;
1.623 - case GLP_LO:
1.624 - return LOWER;
1.625 - case GLP_NF:
1.626 - return FREE;
1.627 - case GLP_NS:
1.628 - return FIXED;
1.629 - default:
1.630 - LEMON_ASSERT(false, "Wrong column status");
1.631 - return LpGlpk::VarStatus();
1.632 - }
1.633 - }
1.634 -
1.635 - LpGlpk::VarStatus LpGlpk::_getRowStatus(int i) const {
1.636 - switch (glp_get_row_stat(lp, i)) {
1.637 - case GLP_BS:
1.638 - return BASIC;
1.639 - case GLP_UP:
1.640 - return UPPER;
1.641 - case GLP_LO:
1.642 - return LOWER;
1.643 - case GLP_NF:
1.644 - return FREE;
1.645 - case GLP_NS:
1.646 - return FIXED;
1.647 - default:
1.648 - LEMON_ASSERT(false, "Wrong row status");
1.649 - return LpGlpk::VarStatus();
1.650 - }
1.651 - }
1.652 -
1.653 - LpGlpk::Value LpGlpk::_getPrimalRay(int i) const {
1.654 - if (_primal_ray.empty()) {
1.655 - int row_num = glp_get_num_rows(lp);
1.656 - int col_num = glp_get_num_cols(lp);
1.657 -
1.658 - _primal_ray.resize(col_num + 1, 0.0);
1.659 -
1.660 - int index = glp_get_unbnd_ray(lp);
1.661 - if (index != 0) {
1.662 - // The primal ray is found in primal simplex second phase
1.663 - LEMON_ASSERT((index <= row_num ? glp_get_row_stat(lp, index) :
1.664 - glp_get_col_stat(lp, index - row_num)) != GLP_BS,
1.665 - "Wrong primal ray");
1.666 -
1.667 - bool negate = glp_get_obj_dir(lp) == GLP_MAX;
1.668 -
1.669 - if (index > row_num) {
1.670 - _primal_ray[index - row_num] = 1.0;
1.671 - if (glp_get_col_dual(lp, index - row_num) > 0) {
1.672 - negate = !negate;
1.673 - }
1.674 - } else {
1.675 - if (glp_get_row_dual(lp, index) > 0) {
1.676 - negate = !negate;
1.677 - }
1.678 - }
1.679 -
1.680 - std::vector<int> ray_indexes(row_num + 1);
1.681 - std::vector<Value> ray_values(row_num + 1);
1.682 - int ray_length = glp_eval_tab_col(lp, index, &ray_indexes.front(),
1.683 - &ray_values.front());
1.684 -
1.685 - for (int i = 1; i <= ray_length; ++i) {
1.686 - if (ray_indexes[i] > row_num) {
1.687 - _primal_ray[ray_indexes[i] - row_num] = ray_values[i];
1.688 - }
1.689 - }
1.690 -
1.691 - if (negate) {
1.692 - for (int i = 1; i <= col_num; ++i) {
1.693 - _primal_ray[i] = - _primal_ray[i];
1.694 - }
1.695 - }
1.696 - } else {
1.697 - for (int i = 1; i <= col_num; ++i) {
1.698 - _primal_ray[i] = glp_get_col_prim(lp, i);
1.699 - }
1.700 - }
1.701 - }
1.702 - return _primal_ray[i];
1.703 - }
1.704 -
1.705 - LpGlpk::Value LpGlpk::_getDualRay(int i) const {
1.706 - if (_dual_ray.empty()) {
1.707 - int row_num = glp_get_num_rows(lp);
1.708 -
1.709 - _dual_ray.resize(row_num + 1, 0.0);
1.710 -
1.711 - int index = glp_get_unbnd_ray(lp);
1.712 - if (index != 0) {
1.713 - // The dual ray is found in dual simplex second phase
1.714 - LEMON_ASSERT((index <= row_num ? glp_get_row_stat(lp, index) :
1.715 - glp_get_col_stat(lp, index - row_num)) == GLP_BS,
1.716 -
1.717 - "Wrong dual ray");
1.718 -
1.719 - int idx;
1.720 - bool negate = false;
1.721 -
1.722 - if (index > row_num) {
1.723 - idx = glp_get_col_bind(lp, index - row_num);
1.724 - if (glp_get_col_prim(lp, index - row_num) >
1.725 - glp_get_col_ub(lp, index - row_num)) {
1.726 - negate = true;
1.727 - }
1.728 - } else {
1.729 - idx = glp_get_row_bind(lp, index);
1.730 - if (glp_get_row_prim(lp, index) > glp_get_row_ub(lp, index)) {
1.731 - negate = true;
1.732 - }
1.733 - }
1.734 -
1.735 - _dual_ray[idx] = negate ? - 1.0 : 1.0;
1.736 -
1.737 - glp_btran(lp, &_dual_ray.front());
1.738 - } else {
1.739 - double eps = 1e-7;
1.740 - // The dual ray is found in primal simplex first phase
1.741 - // We assume that the glpk minimizes the slack to get feasible solution
1.742 - for (int i = 1; i <= row_num; ++i) {
1.743 - int index = glp_get_bhead(lp, i);
1.744 - if (index <= row_num) {
1.745 - double res = glp_get_row_prim(lp, index);
1.746 - if (res > glp_get_row_ub(lp, index) + eps) {
1.747 - _dual_ray[i] = -1;
1.748 - } else if (res < glp_get_row_lb(lp, index) - eps) {
1.749 - _dual_ray[i] = 1;
1.750 - } else {
1.751 - _dual_ray[i] = 0;
1.752 - }
1.753 - _dual_ray[i] *= glp_get_rii(lp, index);
1.754 - } else {
1.755 - double res = glp_get_col_prim(lp, index - row_num);
1.756 - if (res > glp_get_col_ub(lp, index - row_num) + eps) {
1.757 - _dual_ray[i] = -1;
1.758 - } else if (res < glp_get_col_lb(lp, index - row_num) - eps) {
1.759 - _dual_ray[i] = 1;
1.760 - } else {
1.761 - _dual_ray[i] = 0;
1.762 - }
1.763 - _dual_ray[i] /= glp_get_sjj(lp, index - row_num);
1.764 - }
1.765 - }
1.766 -
1.767 - glp_btran(lp, &_dual_ray.front());
1.768 -
1.769 - for (int i = 1; i <= row_num; ++i) {
1.770 - _dual_ray[i] /= glp_get_rii(lp, i);
1.771 - }
1.772 - }
1.773 - }
1.774 - return _dual_ray[i];
1.775 - }
1.776 -
1.777 - LpGlpk::ProblemType LpGlpk::_getPrimalType() const {
1.778 - if (glp_get_status(lp) == GLP_OPT)
1.779 - return OPTIMAL;
1.780 - switch (glp_get_prim_stat(lp)) {
1.781 - case GLP_UNDEF:
1.782 - return UNDEFINED;
1.783 - case GLP_FEAS:
1.784 - case GLP_INFEAS:
1.785 - if (glp_get_dual_stat(lp) == GLP_NOFEAS) {
1.786 - return UNBOUNDED;
1.787 - } else {
1.788 - return UNDEFINED;
1.789 - }
1.790 - case GLP_NOFEAS:
1.791 - return INFEASIBLE;
1.792 - default:
1.793 - LEMON_ASSERT(false, "Wrong primal type");
1.794 - return LpGlpk::ProblemType();
1.795 - }
1.796 - }
1.797 -
1.798 - LpGlpk::ProblemType LpGlpk::_getDualType() const {
1.799 - if (glp_get_status(lp) == GLP_OPT)
1.800 - return OPTIMAL;
1.801 - switch (glp_get_dual_stat(lp)) {
1.802 - case GLP_UNDEF:
1.803 - return UNDEFINED;
1.804 - case GLP_FEAS:
1.805 - case GLP_INFEAS:
1.806 - if (glp_get_prim_stat(lp) == GLP_NOFEAS) {
1.807 - return UNBOUNDED;
1.808 - } else {
1.809 - return UNDEFINED;
1.810 - }
1.811 - case GLP_NOFEAS:
1.812 - return INFEASIBLE;
1.813 - default:
1.814 - LEMON_ASSERT(false, "Wrong primal type");
1.815 - return LpGlpk::ProblemType();
1.816 - }
1.817 - }
1.818 -
1.819 - void LpGlpk::presolver(bool b) {
1.820 - lpx_set_int_parm(lp, LPX_K_PRESOL, b ? 1 : 0);
1.821 - }
1.822 -
1.823 - void LpGlpk::messageLevel(MessageLevel m) {
1.824 - _message_level = m;
1.825 - }
1.826 -
1.827 - // MipGlpk members
1.828 -
1.829 - MipGlpk::MipGlpk()
1.830 - : LpBase(), GlpkBase(), MipSolver() {
1.831 - messageLevel(MESSAGE_NO_OUTPUT);
1.832 - }
1.833 -
1.834 - MipGlpk::MipGlpk(const MipGlpk& other)
1.835 - : LpBase(), GlpkBase(other), MipSolver() {
1.836 - messageLevel(MESSAGE_NO_OUTPUT);
1.837 - }
1.838 -
1.839 - void MipGlpk::_setColType(int i, MipGlpk::ColTypes col_type) {
1.840 - switch (col_type) {
1.841 - case INTEGER:
1.842 - glp_set_col_kind(lp, i, GLP_IV);
1.843 - break;
1.844 - case REAL:
1.845 - glp_set_col_kind(lp, i, GLP_CV);
1.846 - break;
1.847 - }
1.848 - }
1.849 -
1.850 - MipGlpk::ColTypes MipGlpk::_getColType(int i) const {
1.851 - switch (glp_get_col_kind(lp, i)) {
1.852 - case GLP_IV:
1.853 - case GLP_BV:
1.854 - return INTEGER;
1.855 - default:
1.856 - return REAL;
1.857 - }
1.858 -
1.859 - }
1.860 -
1.861 - MipGlpk::SolveExitStatus MipGlpk::_solve() {
1.862 - glp_smcp smcp;
1.863 - glp_init_smcp(&smcp);
1.864 -
1.865 - switch (_message_level) {
1.866 - case MESSAGE_NO_OUTPUT:
1.867 - smcp.msg_lev = GLP_MSG_OFF;
1.868 - break;
1.869 - case MESSAGE_ERROR_MESSAGE:
1.870 - smcp.msg_lev = GLP_MSG_ERR;
1.871 - break;
1.872 - case MESSAGE_NORMAL_OUTPUT:
1.873 - smcp.msg_lev = GLP_MSG_ON;
1.874 - break;
1.875 - case MESSAGE_FULL_OUTPUT:
1.876 - smcp.msg_lev = GLP_MSG_ALL;
1.877 - break;
1.878 - }
1.879 - smcp.meth = GLP_DUAL;
1.880 -
1.881 - if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
1.882 - if (glp_get_status(lp) != GLP_OPT) return SOLVED;
1.883 -
1.884 - glp_iocp iocp;
1.885 - glp_init_iocp(&iocp);
1.886 -
1.887 - switch (_message_level) {
1.888 - case MESSAGE_NO_OUTPUT:
1.889 - iocp.msg_lev = GLP_MSG_OFF;
1.890 - break;
1.891 - case MESSAGE_ERROR_MESSAGE:
1.892 - iocp.msg_lev = GLP_MSG_ERR;
1.893 - break;
1.894 - case MESSAGE_NORMAL_OUTPUT:
1.895 - iocp.msg_lev = GLP_MSG_ON;
1.896 - break;
1.897 - case MESSAGE_FULL_OUTPUT:
1.898 - iocp.msg_lev = GLP_MSG_ALL;
1.899 - break;
1.900 - }
1.901 -
1.902 - if (glp_intopt(lp, &iocp) != 0) return UNSOLVED;
1.903 - return SOLVED;
1.904 - }
1.905 -
1.906 -
1.907 - MipGlpk::ProblemType MipGlpk::_getType() const {
1.908 - switch (glp_get_status(lp)) {
1.909 - case GLP_OPT:
1.910 - switch (glp_mip_status(lp)) {
1.911 - case GLP_UNDEF:
1.912 - return UNDEFINED;
1.913 - case GLP_NOFEAS:
1.914 - return INFEASIBLE;
1.915 - case GLP_FEAS:
1.916 - return FEASIBLE;
1.917 - case GLP_OPT:
1.918 - return OPTIMAL;
1.919 - default:
1.920 - LEMON_ASSERT(false, "Wrong problem type.");
1.921 - return MipGlpk::ProblemType();
1.922 - }
1.923 - case GLP_NOFEAS:
1.924 - return INFEASIBLE;
1.925 - case GLP_INFEAS:
1.926 - case GLP_FEAS:
1.927 - if (glp_get_dual_stat(lp) == GLP_NOFEAS) {
1.928 - return UNBOUNDED;
1.929 - } else {
1.930 - return UNDEFINED;
1.931 - }
1.932 - default:
1.933 - LEMON_ASSERT(false, "Wrong problem type.");
1.934 - return MipGlpk::ProblemType();
1.935 - }
1.936 - }
1.937 -
1.938 - MipGlpk::Value MipGlpk::_getSol(int i) const {
1.939 - return glp_mip_col_val(lp, i);
1.940 - }
1.941 -
1.942 - MipGlpk::Value MipGlpk::_getSolValue() const {
1.943 - return glp_mip_obj_val(lp);
1.944 - }
1.945 -
1.946 - MipGlpk* MipGlpk::_newSolver() const { return new MipGlpk; }
1.947 - MipGlpk* MipGlpk::_cloneSolver() const {return new MipGlpk(*this); }
1.948 -
1.949 - const char* MipGlpk::_solverName() const { return "MipGlpk"; }
1.950 -
1.951 - void MipGlpk::messageLevel(MessageLevel m) {
1.952 - _message_level = m;
1.953 - }
1.954 -
1.955 -} //END OF NAMESPACE LEMON