source: lemon-0.x/lemon/lp_glpk.cc @ 2323:8b18b6fed090

/* -*- C++ -*-
 *
 * This file is a part of LEMON, a generic C++ optimization library
 *
 * Copyright (C) 2003-2006
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
 *
 * Permission to use, modify and distribute this software is granted
 * provided that this copyright notice appears in all copies. For
 * precise terms see the accompanying LICENSE file.
 *
 * This software is provided "AS IS" with no warranty of any kind,
 * express or implied, and with no claim as to its suitability for any
 * purpose.
 *
 */

///\file
///\brief Implementation of the LEMON-GLPK lp solver interface.

#include <lemon/lp_glpk.h>
//#include <iostream>
namespace lemon {


  LpGlpk::LpGlpk() : Parent(), 
                     lp(lpx_create_prob()) {
    ///\todo constrol function for this:
    lpx_set_int_parm(lp, LPX_K_DUAL, 1);
    messageLevel(0);
  }
  
  LpGlpk::LpGlpk(const LpGlpk &glp) : Parent(glp), 
                                      lp(lpx_create_prob()) {
    ///\todo constrol function for this:
    lpx_set_int_parm(lp, LPX_K_DUAL, 1);
    messageLevel(0);
    //Coefficient matrix, row bounds
    lpx_add_rows(lp, lpx_get_num_rows(glp.lp));
    lpx_add_cols(lp, lpx_get_num_cols(glp.lp));
    int len;
    int ind[1+lpx_get_num_cols(glp.lp)];
    Value val[1+lpx_get_num_cols(glp.lp)];
    for (int i=1;i<=lpx_get_num_rows(glp.lp);++i)
      {
        len=lpx_get_mat_row(glp.lp,i,ind,val);
        lpx_set_mat_row(lp, i,len,ind,val);
        lpx_set_row_bnds(lp,i,lpx_get_row_type(glp.lp,i),
                         lpx_get_row_lb(glp.lp,i),lpx_get_row_ub(glp.lp,i));
      }

    //Objective function, coloumn bounds
    lpx_set_obj_dir(lp, lpx_get_obj_dir(glp.lp));
    //Objectif function's constant term treated separately
    lpx_set_obj_coef(lp,0,lpx_get_obj_coef(glp.lp,0));
    for (int i=1;i<=lpx_get_num_cols(glp.lp);++i)
      {
        lpx_set_obj_coef(lp,i,lpx_get_obj_coef(glp.lp,i));
        lpx_set_col_bnds(lp,i,lpx_get_col_type(glp.lp,i),
                         lpx_get_col_lb(glp.lp,i),lpx_get_col_ub(glp.lp,i));
      }
  }
  
  LpGlpk::~LpGlpk() {
    lpx_delete_prob(lp);
  }
  
  int LpGlpk::_addCol() { 
    int i=lpx_add_cols(lp, 1);
    _setColLowerBound(i, -INF);
    _setColUpperBound(i, INF);
    return i;
  }

  ///\e


  LpSolverBase &LpGlpk::_newLp()
  {
    LpGlpk* newlp=new LpGlpk;
    return *newlp;
  }
  
  ///\e

  LpSolverBase &LpGlpk::_copyLp()
  {
    LpGlpk* newlp=new LpGlpk(*this);
    return *newlp;
  }

  int LpGlpk::_addRow() { 
    int i=lpx_add_rows(lp, 1);
    return i;
  }

  
  void LpGlpk::_eraseCol(int i) {
    int cols[2];
    cols[1]=i;
    lpx_del_cols(lp, 1, cols);
  }
  
  void LpGlpk::_eraseRow(int i) {
    int rows[2];
    rows[1]=i;
    lpx_del_rows(lp, 1, rows);
  }

  void LpGlpk::_getColName(int col, std::string & name)
  {
    
    char *n = lpx_get_col_name(lp,col);
    name = n?n:"";
  }
  
  
  void LpGlpk::_setColName(int col, const std::string & name)
  {
    lpx_set_col_name(lp,col,const_cast<char*>(name.c_str()));
  }
  
