lemon/lp_glpk.cc
author kpeter
Thu, 17 Apr 2008 21:46:06 +0000
changeset 2607 78e8de179fe2
parent 2591 3b4d5bc3b4fb
child 2622 fa2877651022
permissions -rw-r--r--
Remove SspMinCostFlow, since it is fully replaced by other classes.
     1 /* -*- C++ -*-
     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 ///\file
    20 ///\brief Implementation of the LEMON-GLPK lp solver interface.
    21 
    22 #include <lemon/lp_glpk.h>
    23 //#include <iostream>
    24 
    25 #if GLP_MAJOR_VERSION > 4 || (GLP_MAJOR_VERSION == 4 && GLP_MINOR_VERSION > 15)
    26 #define LEMON_glp(func) (glp_##func)
    27 #define LEMON_lpx(func) (lpx_##func)
    28 
    29 #define LEMON_GLP(def) (GLP_##def)
    30 #define LEMON_LPX(def) (LPX_##def)
    31 
    32 #else
    33 
    34 #define LEMON_glp(func) (lpx_##func)
    35 #define LEMON_lpx(func) (lpx_##func)
    36 
    37 #define LEMON_GLP(def) (LPX_##def)
    38 #define LEMON_LPX(def) (LPX_##def)
    39 
    40 #endif
    41 
    42 namespace lemon {
    43 
    44   LpGlpk::LpGlpk() : Parent() {
    45     solved = false;
    46     rows = _lp_bits::LpId(1);
    47     cols = _lp_bits::LpId(1);
    48     lp = LEMON_glp(create_prob)();
    49     LEMON_glp(create_index)(lp);
    50     LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_DUAL), 1);
    51     messageLevel(0);
    52   }
    53   
    54   LpGlpk::LpGlpk(const LpGlpk &glp) : Parent() {
    55     solved = false;
    56     rows = _lp_bits::LpId(1);
    57     cols = _lp_bits::LpId(1);
    58     lp = LEMON_glp(create_prob)();
    59     LEMON_glp(create_index)(lp);
    60     ///\todo control function for this:
    61     LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_DUAL), 1);
    62     messageLevel(0);
    63     //Coefficient matrix, row bounds
    64     LEMON_glp(add_rows)(lp, LEMON_glp(get_num_rows)(glp.lp));
    65     LEMON_glp(add_cols)(lp, LEMON_glp(get_num_cols)(glp.lp));
    66     int len;
    67     std::vector<int> ind(1+LEMON_glp(get_num_cols)(glp.lp));
    68     std::vector<Value> val(1+LEMON_glp(get_num_cols)(glp.lp));
    69     for (int i=1;i<=LEMON_glp(get_num_rows)(glp.lp);++i)
    70       {
    71 	len=LEMON_glp(get_mat_row)(glp.lp,i,&*ind.begin(),&*val.begin());
    72 	LEMON_glp(set_mat_row)(lp, i,len,&*ind.begin(),&*val.begin());
    73 	LEMON_glp(set_row_bnds)(lp,i,
    74 				LEMON_glp(get_row_type)(glp.lp,i),
    75 				LEMON_glp(get_row_lb)(glp.lp,i),
    76 				LEMON_glp(get_row_ub)(glp.lp,i));
    77       }
    78 
    79     //Objective function, coloumn bounds
    80     LEMON_glp(set_obj_dir)(lp, LEMON_glp(get_obj_dir)(glp.lp));
    81     //Objectif function's constant term treated separately
    82     LEMON_glp(set_obj_coef)(lp,0,LEMON_glp(get_obj_coef)(glp.lp,0));
    83     for (int i=1;i<=LEMON_glp(get_num_cols)(glp.lp);++i)
    84       {
    85 	LEMON_glp(set_obj_coef)(lp,i,
    86 				LEMON_glp(get_obj_coef)(glp.lp,i));
    87 	LEMON_glp(set_col_bnds)(lp,i,
    88 				LEMON_glp(get_col_type)(glp.lp,i),
    89 				LEMON_glp(get_col_lb)(glp.lp,i),
    90 				LEMON_glp(get_col_ub)(glp.lp,i));
    91       }
    92     rows = glp.rows;
    93     cols = glp.cols;
    94   }
    95   
    96   LpGlpk::~LpGlpk() {
    97     LEMON_glp(delete_prob)(lp);
    98   }
    99   
   100   int LpGlpk::_addCol() { 
   101     int i=LEMON_glp(add_cols)(lp, 1);
   102     LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), 0.0, 0.0);
   103     solved = false;
   104     return i;
   105   }
   106 
   107   ///\e
   108 
   109 
   110   LpSolverBase* LpGlpk::_newLp()
   111   {
   112     LpGlpk* newlp = new LpGlpk;
   113     return newlp;
   114   }
   115   
   116   ///\e
   117 
   118   LpSolverBase* LpGlpk::_copyLp()
   119   {
   120     LpGlpk *newlp = new LpGlpk(*this);
   121     return newlp;
   122   }
   123 
   124   int LpGlpk::_addRow() { 
   125     int i=LEMON_glp(add_rows)(lp, 1);
   126     solved = false;
   127     return i;
   128   }
   129 
   130   
   131   void LpGlpk::_eraseCol(int i) {
   132     int ca[2];
   133     ca[1]=i;
   134     LEMON_glp(del_cols)(lp, 1, ca);
   135     solved = false;
   136   }
   137   
   138   void LpGlpk::_eraseRow(int i) {
   139     int ra[2];
   140     ra[1]=i;
   141     LEMON_glp(del_rows)(lp, 1, ra);
   142     solved = false;
   143   }
   144 
   145   void LpGlpk::_getColName(int c, std::string & name) const
   146   {
   147     
   148     const char *n = LEMON_glp(get_col_name)(lp,c);
   149     name = n?n:"";
   150   }
   151   
   152   
   153   void LpGlpk::_setColName(int c, const std::string & name)
   154   {
   155     LEMON_glp(set_col_name)(lp,c,const_cast<char*>(name.c_str()));
   156 
   157   }
   158 
   159   int LpGlpk::_colByName(const std::string& name) const
   160   {
   161     int k = LEMON_glp(find_col)(lp, const_cast<char*>(name.c_str()));
   162     return k > 0 ? k : -1; 
   163   }
   164 
   165   
   166   void LpGlpk::_setRowCoeffs(int i, ConstRowIterator b, ConstRowIterator e) 
   167   {
   168     std::vector<int> indices;
   169     std::vector<Value> values;
   170 
   171     indices.push_back(0);
   172     values.push_back(0);
   173 
   174     for(ConstRowIterator it=b; it!=e; ++it) {
   175       indices.push_back(it->first);
   176       values.push_back(it->second);
   177     }
   178 
   179     LEMON_glp(set_mat_row)(lp, i, values.size() - 1, 
   180 				&indices[0], &values[0]);
   181 
   182     solved = false;
   183   }
   184 
   185   void LpGlpk::_getRowCoeffs(int ix, RowIterator b) const
   186   {
   187     int length = LEMON_glp(get_mat_row)(lp, ix, 0, 0);
   188     
   189     std::vector<int> indices(length + 1);
   190     std::vector<Value> values(length + 1);
   191     
   192     LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);
   193     
   194     for (int i = 1; i <= length; ++i) {
   195       *b = std::make_pair(indices[i], values[i]);
   196       ++b;
   197     }
   198   }
   199   
   200   void LpGlpk::_setColCoeffs(int ix, ConstColIterator b, ConstColIterator e) {
   201 
   202     std::vector<int> indices;
   203     std::vector<Value> values;
   204 
   205     indices.push_back(0);
   206     values.push_back(0);
   207 
   208     for(ConstColIterator it=b; it!=e; ++it) {
   209       indices.push_back(it->first);
   210       values.push_back(it->second);
   211     }
   212     
   213     LEMON_glp(set_mat_col)(lp, ix, values.size() - 1, 
   214 				&indices[0], &values[0]);
   215 
   216     solved = false;
   217   }
   218 
   219   void LpGlpk::_getColCoeffs(int ix, ColIterator b) const
   220   {
   221     int length = LEMON_glp(get_mat_col)(lp, ix, 0, 0);
   222     
   223     std::vector<int> indices(length + 1);
   224     std::vector<Value> values(length + 1);
   225     
   226     LEMON_glp(get_mat_col)(lp, ix, &indices[0], &values[0]);
   227     
   228     for (int i = 1; i <= length; ++i) {
   229       *b = std::make_pair(indices[i], values[i]);
   230       ++b;
   231     }
   232   }
   233 
   234   void LpGlpk::_setCoeff(int ix, int jx, Value value) 
   235   {
   236 
   237     if (LEMON_glp(get_num_cols)(lp) < LEMON_glp(get_num_rows)(lp)) {
   238 
   239       int length=LEMON_glp(get_mat_row)(lp, ix, 0, 0);
   240       
   241       std::vector<int> indices(length + 2);
   242       std::vector<Value> values(length + 2);
   243       
   244       LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);
   245       
   246       //The following code does not suppose that the elements of the
   247       //array indices are sorted
   248       bool found=false;
   249       for (int i = 1; i <= length; ++i) {
   250         if (indices[i]==jx){
   251           found=true;
   252           values[i]=value;
   253           break;
   254         }
   255       }
   256       if (!found){
   257         ++length;
   258         indices[length]=jx;
   259         values[length]=value;
   260       }
   261     
   262       LEMON_glp(set_mat_row)(lp, ix, length, &indices[0], &values[0]);
   263 
   264     } else {
   265 
   266       int length=LEMON_glp(get_mat_col)(lp, jx, 0, 0);
   267       
   268       std::vector<int> indices(length + 2);
   269       std::vector<Value> values(length + 2);
   270       
   271       LEMON_glp(get_mat_col)(lp, jx, &indices[0], &values[0]);
   272       
   273       //The following code does not suppose that the elements of the
   274       //array indices are sorted
   275       bool found=false;
   276       for (int i = 1; i <= length; ++i) {
   277         if (indices[i]==jx){
   278           found=true;
   279           values[i]=value;
   280           break;
   281         }
   282       }
   283       if (!found){
   284         ++length;
   285         indices[length]=ix;
   286         values[length]=value;
   287       }
   288     
   289       LEMON_glp(set_mat_col)(lp, jx, length, &indices[0], &values[0]);
   290     }
   291 
   292     solved = false;
   293   }
   294 
   295   LpGlpk::Value LpGlpk::_getCoeff(int ix, int jx) const
   296   {
   297 
   298     int length=LEMON_glp(get_mat_row)(lp, ix, 0, 0);
   299     
   300     std::vector<int> indices(length + 1);
   301     std::vector<Value> values(length + 1);
   302     
   303     LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);
   304     
   305     //The following code does not suppose that the elements of the
   306     //array indices are sorted
   307     for (int i = 1; i <= length; ++i) {
   308       if (indices[i]==jx){
   309 	return values[i];
   310       }
   311     }
   312     return 0;
   313 
   314   }
   315 
   316 
   317   void LpGlpk::_setColLowerBound(int i, Value lo)
   318   {
   319     if (lo==INF) {
   320       //FIXME error
   321     }
   322     int b=LEMON_glp(get_col_type)(lp, i);
   323     double up=LEMON_glp(get_col_ub)(lp, i);	
   324     if (lo==-INF) {
   325       switch (b) {
   326       case LEMON_GLP(FR):
   327       case LEMON_GLP(LO):
   328 	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), lo, up);
   329 	break;
   330       case LEMON_GLP(UP):
   331 	break;
   332       case LEMON_GLP(DB):
   333       case LEMON_GLP(FX):
   334 	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);
   335 	break;
   336       default: ;
   337 	//FIXME error
   338       }
   339     } else {
   340       switch (b) {
   341       case LEMON_GLP(FR):
   342       case LEMON_GLP(LO):
   343 	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(LO), lo, up);
   344 	break;
   345       case LEMON_GLP(UP):	  
   346       case LEMON_GLP(DB):
   347       case LEMON_GLP(FX):
   348 	if (lo==up) 
   349 	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FX), lo, up);
   350 	else 
   351 	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(DB), lo, up);
   352 	break;
   353       default: ;
   354 	//FIXME error
   355       }
   356     }
   357 
   358     solved = false;
   359   }
   360 
   361   LpGlpk::Value LpGlpk::_getColLowerBound(int i) const
   362   {
   363     int b=LEMON_glp(get_col_type)(lp, i);
   364       switch (b) {
   365       case LEMON_GLP(LO):
   366       case LEMON_GLP(DB):
   367       case LEMON_GLP(FX):
   368 	return LEMON_glp(get_col_lb)(lp, i);	
   369       default: ;
   370 	return -INF;
   371       }
   372   }
   373   
   374   void LpGlpk::_setColUpperBound(int i, Value up)
   375   {
   376     if (up==-INF) {
   377       //FIXME error
   378     }
   379     int b=LEMON_glp(get_col_type)(lp, i);
   380     double lo=LEMON_glp(get_col_lb)(lp, i);
   381     if (up==INF) {
   382       switch (b) {
   383       case LEMON_GLP(FR):
   384       case LEMON_GLP(LO):
   385 	break;
   386       case LEMON_GLP(UP):
   387 	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), lo, up);
   388 	break;
   389       case LEMON_GLP(DB):
   390       case LEMON_GLP(FX):
   391 	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(LO), lo, up);
   392 	break;
   393       default: ;
   394 	//FIXME error
   395       }
   396     } else {
   397       switch (b) {
   398       case LEMON_GLP(FR):
   399 	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);
   400 	break;
   401       case LEMON_GLP(UP):
   402 	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);
   403 	break;
   404       case LEMON_GLP(LO):
   405       case LEMON_GLP(DB):
   406       case LEMON_GLP(FX):
   407 	if (lo==up) 
   408 	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FX), lo, up);
   409 	else 
   410 	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(DB), lo, up);
   411 	break;
   412       default: ;
   413 	//FIXME error
   414       }
   415     }
   416 
   417     solved = false;
   418   }
   419 
   420   LpGlpk::Value LpGlpk::_getColUpperBound(int i) const
   421   {
   422     int b=LEMON_glp(get_col_type)(lp, i);
   423       switch (b) {
   424       case LEMON_GLP(UP):
   425       case LEMON_GLP(DB):
   426       case LEMON_GLP(FX):
   427 	return LEMON_glp(get_col_ub)(lp, i);	
   428       default: ;
   429 	return INF;
   430       }
   431   }
   432   
   433   void LpGlpk::_setRowBounds(int i, Value lb, Value ub)
   434   {
   435     //Bad parameter
   436     if (lb==INF || ub==-INF) {
   437       //FIXME error
   438     }
   439 
   440     if (lb == -INF){
   441       if (ub == INF){
   442 	LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(FR), lb, ub);
   443       }
   444       else{
   445 	LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(UP), lb, ub);
   446       }
   447     }
   448     else{
   449       if (ub==INF){
   450 	LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(LO), lb, ub);
   451 
   452       }
   453       else{
   454 	if (lb == ub){
   455 	  LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(FX), lb, ub);
   456 	}
   457 	else{
   458 	  LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(DB), lb, ub);
   459 	}
   460       }
   461     }
   462 
   463     solved = false;
   464   }
   465 
   466   void LpGlpk::_getRowBounds(int i, Value &lb, Value &ub) const
   467   {
   468 
   469     int b=LEMON_glp(get_row_type)(lp, i);
   470     switch (b) {
   471     case LEMON_GLP(FR):
   472     case LEMON_GLP(UP):
   473       lb = -INF;
   474 	break;
   475     default: 
   476       lb=LEMON_glp(get_row_lb)(lp, i);
   477     }
   478 
   479     switch (b) {
   480     case LEMON_GLP(FR):
   481     case LEMON_GLP(LO):
   482       ub = INF;
   483 	break;
   484     default: 
   485       ub=LEMON_glp(get_row_ub)(lp, i);
   486     }
   487     
   488   }
   489   
   490   void LpGlpk::_setObjCoeff(int i, Value obj_coef)
   491   {
   492     //i=0 means the constant term (shift)
   493     LEMON_glp(set_obj_coef)(lp, i, obj_coef);
   494 
   495     solved = false;
   496   }
   497 
   498   LpGlpk::Value LpGlpk::_getObjCoeff(int i) const {
   499     //i=0 means the constant term (shift)
   500     return LEMON_glp(get_obj_coef)(lp, i);
   501   }
   502 
   503   void LpGlpk::_clearObj()
   504   {
   505     for (int i=0;i<=LEMON_glp(get_num_cols)(lp);++i){
   506       LEMON_glp(set_obj_coef)(lp, i, 0);
   507     }
   508 
   509     solved = false;
   510   }
   511 
   512   LpGlpk::SolveExitStatus LpGlpk::_solve()
   513   {
   514     // A way to check the problem to be solved
   515     //LEMON_glp(write_cpxlp(lp,"naittvan.cpx");    
   516 
   517     LEMON_lpx(std_basis)(lp);
   518     int i =  LEMON_lpx(simplex)(lp);
   519     
   520     switch (i) {
   521     case LEMON_LPX(E_OK): 
   522       solved = true;
   523       return SOLVED;
   524     default:
   525       return UNSOLVED;
   526     }
   527   }
   528 
   529   LpGlpk::Value LpGlpk::_getPrimal(int i) const
   530   {
   531     return LEMON_glp(get_col_prim)(lp,i);
   532   }
   533 
   534   LpGlpk::Value LpGlpk::_getDual(int i) const
   535   {
   536     return LEMON_glp(get_row_dual)(lp,i);
   537   }
   538   
   539   LpGlpk::Value LpGlpk::_getPrimalValue() const
   540   {
   541     return LEMON_glp(get_obj_val)(lp);
   542   }
   543   bool LpGlpk::_isBasicCol(int i) const
   544   {
   545     return (LEMON_glp(get_col_stat)(lp, i)==LEMON_GLP(BS));
   546   }
   547   
   548  
   549   LpGlpk::SolutionStatus LpGlpk::_getPrimalStatus() const
   550   {
   551     if (!solved) return UNDEFINED;
   552     int stat=  LEMON_lpx(get_status)(lp);
   553     switch (stat) {
   554     case LEMON_LPX(UNDEF)://Undefined (no solve has been run yet)
   555       return UNDEFINED;
   556     case LEMON_LPX(NOFEAS)://There is no feasible solution (primal, I guess)
   557     case LEMON_LPX(INFEAS)://Infeasible 
   558       return INFEASIBLE;
   559     case LEMON_LPX(UNBND)://Unbounded
   560       return INFINITE;
   561     case LEMON_LPX(FEAS)://Feasible
   562       return FEASIBLE;
   563     case LEMON_LPX(OPT)://Feasible
   564       return OPTIMAL;
   565     default:
   566       return UNDEFINED; //to avoid gcc warning
   567       //FIXME error
   568     }
   569   }
   570 
   571   LpGlpk::SolutionStatus LpGlpk::_getDualStatus() const
   572   {
   573     if (!solved) return UNDEFINED;
   574     switch (LEMON_lpx(get_dual_stat)(lp)) {
   575     case LEMON_LPX(D_UNDEF)://Undefined (no solve has been run yet)
   576       return UNDEFINED;
   577     case LEMON_LPX(D_NOFEAS)://There is no dual feasible solution 
   578 //    case LEMON_LPX(D_INFEAS://Infeasible 
   579       return INFEASIBLE;
   580     case LEMON_LPX(D_FEAS)://Feasible    
   581       switch (LEMON_lpx(get_status)(lp)) {
   582       case LEMON_LPX(NOFEAS):
   583 	return INFINITE;
   584       case LEMON_LPX(OPT):
   585 	return OPTIMAL;
   586       default:
   587 	return FEASIBLE;
   588       }
   589     default:
   590       return UNDEFINED; //to avoid gcc warning
   591       //FIXME error
   592     }
   593   }
   594 
   595   LpGlpk::ProblemTypes LpGlpk::_getProblemType() const
   596   {
   597     if (!solved) return UNKNOWN;
   598       //int stat=  LEMON_glp(get_status(lp);
   599     int statp=  LEMON_lpx(get_prim_stat)(lp);
   600     int statd=  LEMON_lpx(get_dual_stat)(lp);
   601     if (statp==LEMON_LPX(P_FEAS) && statd==LEMON_LPX(D_FEAS))
   602 	return PRIMAL_DUAL_FEASIBLE;
   603     if (statp==LEMON_LPX(P_FEAS) && statd==LEMON_LPX(D_NOFEAS))
   604 	return PRIMAL_FEASIBLE_DUAL_INFEASIBLE;
   605     if (statp==LEMON_LPX(P_NOFEAS) && statd==LEMON_LPX(D_FEAS))
   606 	return PRIMAL_INFEASIBLE_DUAL_FEASIBLE;
   607     if (statp==LEMON_LPX(P_NOFEAS) && statd==LEMON_LPX(D_NOFEAS))
   608 	return PRIMAL_DUAL_INFEASIBLE;
   609     //In all other cases
   610     return UNKNOWN;
   611   }
   612 
   613   void LpGlpk::_setMax()
   614   {
   615     solved = false;
   616     LEMON_glp(set_obj_dir)(lp, LEMON_GLP(MAX));
   617   }
   618 
   619   void LpGlpk::_setMin()
   620   {
   621     solved = false;
   622     LEMON_glp(set_obj_dir)(lp, LEMON_GLP(MIN));
   623   }
   624 
   625   bool LpGlpk::_isMax() const
   626   {
   627     return (LEMON_glp(get_obj_dir)(lp)==LEMON_GLP(MAX));
   628   }
   629 
   630  
   631 
   632   void LpGlpk::messageLevel(int m)
   633   {
   634     LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_MSGLEV), m);
   635   }
   636 
   637   void LpGlpk::presolver(bool b)
   638   {
   639     LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_PRESOL), b);
   640   }
   641 
   642  
   643 } //END OF NAMESPACE LEMON