lemon-0.x: src/work/marci/lp/lp_solver_wrapper

     1 // -*- c++ -*-

     2 #ifndef LEMON_LP_SOLVER_WRAPPER_H

     3 #define LEMON_LP_SOLVER_WRAPPER_H

     5 ///\ingroup misc

     6 ///\file

     7 ///\brief Dijkstra algorithm.

     9 // #include <stdio.h>

    10 #include <stdlib.h>

    11 // #include <stdio>

    12 //#include <stdlib>

    13 extern "C" {

    14 #include "glpk.h"

    15 }

    17 #include <iostream>

    18 #include <vector>

    19 #include <string>

    20 #include <list>

    21 #include <memory>

    22 #include <utility>

    24 //#include <sage_graph.h>

    25 //#include <lemon/list_graph.h>

    26 //#include <lemon/graph_wrapper.h>

    27 #include <lemon/invalid.h>

    28 //#include <bfs_dfs.h>

    29 //#include <stp.h>

    30 //#include <lemon/max_flow.h>

    31 //#include <augmenting_flow.h>

    32 //#include <iter_map.h>

    34 using std::cout;

    35 using std::cin;

    36 using std::endl;

    38 namespace lemon {

    41   /// \addtogroup misc

    42   /// @{

    44   /// \brief A partitioned vector with iterable classes.

    45   ///

    46   /// This class implements a container in which the data is stored in an

    47   /// stl vector, the range is partitioned into sets and each set is

    48   /// doubly linked in a list.

    49   /// That is, each class is iterable by lemon iterators, and any member of

    50   /// the vector can bo moved to an other class.

    51   template <typename T>

    52   class IterablePartition {

    53   protected:

    54     struct Node {

    55       T data;

    56       int prev; //invalid az -1

    57       int next;

    58     };

    59     std::vector<Node> nodes;

    60     struct Tip {

    61       int first;

    62       int last;

    63     };

    64     std::vector<Tip> tips;

    65   public:

    66     /// The classes are indexed by integers from \c 0 to \c classNum()-1.

    67     int classNum() const { return tips.size(); }

    68     /// This lemon style iterator iterates through a class.

    69     class ClassIt;

    70     /// Constructor. The number of classes is to be given which is fixed

    71     /// over the life of the container.

    72     /// The partition classes are indexed from 0 to class_num-1.

    73     IterablePartition(int class_num) {

    74       for (int i=0; i<class_num; ++i) {

    75 	Tip t;

    76 	t.first=t.last=-1;

    77 	tips.push_back(t);

    78       }

    79     }

    80   protected:

    81     void befuz(ClassIt it, int class_id) {

    82       if (tips[class_id].first==-1) {

    83 	if (tips[class_id].last==-1) {

    84 	  nodes[it.i].prev=nodes[it.i].next=-1;

    85 	  tips[class_id].first=tips[class_id].last=it.i;

    86 	}

    87       } else {

    88 	nodes[it.i].prev=tips[class_id].last;

    89 	nodes[it.i].next=-1;

    90 	nodes[tips[class_id].last].next=it.i;

    91 	tips[class_id].last=it.i;

    92       }

    93     }

    94     void kifuz(ClassIt it, int class_id) {

    95       if (tips[class_id].first==it.i) {

    96 	if (tips[class_id].last==it.i) {

    97 	  tips[class_id].first=tips[class_id].last=-1;

    98 	} else {

    99 	  tips[class_id].first=nodes[it.i].next;

   100 	  nodes[nodes[it.i].next].prev=-1;

   101 	}

   102       } else {

   103 	if (tips[class_id].last==it.i) {

   104 	  tips[class_id].last=nodes[it.i].prev;

   105 	  nodes[nodes[it.i].prev].next=-1;

   106 	} else {

   107 	  nodes[nodes[it.i].next].prev=nodes[it.i].prev;

   108 	  nodes[nodes[it.i].prev].next=nodes[it.i].next;

   109 	}

   110       }

   111     }

   112   public:

   113     /// A new element with data \c t is pushed into the vector and into class

   114     /// \c class_id.

   115     ClassIt push_back(const T& t, int class_id) {

   116       Node n;

   117       n.data=t;

   118       nodes.push_back(n);

   119       int i=nodes.size()-1;

   120       befuz(i, class_id);

   121       return i;

   122     }

   123     /// A member is moved to an other class.

   124     void set(ClassIt it, int old_class_id, int new_class_id) {

   125       kifuz(it.i, old_class_id);

   126       befuz(it.i, new_class_id);

   127     }

   128     /// Returns the data pointed by \c it.

   129     T& operator[](ClassIt it) { return nodes[it.i].data; }

   130     /// Returns the data pointed by \c it.

   131     const T& operator[](ClassIt it) const { return nodes[it.i].data; }

   132     ///.

   133     class ClassIt {

   134       friend class IterablePartition;

   135     protected:

   136       int i;

   137     public:

   138       /// Default constructor.

   139       ClassIt() { }

   140       /// This constructor constructs an iterator which points

   141       /// to the member of th container indexed by the integer _i.

   142       ClassIt(const int& _i) : i(_i) { }

   143       /// Invalid constructor.

   144       ClassIt(const Invalid&) : i(-1) { }

   145     };

   146     /// First member of class \c class_id.

   147     ClassIt& first(ClassIt& it, int class_id) const {

   148       it.i=tips[class_id].first;

   149       return it;

   150     }

   151     /// Next member.

   152     ClassIt& next(ClassIt& it) const {

   153       it.i=nodes[it.i].next;

   154       return it;

   155     }

   156     /// True iff the iterator is valid.

   157     bool valid(const ClassIt& it) const { return it.i!=-1; }

   158   };

   160   /*! \e

   161    */

   162   template <typename _Value>

   163   class LPSolverBase {

   164   public:

   165     /// \e

   166     typedef _Value Value;

   167     /// \e

   168     typedef IterablePartition<int>::ClassIt RowIt;

   169     /// \e

   170     typedef IterablePartition<int>::ClassIt ColIt;

   171   protected:

   172     /// \e

   173     IterablePartition<int> row_iter_map;

   174     /// \e

   175     IterablePartition<int> col_iter_map;

   176     /// \e

   177     const int VALID_ID;

   178     /// \e

   179     const int INVALID_ID;

   180   public:

   181     /// \e

   182     LPSolverBase() : row_iter_map(2),

   183 		     col_iter_map(2),

   184 		     VALID_ID(0), INVALID_ID(1) { }

   185     /// \e

   186     virtual ~LPSolverBase() { }

   187     /// \e

   188     virtual void setMinimize() = 0;

   189     /// \e

   190     virtual void setMaximize() = 0;

   191     /// \e

   192     virtual RowIt addRow() = 0;

   193     /// \e

   194     virtual ColIt addCol() = 0;

   195     /// temporally, glpk style indexing

   196     virtual void setRowCoeffs(RowIt row_it, int num,

   197 			      int* indices, _Value* doubles) = 0;

   198     //pair<RowIt, _Value>-bol kell megadni egy std range-et

   199     /// \e

   200     template <typename Begin, typename End>

   201     void setRowCoeffs(RowIt row_it, Begin begin, End end) {

   202       int mem_length=1+colNum();

   203       int* indices = new int[mem_length];

   204       _Value* doubles = new _Value[mem_length];

   205       int length=0;

   206       for ( ; begin!=end; ++begin) {

   207 	++length;

   208 	indices[length]=col_iter_map[begin->first];

   209 	doubles[length]=begin->second;

   210       }

   211       setRowCoeffs(row_it, length, indices, doubles);

   212       delete [] indices;

   213       delete [] doubles;

   214     }

   215     /// temporally, glpk style indexing

   216     virtual void setColCoeffs(ColIt col_it, int num,

   217 			      int* indices, _Value* doubles) = 0;

   218     //pair<ColIt, _Value>-bol kell megadni egy std range-et

   219     /// \e

   220     template <typename Begin, typename End>

   221     void setColCoeffs(ColIt col_it, Begin begin, End end) {

   222       int mem_length=1+rowNum();

   223       int* indices = new int[mem_length];

   224       _Value* doubles = new _Value[mem_length];

   225       int length=0;

   226       for ( ; begin!=end; ++begin) {

   227 	++length;

   228 	indices[length]=row_iter_map[begin->first];

   229 	doubles[length]=begin->second;

   230       }

   231       setColCoeffs(col_it, length, indices, doubles);

   232       delete [] indices;

   233       delete [] doubles;

   234     }

   235     /// \e

   236     virtual void eraseCol(const ColIt& col_it) = 0;

   237     /// \e

   238     virtual void eraseRow(const RowIt& row_it) = 0;

   239     /// \e

   240     virtual void setColBounds(const ColIt& col_it, int bound_type,

   241 			      _Value lo, _Value up) =0;

   242     /// \e

   243     virtual _Value getObjCoef(const ColIt& col_it) = 0;

   244     /// \e

   245     virtual void setRowBounds(const RowIt& row_it, int bound_type,

   246 			      _Value lo, _Value up) = 0;

   247     /// \e

   248     virtual void setObjCoef(const ColIt& col_it, _Value obj_coef) = 0;

   249     /// \e

   250     virtual void solveSimplex() = 0;

   251     /// \e

   252     virtual void solvePrimalSimplex() = 0;

   253     /// \e

   254     virtual void solveDualSimplex() = 0;

   255     /// \e

   256     virtual _Value getPrimal(const ColIt& col_it) = 0;

   257     /// \e

   258     virtual _Value getObjVal() = 0;

   259     /// \e

   260     virtual int rowNum() const = 0;

   261     /// \e

   262     virtual int colNum() const = 0;

   263     /// \e

   264     virtual int warmUp() = 0;

   265     /// \e

   266     virtual void printWarmUpStatus(int i) = 0;

   267     /// \e

   268     virtual int getPrimalStatus() = 0;

   269     /// \e

   270     virtual void printPrimalStatus(int i) = 0;

   271     /// \e

   272     virtual int getDualStatus() = 0;

   273     /// \e

   274     virtual void printDualStatus(int i) = 0;

   275     /// Returns the status of the slack variable assigned to row \c row_it.

   276     virtual int getRowStat(const RowIt& row_it) = 0;

   277     /// \e

   278     virtual void printRowStatus(int i) = 0;

   279     /// Returns the status of the variable assigned to column \c col_it.

   280     virtual int getColStat(const ColIt& col_it) = 0;

   281     /// \e

   282     virtual void printColStatus(int i) = 0;

   283   };

   286   /// \brief Wrappers for LP solvers

   287   ///

   288   /// This class implements a lemon wrapper for glpk.

   289   /// Later other LP-solvers will be wrapped into lemon.

   290   /// The aim of this class is to give a general surface to different

   291   /// solvers, i.e. it makes possible to write algorithms using LP's,

   292   /// in which the solver can be changed to an other one easily.

   293   class LPSolverWrapper : public LPSolverBase<double> {

   294   public:

   295     typedef LPSolverBase<double> Parent;

   297     //   class Row {

   298     //   protected:

   299     //     int i;

   300     //   public:

   301     //     Row() { }

   302     //     Row(const Invalid&) : i(0) { }

   303     //     Row(const int& _i) : i(_i) { }

   304     //     operator int() const { return i; }

   305     //   };

   306     //   class RowIt : public Row {

   307     //   public:

   308     //     RowIt(const Row& row) : Row(row) { }

   309     //   };

   311     //   class Col {

   312     //   protected:

   313     //     int i;

   314     //   public:

   315     //     Col() { }

   316     //     Col(const Invalid&) : i(0) { }

   317     //     Col(const int& _i) : i(_i) { }

   318     //     operator int() const { return i; }

   319     //   };

   320     //   class ColIt : public Col {

   321     //     ColIt(const Col& col) : Col(col) { }

   322     //   };

   324   public:

   325     /// \e

   326     LPX* lp;

   328   public:

   329     /// \e

   330     LPSolverWrapper() : Parent(),

   331 			lp(lpx_create_prob()) {

   332       lpx_set_int_parm(lp, LPX_K_DUAL, 1);

   333     }

   334     /// \e

   335     ~LPSolverWrapper() {

   336       lpx_delete_prob(lp);

   337     }

   338     /// \e

   339     void setMinimize() {

   340       lpx_set_obj_dir(lp, LPX_MIN);

   341     }

   342     /// \e

   343     void setMaximize() {

   344       lpx_set_obj_dir(lp, LPX_MAX);

   345     }

   346     /// \e

   347     ColIt addCol() {

   348       int i=lpx_add_cols(lp, 1);

   349       ColIt col_it;

   350       col_iter_map.first(col_it, INVALID_ID);

   351       if (col_iter_map.valid(col_it)) { //van hasznalhato hely

   352 	col_iter_map.set(col_it, INVALID_ID, VALID_ID);

   353 	col_iter_map[col_it]=i;

   354 	//col_id_to_lp_col_id[col_iter_map[col_it]]=i;

   355       } else { //a cucc vegere kell inzertalni mert nincs szabad hely

   356 	//col_id_to_lp_col_id.push_back(i);

   357 	//int j=col_id_to_lp_col_id.size()-1;

   358 	col_it=col_iter_map.push_back(i, VALID_ID);

   359       }

   360       //    edge_index_map.set(e, i);

   361       //    lpx_set_col_bnds(lp, i, LPX_DB, 0.0, 1.0);

   362       //    lpx_set_obj_coef(lp, i, cost[e]);

   363       return col_it;

   364     }

   365     /// \e

   366     RowIt addRow() {

   367       int i=lpx_add_rows(lp, 1);

   368       RowIt row_it;

   369       row_iter_map.first(row_it, INVALID_ID);

   370       if (row_iter_map.valid(row_it)) { //van hasznalhato hely

   371 	row_iter_map.set(row_it, INVALID_ID, VALID_ID);

   372 	row_iter_map[row_it]=i;

   373       } else { //a cucc vegere kell inzertalni mert nincs szabad hely

   374 	row_it=row_iter_map.push_back(i, VALID_ID);

   375       }

   376       return row_it;

   377     }

   378     using Parent::setRowCoeffs;

   379     void setRowCoeffs(RowIt row_it, int length,

   380 		      int* indices, double* doubles) {

   381       lpx_set_mat_row(lp, row_iter_map[row_it], length, indices, doubles);

   382     }

   383     using Parent::setColCoeffs;

   384     void setColCoeffs(ColIt col_it, int length,

   385 		      int* indices, double* doubles) {

   386       lpx_set_mat_col(lp, col_iter_map[col_it], length, indices, doubles);

   387     }

   388     //     //pair<RowIt, double>-bol kell megadni egy std range-et

   389     //     /// \e

   390     //     template <typename Begin, typename End>

   391     //     void setColCoeffs(const ColIt& col_it,

   392     // 		      Begin begin, End end) {

   393     //       int mem_length=1+lpx_get_num_rows(lp);

   394     //       int* indices = new int[mem_length];

   395     //       double* doubles = new double[mem_length];

   396     //       int length=0;

   397     //       for ( ; begin!=end; ++begin) {

   398     // 	++length;

   399     // 	indices[length]=row_iter_map[begin->first];

   400     // 	doubles[length]=begin->second;

   401     //       }

   402     //       lpx_set_mat_col(lp, col_iter_map[col_it], length, indices, doubles);

   403     //       delete [] indices;

   404     //       delete [] doubles;

   405     //     }

   406     //     //pair<ColIt, double>-bol kell megadni egy std range-et

   407     //     /// \e

   408     //     template <typename Begin, typename End>

   409     //     void setRowCoeffs(const RowIt& row_it,

   410     // 		      Begin begin, End end) {

   411     //       int mem_length=1+lpx_get_num_cols(lp);

   412     //       int* indices = new int[mem_length];

   413     //       double* doubles = new double[mem_length];

   414     //       int length=0;

   415     //       for ( ; begin!=end; ++begin) {

   416     // 	++length;

   417     // 	indices[length]=col_iter_map[begin->first];

   418     // 	doubles[length]=begin->second;

   419     //       }

   420     //       lpx_set_mat_row(lp, row_iter_map[row_it], length, indices, doubles);

   421     //       delete [] indices;

   422     //       delete [] doubles;

   423     //     }

   424     /// \e

   425     void eraseCol(const ColIt& col_it) {

   426       col_iter_map.set(col_it, VALID_ID, INVALID_ID);

   427       int cols[2];

   428       cols[1]=col_iter_map[col_it];

   429       lpx_del_cols(lp, 1, cols);

   430       col_iter_map[col_it]=0; //glpk specifikus

   431       ColIt it;

   432       for (col_iter_map.first(it, VALID_ID);

   433 	   col_iter_map.valid(it); col_iter_map.next(it)) {

   434 	if (col_iter_map[it]>cols[1]) --col_iter_map[it];

   435       }

   436     }

   437     /// \e

   438     void eraseRow(const RowIt& row_it) {

   439       row_iter_map.set(row_it, VALID_ID, INVALID_ID);

   440       int rows[2];

   441       rows[1]=row_iter_map[row_it];

   442       lpx_del_rows(lp, 1, rows);

   443       row_iter_map[row_it]=0; //glpk specifikus

   444       RowIt it;

   445       for (row_iter_map.first(it, VALID_ID);

   446 	   row_iter_map.valid(it); row_iter_map.next(it)) {

   447 	if (row_iter_map[it]>rows[1]) --row_iter_map[it];

   448       }

   449     }

   450     /// \e

   451     void setColBounds(const ColIt& col_it, int bound_type,

   452 		      double lo, double up) {

   453       lpx_set_col_bnds(lp, col_iter_map[col_it], bound_type, lo, up);

   454     }

   455     /// \e

   456     double getObjCoef(const ColIt& col_it) {

   457       return lpx_get_obj_coef(lp, col_iter_map[col_it]);

   458     }

   459     /// \e

   460     void setRowBounds(const RowIt& row_it, int bound_type,

   461 		      double lo, double up) {

   462       lpx_set_row_bnds(lp, row_iter_map[row_it], bound_type, lo, up);

   463     }

   464     /// \e

   465     void setObjCoef(const ColIt& col_it, double obj_coef) {

   466       lpx_set_obj_coef(lp, col_iter_map[col_it], obj_coef);

   467     }

   468     /// \e

   469     void solveSimplex() { lpx_simplex(lp); }

   470     /// \e

   471     void solvePrimalSimplex() { lpx_simplex(lp); }

   472     /// \e

   473     void solveDualSimplex() { lpx_simplex(lp); }

   474     /// \e

   475     double getPrimal(const ColIt& col_it) {

   476       return lpx_get_col_prim(lp, col_iter_map[col_it]);

   477     }

   478     /// \e

   479     double getObjVal() { return lpx_get_obj_val(lp); }

   480     /// \e

   481     int rowNum() const { return lpx_get_num_rows(lp); }

   482     /// \e

   483     int colNum() const { return lpx_get_num_cols(lp); }

   484     /// \e

   485     int warmUp() { return lpx_warm_up(lp); }

   486     /// \e

   487     void printWarmUpStatus(int i) {

   488       switch (i) {

   489       case LPX_E_OK: cout << "LPX_E_OK" << endl; break;

   490       case LPX_E_EMPTY: cout << "LPX_E_EMPTY" << endl; break;

   491       case LPX_E_BADB: cout << "LPX_E_BADB" << endl; break;

   492       case LPX_E_SING: cout << "LPX_E_SING" << endl; break;

   493       }

   494     }

   495     /// \e

   496     int getPrimalStatus() { return lpx_get_prim_stat(lp); }

   497     /// \e

   498     void printPrimalStatus(int i) {

   499       switch (i) {

   500       case LPX_P_UNDEF: cout << "LPX_P_UNDEF" << endl; break;

   501       case LPX_P_FEAS: cout << "LPX_P_FEAS" << endl; break;

   502       case LPX_P_INFEAS: cout << "LPX_P_INFEAS" << endl; break;

   503       case LPX_P_NOFEAS: cout << "LPX_P_NOFEAS" << endl; break;

   504       }

   505     }

   506     /// \e

   507     int getDualStatus() { return lpx_get_dual_stat(lp); }

   508     /// \e

   509     void printDualStatus(int i) {

   510       switch (i) {

   511       case LPX_D_UNDEF: cout << "LPX_D_UNDEF" << endl; break;

   512       case LPX_D_FEAS: cout << "LPX_D_FEAS" << endl; break;

   513       case LPX_D_INFEAS: cout << "LPX_D_INFEAS" << endl; break;

   514       case LPX_D_NOFEAS: cout << "LPX_D_NOFEAS" << endl; break;

   515       }

   516     }

   517     /// Returns the status of the slack variable assigned to row \c row_it.

   518     int getRowStat(const RowIt& row_it) {

   519       return lpx_get_row_stat(lp, row_iter_map[row_it]);

   520     }

   521     /// \e

   522     void printRowStatus(int i) {

   523       switch (i) {

   524       case LPX_BS: cout << "LPX_BS" << endl; break;

   525       case LPX_NL: cout << "LPX_NL" << endl; break;

   526       case LPX_NU: cout << "LPX_NU" << endl; break;

   527       case LPX_NF: cout << "LPX_NF" << endl; break;

   528       case LPX_NS: cout << "LPX_NS" << endl; break;

   529       }

   530     }

   531     /// Returns the status of the variable assigned to column \c col_it.

   532     int getColStat(const ColIt& col_it) {

   533       return lpx_get_col_stat(lp, col_iter_map[col_it]);

   534     }

   535     /// \e

   536     void printColStatus(int i) {

   537       switch (i) {

   538       case LPX_BS: cout << "LPX_BS" << endl; break;

   539       case LPX_NL: cout << "LPX_NL" << endl; break;

   540       case LPX_NU: cout << "LPX_NU" << endl; break;

   541       case LPX_NF: cout << "LPX_NF" << endl; break;

   542       case LPX_NS: cout << "LPX_NS" << endl; break;

   543       }

   544     }

   545   };

   547   /// @}

   549 } //namespace lemon

   551 #endif //LEMON_LP_SOLVER_WRAPPER_H

author	marci
	Tue, 04 Jan 2005 17:54:41 +0000
changeset 1048	38a49245a701
parent 1031	0b7169db694f
permissions	-rw-r--r--