  void LpGlpk::_setRowCoeffs(int i, LpRowIterator b, LpRowIterator e) 
  {
    std::vector<int> indices;
    std::vector<Value> values;

    indices.push_back(0);
    values.push_back(0);

    for(LpRowIterator it=b; it!=e; ++it) {
      indices.push_back(it->first);
      values.push_back(it->second);
    }

    lpx_set_mat_row(lp, i, values.size() - 1, &indices[0], &values[0]);
  }
  
  void LpGlpk::_setColCoeffs(int i, LpColIterator b, LpColIterator e) {

    std::vector<int> indices;
    std::vector<Value> values;

    indices.push_back(0);
    values.push_back(0);

    for(LpColIterator it=b; it!=e; ++it) {
      indices.push_back(it->first);
      values.push_back(it->second);
    }
    
    lpx_set_mat_col(lp, i, values.size() - 1, &indices[0], &values[0]);
  }


  void LpGlpk::_setCoeff(int row, int col, Value value) 
  {

    if (lpx_get_num_cols(lp) < lpx_get_num_rows(lp)) {

      int length=lpx_get_mat_row(lp, row, 0, 0);
      
      std::vector<int> indices(length + 2);
      std::vector<Value> values(length + 2);
      
      lpx_get_mat_row(lp, row, &indices[0], &values[0]);
      
      //The following code does not suppose that the elements of the
      //array indices are sorted
      bool found=false;
      for (int i = 1; i <= length; ++i) {
        if (indices[i]==col){
          found=true;
          values[i]=value;
          break;
        }
      }
      if (!found){
        ++length;
        indices[length]=col;
        values[length]=value;
      }
    
      lpx_set_mat_row(lp, row, length, &indices[0], &values[0]);

    } else {

      int length=lpx_get_mat_col(lp, col, 0, 0);
      
      std::vector<int> indices(length + 2);
      std::vector<Value> values(length + 2);
      
      lpx_get_mat_col(lp, col, &indices[0], &values[0]);
      
      //The following code does not suppose that the elements of the
      //array indices are sorted
      bool found=false;
      for (int i = 1; i <= length; ++i) {
        if (indices[i]==col){
          found=true;
          values[i]=value;
          break;
        }
      }
      if (!found){
        ++length;
        indices[length]=row;
        values[length]=value;
      }
    
      lpx_set_mat_col(lp, col, length, &indices[0], &values[0]);
    }
  }

  void LpGlpk::_setColLowerBound(int i, Value lo)
  {
    if (lo==INF) {
      //FIXME error
    }
    int b=lpx_get_col_type(lp, i);
    double up=lpx_get_col_ub(lp, i);    
    if (lo==-INF) {
      switch (b) {
      case LPX_FR:
      case LPX_LO:
        lpx_set_col_bnds(lp, i, LPX_FR, lo, up);
        break;
      case LPX_UP:
        break;
      case LPX_DB:
      case LPX_FX:
        lpx_set_col_bnds(lp, i, LPX_UP, lo, up);
        break;
      default: ;
        //FIXME error
      }
    } else {
      switch (b) {
      case LPX_FR:
      case LPX_LO:
        lpx_set_col_bnds(lp, i, LPX_LO, lo, up);
        break;
      case LPX_UP:        
      case LPX_DB:
      case LPX_FX:
        if (lo==up) 
          lpx_set_col_bnds(lp, i, LPX_FX, lo, up);
        else 
          lpx_set_col_bnds(lp, i, LPX_DB, lo, up);
        break;
      default: ;
        //FIXME error
      }
    }

  }
  
