src/work/marci/lp/lp_solver_base.h
author hegyi
Fri, 01 Apr 2005 09:43:52 +0000
changeset 1288 6cc7b573b7b5
parent 1153 4b0468de3a31
permissions -rw-r--r--
Graph displayer is now displaying nodes. Edges remain still undisplayed yet.
     1 // -*- c++ -*-
     2 #ifndef LEMON_LP_SOLVER_BASE_H
     3 #define LEMON_LP_SOLVER_BASE_H
     4 
     5 ///\ingroup misc
     6 ///\file
     7 
     8 // #include <stdio.h>
     9 #include <stdlib.h>
    10 #include <iostream>
    11 #include <map>
    12 #include <limits>
    13 // #include <stdio>
    14 //#include <stdlib>
    15 extern "C" {
    16 #include "glpk.h"
    17 }
    18 
    19 #include <iostream>
    20 #include <vector>
    21 #include <string>
    22 #include <list>
    23 #include <memory>
    24 #include <utility>
    25 
    26 #include <lemon/invalid.h>
    27 #include <expression.h>
    28 //#include <stp.h>
    29 //#include <lemon/max_flow.h>
    30 //#include <augmenting_flow.h>
    31 //#include <iter_map.h>
    32 
    33 using std::cout;
    34 using std::cin;
    35 using std::endl;
    36 
    37 namespace lemon {
    38   
    39   /// \addtogroup misc
    40   /// @{
    41 
    42   /// \brief A partitioned vector with iterable classes.
    43   ///
    44   /// This class implements a container in which the data is stored in an 
    45   /// stl vector, the range is partitioned into sets and each set is 
    46   /// doubly linked in a list. 
    47   /// That is, each class is iterable by lemon iterators, and any member of 
    48   /// the vector can bo moved to an other class.
    49   template <typename T>
    50   class IterablePartition {
    51   protected:
    52     struct Node {
    53       T data;
    54       int prev; //invalid az -1
    55       int next; 
    56     };
    57     std::vector<Node> nodes;
    58     struct Tip {
    59       int first;
    60       int last;
    61     };
    62     std::vector<Tip> tips;
    63   public:
    64     /// The classes are indexed by integers from \c 0 to \c classNum()-1.
    65     int classNum() const { return tips.size(); }
    66     /// This lemon style iterator iterates through a class. 
    67     class Class;
    68     /// Constructor. The number of classes is to be given which is fixed 
    69     /// over the life of the container. 
    70     /// The partition classes are indexed from 0 to class_num-1. 
    71     IterablePartition(int class_num) { 
    72       for (int i=0; i<class_num; ++i) {
    73 	Tip t;
    74 	t.first=t.last=-1;
    75 	tips.push_back(t);
    76       }
    77     }
    78   protected:
    79     void befuz(Class it, int class_id) {
    80       if (tips[class_id].first==-1) {
    81 	if (tips[class_id].last==-1) {
    82 	  nodes[it.i].prev=nodes[it.i].next=-1;
    83 	  tips[class_id].first=tips[class_id].last=it.i;
    84 	}
    85       } else {
    86 	nodes[it.i].prev=tips[class_id].last;
    87 	nodes[it.i].next=-1;
    88 	nodes[tips[class_id].last].next=it.i;
    89 	tips[class_id].last=it.i;
    90       }
    91     }
    92     void kifuz(Class it, int class_id) {
    93       if (tips[class_id].first==it.i) {
    94 	if (tips[class_id].last==it.i) {
    95 	  tips[class_id].first=tips[class_id].last=-1;
    96 	} else {
    97 	  tips[class_id].first=nodes[it.i].next;
    98 	  nodes[nodes[it.i].next].prev=-1;
    99 	}
   100       } else {
   101 	if (tips[class_id].last==it.i) {
   102 	  tips[class_id].last=nodes[it.i].prev;
   103 	  nodes[nodes[it.i].prev].next=-1;
   104 	} else {
   105 	  nodes[nodes[it.i].next].prev=nodes[it.i].prev;
   106 	  nodes[nodes[it.i].prev].next=nodes[it.i].next;
   107 	}
   108       }
   109     }
   110   public:
   111     /// A new element with data \c t is pushed into the vector and into class 
   112     /// \c class_id.
   113     Class push_back(const T& t, int class_id) { 
   114       Node n;
   115       n.data=t;
   116       nodes.push_back(n);
   117       int i=nodes.size()-1;
   118       befuz(i, class_id);
   119       return i;
   120     }
   121     /// A member is moved to an other class.
   122     void set(Class it, int old_class_id, int new_class_id) {
   123       kifuz(it.i, old_class_id);
   124       befuz(it.i, new_class_id);
   125     }
   126     /// Returns the data pointed by \c it.
   127     T& operator[](Class it) { return nodes[it.i].data; }
   128     /// Returns the data pointed by \c it.
   129     const T& operator[](Class it) const { return nodes[it.i].data; }
   130     ///.
   131     class Class {
   132       friend class IterablePartition;
   133     protected:
   134       int i;
   135     public:
   136       /// Default constructor.
   137       Class() { }
   138       /// This constructor constructs an iterator which points
   139       /// to the member of th container indexed by the integer _i.
   140       Class(const int& _i) : i(_i) { }
   141       /// Invalid constructor.
   142       Class(const Invalid&) : i(-1) { }
   143       friend bool operator<(const Class& x, const Class& y);
   144       friend std::ostream& operator<<(std::ostream& os, 
   145 				      const Class& it);
   146       bool operator==(const Class& node) const {return i == node.i;}
   147       bool operator!=(const Class& node) const {return i != node.i;}
   148     };
   149     friend bool operator<(const Class& x, const Class& y) {
   150       return (x.i < y.i);
   151     }
   152     friend std::ostream& operator<<(std::ostream& os, 
   153 				    const Class& it) {
   154       os << it.i;
   155       return os;
   156     }
   157     /// First member of class \c class_id.
   158     Class& first(Class& it, int class_id) const {
   159       it.i=tips[class_id].first;
   160       return it;
   161     }
   162     /// Next member.
   163     Class& next(Class& it) const {
   164       it.i=nodes[it.i].next;
   165       return it;
   166     }
   167     /// True iff the iterator is valid.
   168     bool valid(const Class& it) const { return it.i!=-1; }
   169 
   170     class ClassIt : public Class {
   171       const IterablePartition* iterable_partition;
   172     public:
   173       ClassIt() { }
   174       ClassIt(Invalid i) : Class(i) { }
   175       ClassIt(const IterablePartition& _iterable_partition, 
   176 	      const int& i) : iterable_partition(&_iterable_partition) {
   177         _iterable_partition.first(*this, i);
   178       }
   179       ClassIt(const IterablePartition& _iterable_partition, 
   180 	      const Class& _class) : 
   181 	Class(_class), iterable_partition(&_iterable_partition) { }
   182       ClassIt& operator++() {
   183         iterable_partition->next(*this);
   184         return *this;
   185       }
   186     };
   187 
   188   };
   189 
   190 
   191   /*! \e
   192     \todo kellenene uj iterable structure bele, mert ez nem az igazi
   193     \todo A[x,y]-t cserel. Jobboldal, baloldal csere.
   194     \todo LEKERDEZESEK!!!
   195     \todo DOKSI!!!! Doxygen group!!!
   196     The aim of this class is to give a general surface to different 
   197     solvers, i.e. it makes possible to write algorithms using LP's, 
   198     in which the solver can be changed to an other one easily.
   199     \nosubgrouping
   200   */
   201   template <typename _Value>
   202   class LPSolverBase {
   203     
   204     /*! @name Uncategorized functions and types (public members)
   205     */
   206     //@{
   207   public:
   208 
   209     //UNCATEGORIZED
   210 
   211     /// \e
   212     typedef IterablePartition<int> Rows;
   213     /// \e
   214     typedef IterablePartition<int> Cols;
   215     /// \e
   216     typedef _Value Value;
   217     /// \e
   218     typedef Rows::Class Row;
   219     /// \e
   220     typedef Cols::Class Col;
   221   public:
   222     /// \e
   223     IterablePartition<int> row_iter_map;
   224     /// \e
   225     IterablePartition<int> col_iter_map;
   226     /// \e
   227     std::vector<Row> int_row_map;
   228     /// \e
   229     std::vector<Col> int_col_map;
   230     /// \e
   231     const int VALID_CLASS;
   232     /// \e
   233     const int INVALID_CLASS;
   234     /// \e 
   235     static const _Value INF;
   236   public:
   237     /// \e
   238     LPSolverBase() : row_iter_map(2), 
   239 		     col_iter_map(2), 
   240 		     VALID_CLASS(0), INVALID_CLASS(1) { }
   241     /// \e
   242     virtual ~LPSolverBase() { }
   243     //@}
   244 
   245     /*! @name Medium level interface (public members)
   246       These functions appear in the low level and also in the high level 
   247       interfaces thus these each of these functions have to be implemented 
   248       only once in the different interfaces.
   249       This means that these functions have to be reimplemented for all of the 
   250       different lp solvers. These are basic functions, and have the same 
   251       parameter lists in the low and high level interfaces. 
   252     */
   253     //@{
   254   public:
   255 
   256     //UNCATEGORIZED FUNCTIONS
   257 
   258     /// \e
   259     virtual void setMinimize() = 0;
   260     /// \e
   261     virtual void setMaximize() = 0;
   262 
   263     //SOLVER FUNCTIONS
   264 
   265     /// \e
   266     virtual void solveSimplex() = 0;
   267     /// \e
   268     virtual void solvePrimalSimplex() = 0;
   269     /// \e
   270     virtual void solveDualSimplex() = 0;
   271 
   272     //SOLUTION RETRIEVING
   273 
   274     /// \e
   275     virtual _Value getObjVal() = 0;
   276 
   277     //OTHER FUNCTIONS
   278 
   279     /// \e
   280     virtual int rowNum() const = 0;
   281     /// \e
   282     virtual int colNum() const = 0;
   283     /// \e
   284     virtual int warmUp() = 0;
   285     /// \e
   286     virtual void printWarmUpStatus(int i) = 0;
   287     /// \e
   288     virtual int getPrimalStatus() = 0;
   289     /// \e
   290     virtual void printPrimalStatus(int i) = 0;
   291     /// \e
   292     virtual int getDualStatus() = 0;
   293     /// \e
   294     virtual void printDualStatus(int i) = 0;
   295     /// Returns the status of the slack variable assigned to row \c row.
   296     virtual int getRowStat(const Row& row) = 0;
   297     /// \e
   298     virtual void printRowStatus(int i) = 0;
   299     /// Returns the status of the variable assigned to column \c col.
   300     virtual int getColStat(const Col& col) = 0;
   301     /// \e
   302     virtual void printColStatus(int i) = 0;
   303 
   304     //@}
   305 
   306     /*! @name Low level interface (protected members)
   307       Problem manipulating functions in the low level interface
   308     */
   309     //@{
   310   protected:
   311 
   312     //MATRIX MANIPULATING FUNCTIONS
   313 
   314     /// \e
   315     virtual int _addCol() = 0;
   316     /// \e
   317     virtual int _addRow() = 0;
   318     /// \e
   319     virtual void _eraseCol(int i) = 0;
   320     /// \e
   321     virtual void _eraseRow(int i) = 0;
   322     /// \e
   323     virtual void _setRowCoeffs(int i, 
   324 			       const std::vector<std::pair<int, _Value> >& coeffs) = 0;
   325     /// \e
   326     /// This routine modifies \c coeffs only by the \c push_back method.
   327     virtual void _getRowCoeffs(int i, 
   328 			       std::vector<std::pair<int, _Value> >& coeffs) = 0;
   329     /// \e
   330     virtual void _setColCoeffs(int i, 
   331 			       const std::vector<std::pair<int, _Value> >& coeffs) = 0;
   332     /// \e
   333     /// This routine modifies \c coeffs only by the \c push_back method.
   334     virtual void _getColCoeffs(int i, 
   335 			       std::vector<std::pair<int, _Value> >& coeffs) = 0;
   336     /// \e
   337     virtual void _setCoeff(int col, int row, _Value value) = 0;
   338     /// \e
   339     virtual _Value _getCoeff(int col, int row) = 0;
   340     //  public:
   341     //    /// \e
   342     //    enum Bound { FREE, LOWER, UPPER, DOUBLE, FIXED };
   343   protected:
   344     /// \e
   345     /// The lower bound of a variable (column) have to be given by an 
   346     /// extended number of type _Value, i.e. a finite number of type 
   347     /// _Value or -INF.
   348     virtual void _setColLowerBound(int i, _Value value) = 0;
   349     /// \e
   350     /// The lower bound of a variable (column) is an 
   351     /// extended number of type _Value, i.e. a finite number of type 
   352     /// _Value or -INF.
   353     virtual _Value _getColLowerBound(int i) = 0;
   354     /// \e
   355     /// The upper bound of a variable (column) have to be given by an 
   356     /// extended number of type _Value, i.e. a finite number of type 
   357     /// _Value or INF.
   358     virtual void _setColUpperBound(int i, _Value value) = 0;
   359     /// \e
   360     /// The upper bound of a variable (column) is an 
   361     /// extended number of type _Value, i.e. a finite number of type 
   362     /// _Value or INF.
   363     virtual _Value _getColUpperBound(int i) = 0;
   364     /// \e
   365     /// The lower bound of a linear expression (row) have to be given by an 
   366     /// extended number of type _Value, i.e. a finite number of type 
   367     /// _Value or -INF.
   368     virtual void _setRowLowerBound(int i, _Value value) = 0;
   369     /// \e
   370     /// The lower bound of a linear expression (row) is an 
   371     /// extended number of type _Value, i.e. a finite number of type 
   372     /// _Value or -INF.
   373     virtual _Value _getRowLowerBound(int i) = 0;
   374     /// \e
   375     /// The upper bound of a linear expression (row) have to be given by an 
   376     /// extended number of type _Value, i.e. a finite number of type 
   377     /// _Value or INF.
   378     virtual void _setRowUpperBound(int i, _Value value) = 0;
   379     /// \e
   380     /// The upper bound of a linear expression (row) is an 
   381     /// extended number of type _Value, i.e. a finite number of type 
   382     /// _Value or INF.
   383     virtual _Value _getRowUpperBound(int i) = 0;
   384     /// \e
   385     virtual void _setObjCoeff(int i, _Value obj_coef) = 0;
   386     /// \e
   387     virtual _Value _getObjCoeff(int i) = 0;
   388     
   389     //SOLUTION RETRIEVING
   390 
   391     /// \e
   392     virtual _Value _getPrimal(int i) = 0;
   393     //@}
   394     
   395     /*! @name High level interface (public members)
   396       Problem manipulating functions in the high level interface
   397     */
   398     //@{
   399   public:
   400 
   401     //MATRIX MANIPULATING FUNCTIONS
   402 
   403     /// \e
   404     Col addCol() {
   405       int i=_addCol();  
   406       Col col;
   407       col_iter_map.first(col, INVALID_CLASS);
   408       if (col_iter_map.valid(col)) { //van hasznalhato hely
   409 	col_iter_map.set(col, INVALID_CLASS, VALID_CLASS);
   410 	col_iter_map[col]=i;
   411       } else { //a cucc vegere kell inzertalni mert nincs szabad hely
   412 	col=col_iter_map.push_back(i, VALID_CLASS);
   413       }
   414       int_col_map.push_back(col);
   415       return col;
   416     }
   417     /// \e
   418     Row addRow() {
   419       int i=_addRow();
   420       Row row;
   421       row_iter_map.first(row, INVALID_CLASS);
   422       if (row_iter_map.valid(row)) { //van hasznalhato hely
   423 	row_iter_map.set(row, INVALID_CLASS, VALID_CLASS);
   424 	row_iter_map[row]=i;
   425       } else { //a cucc vegere kell inzertalni mert nincs szabad hely
   426 	row=row_iter_map.push_back(i, VALID_CLASS);
   427       }
   428       int_row_map.push_back(row);
   429       return row;
   430     }
   431     /// \e
   432     void eraseCol(const Col& col) {
   433       col_iter_map.set(col, VALID_CLASS, INVALID_CLASS);
   434       int cols[2];
   435       cols[1]=col_iter_map[col];
   436       _eraseCol(cols[1]);
   437       col_iter_map[col]=0; //glpk specifikus, de kell ez??
   438       Col it;
   439       for (col_iter_map.first(it, VALID_CLASS); 
   440 	   col_iter_map.valid(it); col_iter_map.next(it)) {
   441 	if (col_iter_map[it]>cols[1]) --col_iter_map[it];
   442       }
   443       int_col_map.erase(int_col_map.begin()+cols[1]);
   444     }
   445     /// \e
   446     void eraseRow(const Row& row) {
   447       row_iter_map.set(row, VALID_CLASS, INVALID_CLASS);
   448       int rows[2];
   449       rows[1]=row_iter_map[row];
   450       _eraseRow(rows[1]);
   451       row_iter_map[row]=0; //glpk specifikus, de kell ez??
   452       Row it;
   453       for (row_iter_map.first(it, VALID_CLASS); 
   454 	   row_iter_map.valid(it); row_iter_map.next(it)) {
   455 	if (row_iter_map[it]>rows[1]) --row_iter_map[it];
   456       }
   457       int_row_map.erase(int_row_map.begin()+rows[1]);
   458     }
   459     /// \e
   460     void setCoeff(Col col, Row row, _Value value) {
   461       _setCoeff(col_iter_map[col], row_iter_map[row], value);
   462     }
   463     /// \e
   464     _Value getCoeff(Col col, Row row) {
   465       return _getCoeff(col_iter_map[col], row_iter_map[row], value);
   466     }
   467     /// \e
   468     void setColLowerBound(Col col, _Value lo) {
   469       _setColLowerBound(col_iter_map[col], lo);
   470     }
   471     /// \e
   472     _Value getColLowerBound(Col col) {
   473       return _getColLowerBound(col_iter_map[col]);
   474     }
   475     /// \e
   476     void setColUpperBound(Col col, _Value up) {
   477       _setColUpperBound(col_iter_map[col], up);
   478     }
   479     /// \e
   480     _Value getColUpperBound(Col col) {      
   481       return _getColUpperBound(col_iter_map[col]);
   482     }
   483     /// \e
   484     void setRowLowerBound(Row row, _Value lo) {
   485       _setRowLowerBound(row_iter_map[row], lo);
   486     }
   487     /// \e
   488     _Value getRowLowerBound(Row row) {
   489       return _getRowLowerBound(row_iter_map[row]);
   490     }
   491     /// \e
   492     void setRowUpperBound(Row row, _Value up) {
   493       _setRowUpperBound(row_iter_map[row], up);
   494     }
   495     /// \e
   496     _Value getRowUpperBound(Row row) {      
   497       return _getRowUpperBound(row_iter_map[row]);
   498     }
   499     /// \e
   500     void setObjCoeff(const Col& col, _Value obj_coef) {
   501       _setObjCoeff(col_iter_map[col], obj_coef);
   502     }
   503     /// \e
   504     _Value getObjCoeff(const Col& col) {
   505       return _getObjCoeff(col_iter_map[col]);
   506     }
   507 
   508     //SOLUTION RETRIEVING FUNCTIONS
   509 
   510     /// \e
   511     _Value getPrimal(const Col& col) {
   512       return _getPrimal(col_iter_map[col]);
   513     }    
   514 
   515     //@}
   516 
   517     /*! @name User friend interface
   518       Problem manipulating functions in the user friend interface
   519     */
   520     //@{
   521 
   522     //EXPRESSION TYPES
   523 
   524     /// \e
   525     typedef Expr<Col, _Value> Expression;
   526     /// \e
   527     typedef Expr<Row, _Value> DualExpression;
   528     /// \e
   529     typedef Constr<Col, _Value> Constraint;
   530 
   531     //MATRIX MANIPULATING FUNCTIONS
   532 
   533     /// \e
   534     void setRowCoeffs(Row row, const Expression& expr) {
   535       std::vector<std::pair<int, _Value> > row_coeffs;
   536       for(typename Expression::Data::const_iterator i=expr.data.begin(); 
   537 	  i!=expr.data.end(); ++i) {
   538 	row_coeffs.push_back(std::make_pair
   539 			     (col_iter_map[(*i).first], (*i).second));
   540       }
   541       _setRowCoeffs(row_iter_map[row], row_coeffs);
   542     }
   543     /// \e 
   544     void setRow(Row row, const Constraint& constr) {
   545       setRowCoeffs(row, constr.expr);
   546       setRowLowerBound(row, constr.lo);
   547       setRowUpperBound(row, constr.up);
   548     }
   549     /// \e 
   550     Row addRow(const Constraint& constr) {
   551       Row row=addRow();
   552       setRowCoeffs(row, constr.expr);
   553       setRowLowerBound(row, constr.lo);
   554       setRowUpperBound(row, constr.up);
   555       return row;
   556     }
   557     /// \e
   558     /// This routine modifies \c expr by only adding to it.
   559     void getRowCoeffs(Row row, Expression& expr) {
   560       std::vector<std::pair<int, _Value> > row_coeffs;
   561       _getRowCoeffs(row_iter_map[row], row_coeffs);
   562       for(typename std::vector<std::pair<int, _Value> >::const_iterator 
   563  	    i=row_coeffs.begin(); i!=row_coeffs.end(); ++i) {
   564  	expr+= (*i).second*int_col_map[(*i).first];
   565       }
   566     }
   567     /// \e
   568     void setColCoeffs(Col col, const DualExpression& expr) {
   569       std::vector<std::pair<int, _Value> > col_coeffs;
   570       for(typename DualExpression::Data::const_iterator i=expr.data.begin(); 
   571 	  i!=expr.data.end(); ++i) {
   572 	col_coeffs.push_back(std::make_pair
   573 			     (row_iter_map[(*i).first], (*i).second));
   574       }
   575       _setColCoeffs(col_iter_map[col], col_coeffs);
   576     }
   577     /// \e
   578     /// This routine modifies \c expr by only adding to it.
   579     void getColCoeffs(Col col, DualExpression& expr) {
   580       std::vector<std::pair<int, _Value> > col_coeffs;
   581       _getColCoeffs(col_iter_map[col], col_coeffs);
   582       for(typename std::vector<std::pair<int, _Value> >::const_iterator 
   583  	    i=col_coeffs.begin(); i!=col_coeffs.end(); ++i) {
   584  	expr+= (*i).second*int_row_map[(*i).first];
   585       }
   586     }
   587     /// \e
   588     void setObjCoeffs(const Expression& expr) {
   589       // writing zero everywhere
   590       for(Cols::ClassIt it(col_iter_map, VALID_CLASS); it!=INVALID; ++it)
   591 	setObjCoeff(it, 0.0);
   592       // writing the data needed
   593       for(typename Expression::Data::const_iterator i=expr.data.begin(); 
   594 	  i!=expr.data.end(); ++i) {
   595 	setObjCoeff((*i).first, (*i).second);
   596       }
   597     }
   598     /// \e
   599     /// This routine modifies \c expr by only adding to it.
   600     void getObjCoeffs(Expression& expr) {
   601       for(Cols::ClassIt it(col_iter_map, VALID_CLASS); it!=INVALID; ++it)
   602 	expr+=getObjCoeff(it)*it;
   603     }
   604     //@}
   605 
   606 
   607     /*! @name MIP functions and types (public members)
   608     */
   609     //@{
   610   public:
   611     /// \e
   612     virtual void solveBandB() = 0;
   613     /// \e
   614     virtual void setLP() = 0;
   615     /// \e
   616     virtual void setMIP() = 0;
   617   protected:
   618    /// \e
   619     virtual void _setColCont(int i) = 0;
   620     /// \e
   621     virtual void _setColInt(int i) = 0;
   622     /// \e
   623     virtual _Value _getMIPPrimal(int i) = 0;
   624   public:
   625     /// \e
   626     void setColCont(Col col) {
   627       _setColCont(col_iter_map[col]);
   628     }
   629     /// \e
   630     void setColInt(Col col) {
   631       _setColInt(col_iter_map[col]);
   632     }
   633     /// \e
   634     _Value getMIPPrimal(Col col) {
   635       return _getMIPPrimal(col_iter_map[col]);
   636     }
   637     //@}
   638   };
   639   
   640   template <typename _Value>
   641   const _Value LPSolverBase<_Value>::INF=std::numeric_limits<_Value>::infinity();
   642 
   643 
   644   /// \brief Wrapper for GLPK solver
   645   /// 
   646   /// This class implements a lemon wrapper for GLPK.
   647   class LPGLPK : public LPSolverBase<double> {
   648   public:
   649     typedef LPSolverBase<double> Parent;
   650 
   651   public:
   652     /// \e
   653     LPX* lp;
   654 
   655   public:
   656     /// \e
   657     LPGLPK() : Parent(), 
   658 			lp(lpx_create_prob()) {
   659       int_row_map.push_back(Row());
   660       int_col_map.push_back(Col());
   661       lpx_set_int_parm(lp, LPX_K_DUAL, 1);
   662     }
   663     /// \e
   664     ~LPGLPK() {
   665       lpx_delete_prob(lp);
   666     }
   667 
   668     //MATRIX INDEPEDENT MANIPULATING FUNCTIONS
   669 
   670     /// \e
   671     void setMinimize() { 
   672       lpx_set_obj_dir(lp, LPX_MIN);
   673     }
   674     /// \e
   675     void setMaximize() { 
   676       lpx_set_obj_dir(lp, LPX_MAX);
   677     }
   678 
   679     //LOW LEVEL INTERFACE, MATRIX MANIPULATING FUNCTIONS
   680 
   681   protected:
   682     /// \e
   683     int _addCol() { 
   684       int i=lpx_add_cols(lp, 1);
   685       _setColLowerBound(i, -INF);
   686       _setColUpperBound(i, INF);
   687       return i;
   688     }
   689     /// \e
   690     int _addRow() { 
   691       int i=lpx_add_rows(lp, 1);
   692       return i;
   693     }
   694     /// \e
   695     virtual void _setRowCoeffs(int i, 
   696 			       const std::vector<std::pair<int, double> >& coeffs) {
   697       int mem_length=1+colNum();
   698       int* indices = new int[mem_length];
   699       double* doubles = new double[mem_length];
   700       int length=0;
   701       for (std::vector<std::pair<int, double> >::
   702 	     const_iterator it=coeffs.begin(); it!=coeffs.end(); ++it) {
   703 	++length;
   704 	indices[length]=it->first;
   705 	doubles[length]=it->second;
   706       }
   707       lpx_set_mat_row(lp, i, length, indices, doubles);
   708       delete [] indices;
   709       delete [] doubles;
   710     }
   711     /// \e
   712     virtual void _getRowCoeffs(int i, 
   713 			       std::vector<std::pair<int, double> >& coeffs) {
   714       int mem_length=1+colNum();
   715       int* indices = new int[mem_length];
   716       double* doubles = new double[mem_length];
   717       int length=lpx_get_mat_row(lp, i, indices, doubles);
   718       for (int i=1; i<=length; ++i) {
   719 	coeffs.push_back(std::make_pair(indices[i], doubles[i]));
   720       }
   721       delete [] indices;
   722       delete [] doubles;
   723     }
   724     /// \e
   725     virtual void _setColCoeffs(int i, 
   726 			       const std::vector<std::pair<int, double> >& coeffs) {
   727       int mem_length=1+rowNum();
   728       int* indices = new int[mem_length];
   729       double* doubles = new double[mem_length];
   730       int length=0;
   731       for (std::vector<std::pair<int, double> >::
   732 	     const_iterator it=coeffs.begin(); it!=coeffs.end(); ++it) {
   733 	++length;
   734 	indices[length]=it->first;
   735 	doubles[length]=it->second;
   736       }
   737       lpx_set_mat_col(lp, i, length, indices, doubles);
   738       delete [] indices;
   739       delete [] doubles;
   740     }
   741     /// \e
   742     virtual void _getColCoeffs(int i, 
   743 			       std::vector<std::pair<int, double> >& coeffs) {
   744       int mem_length=1+rowNum();
   745       int* indices = new int[mem_length];
   746       double* doubles = new double[mem_length];
   747       int length=lpx_get_mat_col(lp, i, indices, doubles);
   748       for (int i=1; i<=length; ++i) {
   749 	coeffs.push_back(std::make_pair(indices[i], doubles[i]));
   750       }
   751       delete [] indices;
   752       delete [] doubles;
   753     }
   754     /// \e
   755     virtual void _eraseCol(int i) {
   756       int cols[2];
   757       cols[1]=i;
   758       lpx_del_cols(lp, 1, cols);
   759     }
   760     virtual void _eraseRow(int i) {
   761       int rows[2];
   762       rows[1]=i;
   763       lpx_del_rows(lp, 1, rows);
   764     }
   765     void _setCoeff(int col, int row, double value) {
   766       /// FIXME not yet implemented
   767     }
   768     double _getCoeff(int col, int row) {
   769       /// FIXME not yet implemented
   770       return 0.0;
   771     }
   772     virtual void _setColLowerBound(int i, double lo) {
   773       if (lo==INF) {
   774 	//FIXME error
   775       }
   776       int b=lpx_get_col_type(lp, i);
   777       double up=lpx_get_col_ub(lp, i);	
   778       if (lo==-INF) {
   779 	switch (b) {
   780 	case LPX_FR:
   781 	case LPX_LO:
   782 	  lpx_set_col_bnds(lp, i, LPX_FR, lo, up);
   783 	  break;
   784 	case LPX_UP:
   785 	  break;
   786 	case LPX_DB:
   787 	case LPX_FX:
   788 	  lpx_set_col_bnds(lp, i, LPX_UP, lo, up);
   789 	  break;
   790 	default: ;
   791 	  //FIXME error
   792 	}
   793       } else {
   794 	switch (b) {
   795 	case LPX_FR:
   796 	case LPX_LO:
   797 	  lpx_set_col_bnds(lp, i, LPX_LO, lo, up);
   798 	  break;
   799 	case LPX_UP:	  
   800 	case LPX_DB:
   801 	case LPX_FX:
   802 	  if (lo==up) 
   803 	    lpx_set_col_bnds(lp, i, LPX_FX, lo, up);
   804 	  else 
   805 	    lpx_set_col_bnds(lp, i, LPX_DB, lo, up);
   806 	  break;
   807 	default: ;
   808 	  //FIXME error
   809 	}
   810       }
   811     }
   812     virtual double _getColLowerBound(int i) {
   813       int b=lpx_get_col_type(lp, i);
   814       switch (b) {
   815       case LPX_FR:
   816 	return -INF;
   817       case LPX_LO:
   818 	return lpx_get_col_lb(lp, i);
   819       case LPX_UP:
   820 	return -INF;
   821       case LPX_DB:
   822       case LPX_FX:
   823 	return lpx_get_col_lb(lp, i);
   824       default: ;
   825 	//FIXME error
   826 	return 0.0;
   827       }
   828     }
   829     virtual void _setColUpperBound(int i, double up) {
   830       if (up==-INF) {
   831 	//FIXME error
   832       }
   833       int b=lpx_get_col_type(lp, i);
   834       double lo=lpx_get_col_lb(lp, i);
   835       if (up==INF) {
   836 	switch (b) {
   837 	case LPX_FR:
   838 	case LPX_LO:
   839 	  break;
   840 	case LPX_UP:
   841 	  lpx_set_col_bnds(lp, i, LPX_FR, lo, up);
   842 	  break;
   843 	case LPX_DB:
   844 	case LPX_FX:
   845 	  lpx_set_col_bnds(lp, i, LPX_LO, lo, up);
   846 	  break;
   847 	default: ;
   848 	  //FIXME error
   849 	}
   850       } else {
   851 	switch (b) {
   852 	case LPX_FR:
   853 	  lpx_set_col_bnds(lp, i, LPX_UP, lo, up);
   854 	  break;
   855 	case LPX_LO:
   856 	  if (lo==up) 
   857 	    lpx_set_col_bnds(lp, i, LPX_FX, lo, up);
   858 	  else
   859 	    lpx_set_col_bnds(lp, i, LPX_DB, lo, up);
   860 	  break;
   861 	case LPX_UP:
   862 	  lpx_set_col_bnds(lp, i, LPX_UP, lo, up);
   863 	  break;
   864 	case LPX_DB:
   865 	case LPX_FX:
   866 	  if (lo==up) 
   867 	    lpx_set_col_bnds(lp, i, LPX_FX, lo, up);
   868 	  else 
   869 	    lpx_set_col_bnds(lp, i, LPX_DB, lo, up);
   870 	  break;
   871 	default: ;
   872 	  //FIXME error
   873 	}
   874       }
   875     }
   876     virtual double _getColUpperBound(int i) {
   877       int b=lpx_get_col_type(lp, i);
   878       switch (b) {
   879       case LPX_FR:
   880       case LPX_LO:
   881 	return INF;
   882       case LPX_UP:
   883       case LPX_DB:
   884       case LPX_FX:
   885 	return lpx_get_col_ub(lp, i);
   886       default: ;
   887 	//FIXME error
   888 	return 0.0;
   889       }
   890     }
   891     virtual void _setRowLowerBound(int i, double lo) {
   892       if (lo==INF) {
   893 	//FIXME error
   894       }
   895       int b=lpx_get_row_type(lp, i);
   896       double up=lpx_get_row_ub(lp, i);	
   897       if (lo==-INF) {
   898 	switch (b) {
   899 	case LPX_FR:
   900 	case LPX_LO:
   901 	  lpx_set_row_bnds(lp, i, LPX_FR, lo, up);
   902 	  break;
   903 	case LPX_UP:
   904 	  break;
   905 	case LPX_DB:
   906 	case LPX_FX:
   907 	  lpx_set_row_bnds(lp, i, LPX_UP, lo, up);
   908 	  break;
   909 	default: ;
   910 	  //FIXME error
   911 	}
   912       } else {
   913 	switch (b) {
   914 	case LPX_FR:
   915 	case LPX_LO:
   916 	  lpx_set_row_bnds(lp, i, LPX_LO, lo, up);
   917 	  break;
   918 	case LPX_UP:	  
   919 	case LPX_DB:
   920 	case LPX_FX:
   921 	  if (lo==up) 
   922 	    lpx_set_row_bnds(lp, i, LPX_FX, lo, up);
   923 	  else 
   924 	    lpx_set_row_bnds(lp, i, LPX_DB, lo, up);
   925 	  break;
   926 	default: ;
   927 	  //FIXME error
   928 	}
   929       }
   930     }
   931     virtual double _getRowLowerBound(int i) {
   932       int b=lpx_get_row_type(lp, i);
   933       switch (b) {
   934       case LPX_FR:
   935 	return -INF;
   936       case LPX_LO:
   937 	return lpx_get_row_lb(lp, i);
   938       case LPX_UP:
   939 	return -INF;
   940       case LPX_DB:
   941       case LPX_FX:
   942 	return lpx_get_row_lb(lp, i);
   943       default: ;
   944 	//FIXME error
   945 	return 0.0;
   946       }
   947     }
   948     virtual void _setRowUpperBound(int i, double up) {
   949       if (up==-INF) {
   950 	//FIXME error
   951       }
   952       int b=lpx_get_row_type(lp, i);
   953       double lo=lpx_get_row_lb(lp, i);
   954       if (up==INF) {
   955 	switch (b) {
   956 	case LPX_FR:
   957 	case LPX_LO:
   958 	  break;
   959 	case LPX_UP:
   960 	  lpx_set_row_bnds(lp, i, LPX_FR, lo, up);
   961 	  break;
   962 	case LPX_DB:
   963 	case LPX_FX:
   964 	  lpx_set_row_bnds(lp, i, LPX_LO, lo, up);
   965 	  break;
   966 	default: ;
   967 	  //FIXME error
   968 	}
   969       } else {
   970 	switch (b) {
   971 	case LPX_FR:
   972 	  lpx_set_row_bnds(lp, i, LPX_UP, lo, up);
   973 	  break;
   974 	case LPX_LO:
   975 	  if (lo==up) 
   976 	    lpx_set_row_bnds(lp, i, LPX_FX, lo, up);
   977 	  else
   978 	    lpx_set_row_bnds(lp, i, LPX_DB, lo, up);
   979 	  break;
   980 	case LPX_UP:
   981 	  lpx_set_row_bnds(lp, i, LPX_UP, lo, up);
   982 	  break;
   983 	case LPX_DB:
   984 	case LPX_FX:
   985 	  if (lo==up) 
   986 	    lpx_set_row_bnds(lp, i, LPX_FX, lo, up);
   987 	  else 
   988 	    lpx_set_row_bnds(lp, i, LPX_DB, lo, up);
   989 	  break;
   990 	default: ;
   991 	  //FIXME error
   992 	}
   993       }
   994     }
   995     virtual double _getRowUpperBound(int i) {
   996       int b=lpx_get_row_type(lp, i);
   997       switch (b) {
   998       case LPX_FR:
   999       case LPX_LO:
  1000 	return INF;
  1001       case LPX_UP:
  1002       case LPX_DB:
  1003       case LPX_FX:
  1004 	return lpx_get_row_ub(lp, i);
  1005       default: ;
  1006 	//FIXME error
  1007 	return 0.0;
  1008       }
  1009     }
  1010     /// \e
  1011     virtual double _getObjCoeff(int i) { 
  1012       return lpx_get_obj_coef(lp, i);
  1013     }
  1014     /// \e
  1015     virtual void _setObjCoeff(int i, double obj_coef) { 
  1016       lpx_set_obj_coef(lp, i, obj_coef);
  1017     }
  1018   public:
  1019     /// \e
  1020     void solveSimplex() { lpx_simplex(lp); }
  1021     /// \e
  1022     void solvePrimalSimplex() { lpx_simplex(lp); }
  1023     /// \e
  1024     void solveDualSimplex() { lpx_simplex(lp); }
  1025   protected:
  1026     virtual double _getPrimal(int i) {
  1027       return lpx_get_col_prim(lp, i);
  1028     }
  1029   public:
  1030     /// \e
  1031     double getObjVal() { return lpx_get_obj_val(lp); }
  1032     /// \e
  1033     int rowNum() const { return lpx_get_num_rows(lp); }
  1034     /// \e
  1035     int colNum() const { return lpx_get_num_cols(lp); }
  1036     /// \e
  1037     int warmUp() { return lpx_warm_up(lp); }
  1038     /// \e
  1039     void printWarmUpStatus(int i) {
  1040       switch (i) {
  1041       case LPX_E_OK: cout << "LPX_E_OK" << endl; break;
  1042       case LPX_E_EMPTY: cout << "LPX_E_EMPTY" << endl; break;	
  1043       case LPX_E_BADB: cout << "LPX_E_BADB" << endl; break;
  1044       case LPX_E_SING: cout << "LPX_E_SING" << endl; break;
  1045       }
  1046     }
  1047     /// \e
  1048     int getPrimalStatus() { return lpx_get_prim_stat(lp); }
  1049     /// \e
  1050     void printPrimalStatus(int i) {
  1051       switch (i) {
  1052       case LPX_P_UNDEF: cout << "LPX_P_UNDEF" << endl; break;
  1053       case LPX_P_FEAS: cout << "LPX_P_FEAS" << endl; break;	
  1054       case LPX_P_INFEAS: cout << "LPX_P_INFEAS" << endl; break;
  1055       case LPX_P_NOFEAS: cout << "LPX_P_NOFEAS" << endl; break;
  1056       }
  1057     }
  1058     /// \e
  1059     int getDualStatus() { return lpx_get_dual_stat(lp); }
  1060     /// \e
  1061     void printDualStatus(int i) {
  1062       switch (i) {
  1063       case LPX_D_UNDEF: cout << "LPX_D_UNDEF" << endl; break;
  1064       case LPX_D_FEAS: cout << "LPX_D_FEAS" << endl; break;	
  1065       case LPX_D_INFEAS: cout << "LPX_D_INFEAS" << endl; break;
  1066       case LPX_D_NOFEAS: cout << "LPX_D_NOFEAS" << endl; break;
  1067       }
  1068     }
  1069     /// Returns the status of the slack variable assigned to row \c row.
  1070     int getRowStat(const Row& row) { 
  1071       return lpx_get_row_stat(lp, row_iter_map[row]); 
  1072     }
  1073     /// \e
  1074     void printRowStatus(int i) {
  1075       switch (i) {
  1076       case LPX_BS: cout << "LPX_BS" << endl; break;
  1077       case LPX_NL: cout << "LPX_NL" << endl; break;	
  1078       case LPX_NU: cout << "LPX_NU" << endl; break;
  1079       case LPX_NF: cout << "LPX_NF" << endl; break;
  1080       case LPX_NS: cout << "LPX_NS" << endl; break;
  1081       }
  1082     }
  1083     /// Returns the status of the variable assigned to column \c col.
  1084     int getColStat(const Col& col) { 
  1085       return lpx_get_col_stat(lp, col_iter_map[col]); 
  1086     }
  1087     /// \e
  1088     void printColStatus(int i) {
  1089       switch (i) {
  1090       case LPX_BS: cout << "LPX_BS" << endl; break;
  1091       case LPX_NL: cout << "LPX_NL" << endl; break;	
  1092       case LPX_NU: cout << "LPX_NU" << endl; break;
  1093       case LPX_NF: cout << "LPX_NF" << endl; break;
  1094       case LPX_NS: cout << "LPX_NS" << endl; break;
  1095       }
  1096     }
  1097 
  1098     // MIP
  1099     /// \e
  1100     void solveBandB() { lpx_integer(lp); }
  1101     /// \e
  1102     void setLP() { lpx_set_class(lp, LPX_LP); }
  1103     /// \e
  1104     void setMIP() { lpx_set_class(lp, LPX_MIP); }
  1105   protected:
  1106     /// \e
  1107     void _setColCont(int i) { lpx_set_col_kind(lp, i, LPX_CV); }
  1108     /// \e
  1109     void _setColInt(int i) { lpx_set_col_kind(lp, i, LPX_IV); }
  1110     /// \e
  1111     double _getMIPPrimal(int i) { return lpx_mip_col_val(lp, i); }
  1112   };
  1113   
  1114   /// @}
  1115 
  1116 } //namespace lemon
  1117 
  1118 #endif //LEMON_LP_SOLVER_BASE_H