lemon/soplex.cc
author Peter Kovacs <kpeter@inf.elte.hu>
Thu, 12 Nov 2009 23:26:13 +0100
changeset 806 fa6f37d7a25b
parent 576 745e182d0139
child 877 141f9c0db4a3
permissions -rw-r--r--
Entirely rework CapacityScaling (#180)

- Use the new interface similarly to NetworkSimplex.
- Rework the implementation using an efficient internal structure
for handling the residual network. This improvement made the
code much faster (up to 2-5 times faster on large graphs).
- Handle GEQ supply type (LEQ is not supported).
- Handle negative costs for arcs of finite capacity.
(Note that this algorithm cannot handle arcs of negative cost
and infinite upper bound, thus it returns UNBOUNDED if such
an arc exists.)
- Extend the documentation.
     1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
     2  *
     3  * This file is a part of LEMON, a generic C++ optimization library.
     4  *
     5  * Copyright (C) 2003-2008
     6  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7  * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8  *
     9  * Permission to use, modify and distribute this software is granted
    10  * provided that this copyright notice appears in all copies. For
    11  * precise terms see the accompanying LICENSE file.
    12  *
    13  * This software is provided "AS IS" with no warranty of any kind,
    14  * express or implied, and with no claim as to its suitability for any
    15  * purpose.
    16  *
    17  */
    18 
    19 #include <iostream>
    20 #include <lemon/soplex.h>
    21 
    22 #include <soplex.h>
    23 #include <spxout.h>
    24 
    25 
    26 ///\file
    27 ///\brief Implementation of the LEMON-SOPLEX lp solver interface.
    28 namespace lemon {
    29 
    30   SoplexLp::SoplexLp() {
    31     soplex = new soplex::SoPlex;
    32     messageLevel(MESSAGE_NOTHING);
    33   }
    34 
    35   SoplexLp::~SoplexLp() {
    36     delete soplex;
    37   }
    38 
    39   SoplexLp::SoplexLp(const SoplexLp& lp) {
    40     rows = lp.rows;
    41     cols = lp.cols;
    42 
    43     soplex = new soplex::SoPlex;
    44     (*static_cast<soplex::SPxLP*>(soplex)) = *(lp.soplex);
    45 
    46     _col_names = lp._col_names;
    47     _col_names_ref = lp._col_names_ref;
    48 
    49     _row_names = lp._row_names;
    50     _row_names_ref = lp._row_names_ref;
    51 
    52     messageLevel(MESSAGE_NOTHING);
    53   }
    54 
    55   void SoplexLp::_clear_temporals() {
    56     _primal_values.clear();
    57     _dual_values.clear();
    58   }
    59 
    60   SoplexLp* SoplexLp::newSolver() const {
    61     SoplexLp* newlp = new SoplexLp();
    62     return newlp;
    63   }
    64 
    65   SoplexLp* SoplexLp::cloneSolver() const {
    66     SoplexLp* newlp = new SoplexLp(*this);
    67     return newlp;
    68   }
    69 
    70   const char* SoplexLp::_solverName() const { return "SoplexLp"; }
    71 
    72   int SoplexLp::_addCol() {
    73     soplex::LPCol c;
    74     c.setLower(-soplex::infinity);
    75     c.setUpper(soplex::infinity);
    76     soplex->addCol(c);
    77 
    78     _col_names.push_back(std::string());
    79 
    80     return soplex->nCols() - 1;
    81   }
    82 
    83   int SoplexLp::_addRow() {
    84     soplex::LPRow r;
    85     r.setLhs(-soplex::infinity);
    86     r.setRhs(soplex::infinity);
    87     soplex->addRow(r);
    88 
    89     _row_names.push_back(std::string());
    90 
    91     return soplex->nRows() - 1;
    92   }
    93 
    94   int SoplexLp::_addRow(Value l, ExprIterator b, ExprIterator e, Value u) {
    95     soplex::DSVector v;
    96     for (ExprIterator it = b; it != e; ++it) {
    97       v.add(it->first, it->second);
    98     }
    99     soplex::LPRow r(l, v, u);
   100     soplex->addRow(r);
   101 
   102     _row_names.push_back(std::string());
   103 
   104     return soplex->nRows() - 1;
   105   }
   106 
   107 
   108   void SoplexLp::_eraseCol(int i) {
   109     soplex->removeCol(i);
   110     _col_names_ref.erase(_col_names[i]);
   111     _col_names[i] = _col_names.back();
   112     _col_names_ref[_col_names.back()] = i;
   113     _col_names.pop_back();
   114   }
   115 
   116   void SoplexLp::_eraseRow(int i) {
   117     soplex->removeRow(i);
   118     _row_names_ref.erase(_row_names[i]);
   119     _row_names[i] = _row_names.back();
   120     _row_names_ref[_row_names.back()] = i;
   121     _row_names.pop_back();
   122   }
   123 
   124   void SoplexLp::_eraseColId(int i) {
   125     cols.eraseIndex(i);
   126     cols.relocateIndex(i, cols.maxIndex());
   127   }
   128   void SoplexLp::_eraseRowId(int i) {
   129     rows.eraseIndex(i);
   130     rows.relocateIndex(i, rows.maxIndex());
   131   }
   132 
   133   void SoplexLp::_getColName(int c, std::string &name) const {
   134     name = _col_names[c];
   135   }
   136 
   137   void SoplexLp::_setColName(int c, const std::string &name) {
   138     _col_names_ref.erase(_col_names[c]);
   139     _col_names[c] = name;
   140     if (!name.empty()) {
   141       _col_names_ref.insert(std::make_pair(name, c));
   142     }
   143   }
   144 
   145   int SoplexLp::_colByName(const std::string& name) const {
   146     std::map<std::string, int>::const_iterator it =
   147       _col_names_ref.find(name);
   148     if (it != _col_names_ref.end()) {
   149       return it->second;
   150     } else {
   151       return -1;
   152     }
   153   }
   154 
   155   void SoplexLp::_getRowName(int r, std::string &name) const {
   156     name = _row_names[r];
   157   }
   158 
   159   void SoplexLp::_setRowName(int r, const std::string &name) {
   160     _row_names_ref.erase(_row_names[r]);
   161     _row_names[r] = name;
   162     if (!name.empty()) {
   163       _row_names_ref.insert(std::make_pair(name, r));
   164     }
   165   }
   166 
   167   int SoplexLp::_rowByName(const std::string& name) const {
   168     std::map<std::string, int>::const_iterator it =
   169       _row_names_ref.find(name);
   170     if (it != _row_names_ref.end()) {
   171       return it->second;
   172     } else {
   173       return -1;
   174     }
   175   }
   176 
   177 
   178   void SoplexLp::_setRowCoeffs(int i, ExprIterator b, ExprIterator e) {
   179     for (int j = 0; j < soplex->nCols(); ++j) {
   180       soplex->changeElement(i, j, 0.0);
   181     }
   182     for(ExprIterator it = b; it != e; ++it) {
   183       soplex->changeElement(i, it->first, it->second);
   184     }
   185   }
   186 
   187   void SoplexLp::_getRowCoeffs(int i, InsertIterator b) const {
   188     const soplex::SVector& vec = soplex->rowVector(i);
   189     for (int k = 0; k < vec.size(); ++k) {
   190       *b = std::make_pair(vec.index(k), vec.value(k));
   191       ++b;
   192     }
   193   }
   194 
   195   void SoplexLp::_setColCoeffs(int j, ExprIterator b, ExprIterator e) {
   196     for (int i = 0; i < soplex->nRows(); ++i) {
   197       soplex->changeElement(i, j, 0.0);
   198     }
   199     for(ExprIterator it = b; it != e; ++it) {
   200       soplex->changeElement(it->first, j, it->second);
   201     }
   202   }
   203 
   204   void SoplexLp::_getColCoeffs(int i, InsertIterator b) const {
   205     const soplex::SVector& vec = soplex->colVector(i);
   206     for (int k = 0; k < vec.size(); ++k) {
   207       *b = std::make_pair(vec.index(k), vec.value(k));
   208       ++b;
   209     }
   210   }
   211 
   212   void SoplexLp::_setCoeff(int i, int j, Value value) {
   213     soplex->changeElement(i, j, value);
   214   }
   215 
   216   SoplexLp::Value SoplexLp::_getCoeff(int i, int j) const {
   217     return soplex->rowVector(i)[j];
   218   }
   219 
   220   void SoplexLp::_setColLowerBound(int i, Value value) {
   221     LEMON_ASSERT(value != INF, "Invalid bound");
   222     soplex->changeLower(i, value != -INF ? value : -soplex::infinity);
   223   }
   224 
   225   SoplexLp::Value SoplexLp::_getColLowerBound(int i) const {
   226     double value = soplex->lower(i);
   227     return value != -soplex::infinity ? value : -INF;
   228   }
   229 
   230   void SoplexLp::_setColUpperBound(int i, Value value) {
   231     LEMON_ASSERT(value != -INF, "Invalid bound");
   232     soplex->changeUpper(i, value != INF ? value : soplex::infinity);
   233   }
   234 
   235   SoplexLp::Value SoplexLp::_getColUpperBound(int i) const {
   236     double value = soplex->upper(i);
   237     return value != soplex::infinity ? value : INF;
   238   }
   239 
   240   void SoplexLp::_setRowLowerBound(int i, Value lb) {
   241     LEMON_ASSERT(lb != INF, "Invalid bound");
   242     soplex->changeRange(i, lb != -INF ? lb : -soplex::infinity, soplex->rhs(i));
   243   }
   244 
   245   SoplexLp::Value SoplexLp::_getRowLowerBound(int i) const {
   246     double res = soplex->lhs(i);
   247     return res == -soplex::infinity ? -INF : res;
   248   }
   249 
   250   void SoplexLp::_setRowUpperBound(int i, Value ub) {
   251     LEMON_ASSERT(ub != -INF, "Invalid bound");
   252     soplex->changeRange(i, soplex->lhs(i), ub != INF ? ub : soplex::infinity);
   253   }
   254 
   255   SoplexLp::Value SoplexLp::_getRowUpperBound(int i) const {
   256     double res = soplex->rhs(i);
   257     return res == soplex::infinity ? INF : res;
   258   }
   259 
   260   void SoplexLp::_setObjCoeffs(ExprIterator b, ExprIterator e) {
   261     for (int j = 0; j < soplex->nCols(); ++j) {
   262       soplex->changeObj(j, 0.0);
   263     }
   264     for (ExprIterator it = b; it != e; ++it) {
   265       soplex->changeObj(it->first, it->second);
   266     }
   267   }
   268 
   269   void SoplexLp::_getObjCoeffs(InsertIterator b) const {
   270     for (int j = 0; j < soplex->nCols(); ++j) {
   271       Value coef = soplex->obj(j);
   272       if (coef != 0.0) {
   273         *b = std::make_pair(j, coef);
   274         ++b;
   275       }
   276     }
   277   }
   278 
   279   void SoplexLp::_setObjCoeff(int i, Value obj_coef) {
   280     soplex->changeObj(i, obj_coef);
   281   }
   282 
   283   SoplexLp::Value SoplexLp::_getObjCoeff(int i) const {
   284     return soplex->obj(i);
   285   }
   286 
   287   SoplexLp::SolveExitStatus SoplexLp::_solve() {
   288 
   289     _clear_temporals();
   290     
   291     _applyMessageLevel();
   292 
   293     soplex::SPxSolver::Status status = soplex->solve();
   294 
   295     switch (status) {
   296     case soplex::SPxSolver::OPTIMAL:
   297     case soplex::SPxSolver::INFEASIBLE:
   298     case soplex::SPxSolver::UNBOUNDED:
   299       return SOLVED;
   300     default:
   301       return UNSOLVED;
   302     }
   303   }
   304 
   305   SoplexLp::Value SoplexLp::_getPrimal(int i) const {
   306     if (_primal_values.empty()) {
   307       _primal_values.resize(soplex->nCols());
   308       soplex::Vector pv(_primal_values.size(), &_primal_values.front());
   309       soplex->getPrimal(pv);
   310     }
   311     return _primal_values[i];
   312   }
   313 
   314   SoplexLp::Value SoplexLp::_getDual(int i) const {
   315     if (_dual_values.empty()) {
   316       _dual_values.resize(soplex->nRows());
   317       soplex::Vector dv(_dual_values.size(), &_dual_values.front());
   318       soplex->getDual(dv);
   319     }
   320     return _dual_values[i];
   321   }
   322 
   323   SoplexLp::Value SoplexLp::_getPrimalValue() const {
   324     return soplex->objValue();
   325   }
   326 
   327   SoplexLp::VarStatus SoplexLp::_getColStatus(int i) const {
   328     switch (soplex->getBasisColStatus(i)) {
   329     case soplex::SPxSolver::BASIC:
   330       return BASIC;
   331     case soplex::SPxSolver::ON_UPPER:
   332       return UPPER;
   333     case soplex::SPxSolver::ON_LOWER:
   334       return LOWER;
   335     case soplex::SPxSolver::FIXED:
   336       return FIXED;
   337     case soplex::SPxSolver::ZERO:
   338       return FREE;
   339     default:
   340       LEMON_ASSERT(false, "Wrong column status");
   341       return VarStatus();
   342     }
   343   }
   344 
   345   SoplexLp::VarStatus SoplexLp::_getRowStatus(int i) const {
   346     switch (soplex->getBasisRowStatus(i)) {
   347     case soplex::SPxSolver::BASIC:
   348       return BASIC;
   349     case soplex::SPxSolver::ON_UPPER:
   350       return UPPER;
   351     case soplex::SPxSolver::ON_LOWER:
   352       return LOWER;
   353     case soplex::SPxSolver::FIXED:
   354       return FIXED;
   355     case soplex::SPxSolver::ZERO:
   356       return FREE;
   357     default:
   358       LEMON_ASSERT(false, "Wrong row status");
   359       return VarStatus();
   360     }
   361   }
   362 
   363   SoplexLp::Value SoplexLp::_getPrimalRay(int i) const {
   364     if (_primal_ray.empty()) {
   365       _primal_ray.resize(soplex->nCols());
   366       soplex::Vector pv(_primal_ray.size(), &_primal_ray.front());
   367       soplex->getDualfarkas(pv);
   368     }
   369     return _primal_ray[i];
   370   }
   371 
   372   SoplexLp::Value SoplexLp::_getDualRay(int i) const {
   373     if (_dual_ray.empty()) {
   374       _dual_ray.resize(soplex->nRows());
   375       soplex::Vector dv(_dual_ray.size(), &_dual_ray.front());
   376       soplex->getDualfarkas(dv);
   377     }
   378     return _dual_ray[i];
   379   }
   380 
   381   SoplexLp::ProblemType SoplexLp::_getPrimalType() const {
   382     switch (soplex->status()) {
   383     case soplex::SPxSolver::OPTIMAL:
   384       return OPTIMAL;
   385     case soplex::SPxSolver::UNBOUNDED:
   386       return UNBOUNDED;
   387     case soplex::SPxSolver::INFEASIBLE:
   388       return INFEASIBLE;
   389     default:
   390       return UNDEFINED;
   391     }
   392   }
   393 
   394   SoplexLp::ProblemType SoplexLp::_getDualType() const {
   395     switch (soplex->status()) {
   396     case soplex::SPxSolver::OPTIMAL:
   397       return OPTIMAL;
   398     case soplex::SPxSolver::UNBOUNDED:
   399       return UNBOUNDED;
   400     case soplex::SPxSolver::INFEASIBLE:
   401       return INFEASIBLE;
   402     default:
   403       return UNDEFINED;
   404     }
   405   }
   406 
   407   void SoplexLp::_setSense(Sense sense) {
   408     switch (sense) {
   409     case MIN:
   410       soplex->changeSense(soplex::SPxSolver::MINIMIZE);
   411       break;
   412     case MAX:
   413       soplex->changeSense(soplex::SPxSolver::MAXIMIZE);
   414     }
   415   }
   416 
   417   SoplexLp::Sense SoplexLp::_getSense() const {
   418     switch (soplex->spxSense()) {
   419     case soplex::SPxSolver::MAXIMIZE:
   420       return MAX;
   421     case soplex::SPxSolver::MINIMIZE:
   422       return MIN;
   423     default:
   424       LEMON_ASSERT(false, "Wrong sense.");
   425       return SoplexLp::Sense();
   426     }
   427   }
   428 
   429   void SoplexLp::_clear() {
   430     soplex->clear();
   431     _col_names.clear();
   432     _col_names_ref.clear();
   433     _row_names.clear();
   434     _row_names_ref.clear();
   435     cols.clear();
   436     rows.clear();
   437     _clear_temporals();
   438   }
   439 
   440   void SoplexLp::_messageLevel(MessageLevel level) {
   441     switch (level) {
   442     case MESSAGE_NOTHING:
   443       _message_level = -1;
   444       break;
   445     case MESSAGE_ERROR:
   446       _message_level = soplex::SPxOut::ERROR;
   447       break;
   448     case MESSAGE_WARNING:
   449       _message_level = soplex::SPxOut::WARNING;
   450       break;
   451     case MESSAGE_NORMAL:
   452       _message_level = soplex::SPxOut::INFO2;
   453       break;
   454     case MESSAGE_VERBOSE:
   455       _message_level = soplex::SPxOut::DEBUG;
   456       break;
   457     }
   458   }
   459 
   460   void SoplexLp::_applyMessageLevel() {
   461     soplex::Param::setVerbose(_message_level);
   462   }
   463 
   464 } //namespace lemon
   465