  void LpGlpk::_setColUpperBound(int i, Value up)
  {
    if (up==-INF) {
      //FIXME error
    }
    int b=lpx_get_col_type(lp, i);
    double lo=lpx_get_col_lb(lp, i);
    if (up==INF) {
      switch (b) {
      case LPX_FR:
      case LPX_LO:
        break;
      case LPX_UP:
        lpx_set_col_bnds(lp, i, LPX_FR, lo, up);
        break;
      case LPX_DB:
      case LPX_FX:
        lpx_set_col_bnds(lp, i, LPX_LO, lo, up);
        break;
      default: ;
        //FIXME error
      }
    } else {
      switch (b) {
      case LPX_FR:
        lpx_set_col_bnds(lp, i, LPX_UP, lo, up);
        break;
      case LPX_UP:
        lpx_set_col_bnds(lp, i, LPX_UP, lo, up);
        break;
      case LPX_LO:
      case LPX_DB:
      case LPX_FX:
        if (lo==up) 
          lpx_set_col_bnds(lp, i, LPX_FX, lo, up);
        else 
          lpx_set_col_bnds(lp, i, LPX_DB, lo, up);
        break;
      default: ;
        //FIXME error
      }
    }
  }
  
//   void LpGlpk::_setRowLowerBound(int i, Value lo)
//   {
//     if (lo==INF) {
//       //FIXME error
//     }
//     int b=lpx_get_row_type(lp, i);
//     double up=lpx_get_row_ub(lp, i); 
//     if (lo==-INF) {
//       switch (b) {
//       case LPX_FR:
//       case LPX_LO:
//      lpx_set_row_bnds(lp, i, LPX_FR, lo, up);
//      break;
//       case LPX_UP:
//      break;
//       case LPX_DB:
//       case LPX_FX:
//      lpx_set_row_bnds(lp, i, LPX_UP, lo, up);
//      break;
//       default: ;
//      //FIXME error
//       }
//     } else {
//       switch (b) {
//       case LPX_FR:
//       case LPX_LO:
//      lpx_set_row_bnds(lp, i, LPX_LO, lo, up);
//      break;
//       case LPX_UP:     
//       case LPX_DB:
//       case LPX_FX:
//      if (lo==up) 
//        lpx_set_row_bnds(lp, i, LPX_FX, lo, up);
//      else 
//        lpx_set_row_bnds(lp, i, LPX_DB, lo, up);
//      break;
//       default: ;
//      //FIXME error
//       }
//     }
//   }
  
//   void LpGlpk::_setRowUpperBound(int i, Value up)
//   {
//     if (up==-INF) {
//       //FIXME error
//     }
//     int b=lpx_get_row_type(lp, i);
//     double lo=lpx_get_row_lb(lp, i);
//     if (up==INF) {
//       switch (b) {
//       case LPX_FR:
//       case LPX_LO:
//      break;
//       case LPX_UP:
//      lpx_set_row_bnds(lp, i, LPX_FR, lo, up);
//      break;
//       case LPX_DB:
//       case LPX_FX:
//      lpx_set_row_bnds(lp, i, LPX_LO, lo, up);
//      break;
//       default: ;
//      //FIXME error
//       }
//     } else {
//       switch (b) {
//       case LPX_FR:
//      lpx_set_row_bnds(lp, i, LPX_UP, lo, up);
//      break;
//       case LPX_UP:
//      lpx_set_row_bnds(lp, i, LPX_UP, lo, up);
//      break;
//       case LPX_LO:
//       case LPX_DB:
//       case LPX_FX:
//      if (lo==up) 
//        lpx_set_row_bnds(lp, i, LPX_FX, lo, up);
//      else 
//        lpx_set_row_bnds(lp, i, LPX_DB, lo, up);
//      break;
//       default: ;
//      //FIXME error
//       }
//     }
//   }

  void LpGlpk::_setRowBounds(int i, Value lb, Value ub)
  {
    //Bad parameter
    if (lb==INF || ub==-INF) {
      //FIXME error
    }

    if (lb == -INF){
      if (ub == INF){
        lpx_set_row_bnds(lp, i, LPX_FR, lb, ub);
      }
      else{
        lpx_set_row_bnds(lp, i, LPX_UP, lb, ub);
      }
    }
    else{
      if (ub==INF){
        lpx_set_row_bnds(lp, i, LPX_LO, lb, ub);

      }
      else{
        if (lb == ub){
          lpx_set_row_bnds(lp, i, LPX_FX, lb, ub);
        }
        else{
          lpx_set_row_bnds(lp, i, LPX_DB, lb, ub);
        }
      }
    }

  }
  
  void LpGlpk::_setObjCoeff(int i, Value obj_coef)
  {
    //i=0 means the constant term (shift)
    lpx_set_obj_coef(lp, i, obj_coef);
  }

  void LpGlpk::_clearObj()
  {
    for (int i=0;i<=lpx_get_num_cols(lp);++i){
      lpx_set_obj_coef(lp, i, 0);
    }
  }
//   void LpGlpk::_setObj(int length,
//                     int  const * indices, 
//                     Value  const * values )
//   {
//     Value new_values[1+lpx_num_cols()];
//     for (i=0;i<=lpx_num_cols();++i){
//       new_values[i]=0;
//     }
//     for (i=1;i<=length;++i){
//       new_values[indices[i]]=values[i];
//     }
    
//     for (i=0;i<=lpx_num_cols();++i){
//       lpx_set_obj_coef(lp, i, new_values[i]);
//     }
//   }

