src/work/marci/lp/lp_solver_wrapper.h
author marci
Tue, 31 Aug 2004 11:26:59 +0000
changeset 775 e46a1f0623a0
parent 765 4405b6be83bb
child 888 cc3590763f7f
permissions -rw-r--r--
ResGraphWrapper<Graph> is done, so does dimacs.h.
     1 // -*- c++ -*-
     2 #ifndef HUGO_LP_SOLVER_WRAPPER_H
     3 #define HUGO_LP_SOLVER_WRAPPER
     4 
     5 ///\ingroup misc
     6 ///\file
     7 ///\brief Dijkstra algorithm.
     8 
     9 // #include <stdio.h>
    10 #include <stdlib.h>
    11 // #include <stdio>
    12 //#include <stdlib>
    13 #include "glpk.h"
    14 
    15 #include <iostream>
    16 #include <vector>
    17 #include <string>
    18 #include <list>
    19 #include <memory>
    20 #include <utility>
    21 
    22 //#include <sage_graph.h>
    23 //#include <hugo/list_graph.h>
    24 //#include <hugo/graph_wrapper.h>
    25 #include <hugo/invalid.h>
    26 //#include <bfs_dfs.h>
    27 //#include <stp.h>
    28 //#include <hugo/max_flow.h>
    29 //#include <augmenting_flow.h>
    30 //#include <iter_map.h>
    31 
    32 using std::cout;
    33 using std::cin;
    34 using std::endl;
    35 
    36 namespace hugo {
    37 
    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 hugo 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 hugo style iterator iterates through a class. 
    67     class ClassIt;
    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(ClassIt 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(ClassIt 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     ClassIt 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(ClassIt 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[](ClassIt it) { return nodes[it.i].data; }
   128     /// Returns the data pointed by \c it.
   129     const T& operator[](ClassIt it) const { return nodes[it.i].data; }
   130     ///.
   131     class ClassIt {
   132       friend class IterablePartition;
   133     protected:
   134       int i;
   135     public:
   136       /// Default constructor.
   137       ClassIt() { }
   138       /// This constructor constructs an iterator which points
   139       /// to the member of th container indexed by the integer _i.
   140       ClassIt(const int& _i) : i(_i) { }
   141       /// Invalid constructor.
   142       ClassIt(const Invalid&) : i(-1) { }
   143     };
   144     /// First member of class \c class_id.
   145     ClassIt& first(ClassIt& it, int class_id) const {
   146       it.i=tips[class_id].first;
   147       return it;
   148     }
   149     /// Next member.
   150     ClassIt& next(ClassIt& it) const {
   151       it.i=nodes[it.i].next;
   152       return it;
   153     }
   154     /// True iff the iterator is valid.
   155     bool valid(const ClassIt& it) const { return it.i!=-1; }
   156   };
   157   
   158   /// \brief Wrappers for LP solvers
   159   /// 
   160   /// This class implements a hugo wrapper for glpk.
   161   /// Later other LP-solvers will be wrapped into hugo.
   162   /// The aim of this class is to give a general surface to different 
   163   /// solvers, i.e. it makes possible to write algorithms using LP's, 
   164   /// in which the solver can be changed to an other one easily.
   165   class LPSolverWrapper {
   166   public:
   167 
   168 //   class Row {
   169 //   protected:
   170 //     int i;
   171 //   public:
   172 //     Row() { }
   173 //     Row(const Invalid&) : i(0) { }
   174 //     Row(const int& _i) : i(_i) { }
   175 //     operator int() const { return i; }
   176 //   };
   177 //   class RowIt : public Row {
   178 //   public:
   179 //     RowIt(const Row& row) : Row(row) { }
   180 //   };
   181 
   182 //   class Col {
   183 //   protected:
   184 //     int i;
   185 //   public:
   186 //     Col() { }
   187 //     Col(const Invalid&) : i(0) { }
   188 //     Col(const int& _i) : i(_i) { }
   189 //     operator int() const { return i; }
   190 //   };
   191 //   class ColIt : public Col {
   192 //     ColIt(const Col& col) : Col(col) { }
   193 //   };
   194 
   195   public:
   196     ///.
   197     LPX* lp;
   198     ///.
   199     typedef IterablePartition<int>::ClassIt RowIt;
   200     ///.
   201     IterablePartition<int> row_iter_map;
   202     ///.
   203     typedef IterablePartition<int>::ClassIt ColIt;
   204     ///.
   205     IterablePartition<int> col_iter_map;
   206     //std::vector<int> row_id_to_lp_row_id;
   207     //std::vector<int> col_id_to_lp_col_id;
   208     ///.
   209     const int VALID_ID;
   210     ///.
   211     const int INVALID_ID;
   212 
   213   public:
   214     ///.
   215     LPSolverWrapper() : lp(lpx_create_prob()), 
   216 			row_iter_map(2), 
   217 			col_iter_map(2), 
   218 			//row_id_to_lp_row_id(), col_id_to_lp_col_id(), 
   219 			VALID_ID(0), INVALID_ID(1) {
   220       lpx_set_int_parm(lp, LPX_K_DUAL, 1);
   221     }
   222     ///.
   223     ~LPSolverWrapper() {
   224       lpx_delete_prob(lp);
   225     }
   226     ///.
   227     void setMinimize() { 
   228       lpx_set_obj_dir(lp, LPX_MIN);
   229     }
   230     ///.
   231     void setMaximize() { 
   232       lpx_set_obj_dir(lp, LPX_MAX);
   233     }
   234     ///.
   235     ColIt addCol() {
   236       int i=lpx_add_cols(lp, 1);  
   237       ColIt col_it;
   238       col_iter_map.first(col_it, INVALID_ID);
   239       if (col_iter_map.valid(col_it)) { //van hasznalhato hely
   240 	col_iter_map.set(col_it, INVALID_ID, VALID_ID);
   241 	col_iter_map[col_it]=i;
   242 	//col_id_to_lp_col_id[col_iter_map[col_it]]=i;
   243       } else { //a cucc vegere kell inzertalni mert nincs szabad hely
   244 	//col_id_to_lp_col_id.push_back(i);
   245 	//int j=col_id_to_lp_col_id.size()-1;
   246 	col_it=col_iter_map.push_back(i, VALID_ID);
   247       }
   248 //    edge_index_map.set(e, i);
   249 //    lpx_set_col_bnds(lp, i, LPX_DB, 0.0, 1.0);
   250 //    lpx_set_obj_coef(lp, i, cost[e]);    
   251       return col_it;
   252     }
   253     ///.
   254     RowIt addRow() {
   255       int i=lpx_add_rows(lp, 1);  
   256       RowIt row_it;
   257       row_iter_map.first(row_it, INVALID_ID);
   258       if (row_iter_map.valid(row_it)) { //van hasznalhato hely
   259 	row_iter_map.set(row_it, INVALID_ID, VALID_ID);
   260 	row_iter_map[row_it]=i;
   261       } else { //a cucc vegere kell inzertalni mert nincs szabad hely
   262 	row_it=row_iter_map.push_back(i, VALID_ID);
   263       }
   264       return row_it;
   265     }
   266     //pair<RowIt, double>-bol kell megadni egy std range-et
   267     ///.
   268     template <typename Begin, typename End>
   269     void setColCoeffs(const ColIt& col_it, 
   270 		      Begin begin, End end) {
   271       int mem_length=1+lpx_get_num_rows(lp);
   272       int* indices = new int[mem_length];
   273       double* doubles = new double[mem_length];
   274       int length=0;
   275       for ( ; begin!=end; ++begin) {
   276 	++length;
   277 	indices[length]=row_iter_map[begin->first];
   278 	doubles[length]=begin->second;
   279       }
   280       lpx_set_mat_col(lp, col_iter_map[col_it], length, indices, doubles);
   281       delete [] indices;
   282       delete [] doubles;
   283     }
   284     //pair<ColIt, double>-bol kell megadni egy std range-et
   285     ///.
   286     template <typename Begin, typename End>
   287     void setRowCoeffs(const RowIt& row_it, 
   288 		      Begin begin, End end) {
   289       int mem_length=1+lpx_get_num_cols(lp);
   290       int* indices = new int[mem_length];
   291       double* doubles = new double[mem_length];
   292       int length=0;
   293       for ( ; begin!=end; ++begin) {
   294 	++length;
   295 	indices[length]=col_iter_map[begin->first];
   296 	doubles[length]=begin->second;
   297       }
   298       lpx_set_mat_row(lp, row_iter_map[row_it], length, indices, doubles);
   299       delete [] indices;
   300       delete [] doubles;
   301     }
   302     ///.
   303     void eraseCol(const ColIt& col_it) {
   304       col_iter_map.set(col_it, VALID_ID, INVALID_ID);
   305       int cols[2];
   306       cols[1]=col_iter_map[col_it];
   307       lpx_del_cols(lp, 1, cols);
   308       col_iter_map[col_it]=0; //glpk specifikus
   309       ColIt it;
   310       for (col_iter_map.first(it, VALID_ID); 
   311 	   col_iter_map.valid(it); col_iter_map.next(it)) {
   312 	if (col_iter_map[it]>cols[1]) --col_iter_map[it];
   313       }
   314     }
   315     ///.
   316     void eraseRow(const RowIt& row_it) {
   317       row_iter_map.set(row_it, VALID_ID, INVALID_ID);
   318       int rows[2];
   319       rows[1]=row_iter_map[row_it];
   320       lpx_del_rows(lp, 1, rows);
   321       row_iter_map[row_it]=0; //glpk specifikus
   322       RowIt it;
   323       for (row_iter_map.first(it, VALID_ID); 
   324 	   row_iter_map.valid(it); row_iter_map.next(it)) {
   325 	if (row_iter_map[it]>rows[1]) --row_iter_map[it];
   326       }
   327     }
   328     ///.
   329     void setColBounds(const ColIt& col_it, int bound_type, 
   330 		      double lo, double up) {
   331       lpx_set_col_bnds(lp, col_iter_map[col_it], bound_type, lo, up);
   332     }
   333     ///.
   334     double getObjCoef(const ColIt& col_it) { 
   335       return lpx_get_obj_coef(lp, col_iter_map[col_it]);
   336     }
   337     ///.
   338     void setRowBounds(const RowIt& row_it, int bound_type, 
   339 		      double lo, double up) {
   340       lpx_set_row_bnds(lp, row_iter_map[row_it], bound_type, lo, up);
   341     }
   342 //   void setObjCoef(const RowIt& row_it, double obj_coef) { 
   343 //     lpx_set_obj_coef(lp, row_iter_map[row_it], obj_coef);
   344 //   }
   345     ///.
   346     void solveSimplex() { lpx_simplex(lp); }
   347     ///.
   348     void solvePrimalSimplex() { lpx_simplex(lp); }
   349     ///.
   350     void solveDualSimplex() { lpx_simplex(lp); }
   351     ///.
   352     double getPrimal(const ColIt& col_it) {
   353       return lpx_get_col_prim(lp, col_iter_map[col_it]);
   354     }
   355     ///.
   356     double getObjVal() { return lpx_get_obj_val(lp); }
   357     ///.
   358     int rowNum() const { return lpx_get_num_rows(lp); }
   359     ///.
   360     int colNum() const { return lpx_get_num_cols(lp); }
   361     ///.
   362     int warmUp() { return lpx_warm_up(lp); }
   363     ///.
   364     void printWarmUpStatus(int i) {
   365       switch (i) {
   366 	case LPX_E_OK: cout << "LPX_E_OK" << endl; break;
   367 	case LPX_E_EMPTY: cout << "LPX_E_EMPTY" << endl; break;	
   368 	case LPX_E_BADB: cout << "LPX_E_BADB" << endl; break;
   369 	case LPX_E_SING: cout << "LPX_E_SING" << endl; break;
   370       }
   371     }
   372     ///.
   373     int getPrimalStatus() { return lpx_get_prim_stat(lp); }
   374     ///.
   375     void printPrimalStatus(int i) {
   376       switch (i) {
   377 	case LPX_P_UNDEF: cout << "LPX_P_UNDEF" << endl; break;
   378 	case LPX_P_FEAS: cout << "LPX_P_FEAS" << endl; break;	
   379 	case LPX_P_INFEAS: cout << "LPX_P_INFEAS" << endl; break;
   380 	case LPX_P_NOFEAS: cout << "LPX_P_NOFEAS" << endl; break;
   381       }
   382     }
   383     ///.
   384     int getDualStatus() { return lpx_get_dual_stat(lp); }
   385     ///.
   386     void printDualStatus(int i) {
   387       switch (i) {
   388 	case LPX_D_UNDEF: cout << "LPX_D_UNDEF" << endl; break;
   389 	case LPX_D_FEAS: cout << "LPX_D_FEAS" << endl; break;	
   390 	case LPX_D_INFEAS: cout << "LPX_D_INFEAS" << endl; break;
   391 	case LPX_D_NOFEAS: cout << "LPX_D_NOFEAS" << endl; break;
   392       }
   393     }
   394     /// Returns the status of the slack variable assigned to row \c row_it.
   395     int getRowStat(const RowIt& row_it) { 
   396       return lpx_get_row_stat(lp, row_iter_map[row_it]); 
   397     }
   398     ///.
   399     void printRowStatus(int i) {
   400       switch (i) {
   401 	case LPX_BS: cout << "LPX_BS" << endl; break;
   402 	case LPX_NL: cout << "LPX_NL" << endl; break;	
   403 	case LPX_NU: cout << "LPX_NU" << endl; break;
   404 	case LPX_NF: cout << "LPX_NF" << endl; break;
   405 	case LPX_NS: cout << "LPX_NS" << endl; break;
   406       }
   407     }
   408     /// Returns the status of the variable assigned to column \c col_it.
   409     int getColStat(const ColIt& col_it) { 
   410       return lpx_get_col_stat(lp, col_iter_map[col_it]); 
   411     }
   412     ///.
   413     void printColStatus(int i) {
   414       switch (i) {
   415 	case LPX_BS: cout << "LPX_BS" << endl; break;
   416 	case LPX_NL: cout << "LPX_NL" << endl; break;	
   417 	case LPX_NU: cout << "LPX_NU" << endl; break;
   418 	case LPX_NF: cout << "LPX_NF" << endl; break;
   419 	case LPX_NS: cout << "LPX_NS" << endl; break;
   420       }
   421     }
   422   };
   423   
   424   /// @}
   425 
   426 } //namespace hugo
   427 
   428 #endif //HUGO_LP_SOLVER_WRAPPER_H