  LpGlpk::SolveExitStatus LpGlpk::_solve()
  {
    int i =  lpx_simplex(lp);
    switch (i) {
    case LPX_E_OK: 
      return SOLVED;
    default:
      return UNSOLVED;
    }
  }

  LpGlpk::Value LpGlpk::_getPrimal(int i)
  {
    return lpx_get_col_prim(lp,i);
  }

  LpGlpk::Value LpGlpk::_getDual(int i)
  {
    return lpx_get_row_dual(lp,i);
  }
  
  LpGlpk::Value LpGlpk::_getPrimalValue()
  {
    return lpx_get_obj_val(lp);
  }
  bool LpGlpk::_isBasicCol(int i) {
    return (lpx_get_col_stat(lp, i)==LPX_BS);
  }
  
 
  LpGlpk::SolutionStatus LpGlpk::_getPrimalStatus()
  {
    int stat=  lpx_get_status(lp);
    switch (stat) {
    case LPX_UNDEF://Undefined (no solve has been run yet)
      return UNDEFINED;
    case LPX_NOFEAS://There is no feasible solution (primal, I guess)
    case LPX_INFEAS://Infeasible 
      return INFEASIBLE;
    case LPX_UNBND://Unbounded
      return INFINITE;
    case LPX_FEAS://Feasible
      return FEASIBLE;
    case LPX_OPT://Feasible
      return OPTIMAL;
    default:
      return UNDEFINED; //to avoid gcc warning
      //FIXME error
    }
  }

  LpGlpk::SolutionStatus LpGlpk::_getDualStatus()
  {
//     std::cout<<"Itt megy: "<<lpx_get_dual_stat(lp)<<std::endl;
//     std::cout<<"Itt a primal: "<<lpx_get_prim_stat(lp)<<std::endl;

    switch (lpx_get_dual_stat(lp)) {
    case LPX_D_UNDEF://Undefined (no solve has been run yet)
      return UNDEFINED;
    case LPX_D_NOFEAS://There is no dual feasible solution 
//    case LPX_D_INFEAS://Infeasible 
      return INFEASIBLE;
    case LPX_D_FEAS://Feasible    
      switch (lpx_get_status(lp)) {
      case LPX_NOFEAS:
        return INFINITE;
      case LPX_OPT:
        return OPTIMAL;
      default:
        return FEASIBLE;
      }
    default:
      return UNDEFINED; //to avoid gcc warning
      //FIXME error
    }
  }

  LpGlpk::ProblemTypes LpGlpk::_getProblemType()
  {
      //int stat=  lpx_get_status(lp);
    int statp=  lpx_get_prim_stat(lp);
    int statd=  lpx_get_dual_stat(lp);
    if (statp==LPX_P_FEAS && statd==LPX_D_FEAS)
        return PRIMAL_DUAL_FEASIBLE;
    if (statp==LPX_P_FEAS && statd==LPX_D_NOFEAS)
        return PRIMAL_FEASIBLE_DUAL_INFEASIBLE;
    if (statp==LPX_P_NOFEAS && statd==LPX_D_FEAS)
        return PRIMAL_INFEASIBLE_DUAL_FEASIBLE;
    if (statp==LPX_P_NOFEAS && statd==LPX_D_NOFEAS)
        return PRIMAL_DUAL_INFEASIBLE;
    //In all other cases
    return UNKNOWN;
  }

  void LpGlpk::_setMax()
  {
    lpx_set_obj_dir(lp, LPX_MAX);
  }

  void LpGlpk::_setMin()
  {
    lpx_set_obj_dir(lp, LPX_MIN);
  }

 
  void LpGlpk::messageLevel(int m)
  {
    lpx_set_int_parm(lp, LPX_K_MSGLEV, m);
  }

  void LpGlpk::presolver(bool b)
  {
    lpx_set_int_parm(lp, LPX_K_PRESOL, b);
  }

 
} //END OF NAMESPACE LEMON