/* -*- C++ -*-

 *

 * This file is a part of LEMON, a generic C++ optimization library

 *

 * Copyright (C) 2003-2007

 * 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>

#if GLP_MAJOR_VERSION > 4 || (GLP_MAJOR_VERSION == 4 && GLP_MINOR_VERSION > 15)

#define LEMON_glp(func) (glp_##func)

#define LEMON_lpx(func) (lpx_##func)

#define LEMON_GLP(def) (GLP_##def)

#define LEMON_LPX(def) (LPX_##def)

#else

#define LEMON_glp(func) (lpx_##func)

#define LEMON_lpx(func) (lpx_##func)

#define LEMON_GLP(def) (LPX_##def)

#define LEMON_LPX(def) (LPX_##def)

#endif

namespace lemon {

  LpGlpk::LpGlpk() : Parent() {

    solved = false;

    rows = _lp_bits::LpId(1);

    cols = _lp_bits::LpId(1);

    lp = LEMON_glp(create_prob)();

    LEMON_glp(create_index)(lp);

    LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_DUAL), 1);

    messageLevel(0);

  }

  LpGlpk::LpGlpk(const LpGlpk &glp) : Parent() {

    solved = false;

    rows = _lp_bits::LpId(1);

    cols = _lp_bits::LpId(1);

    lp = LEMON_glp(create_prob)();

    LEMON_glp(create_index)(lp);

    ///\todo control function for this:

    LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_DUAL), 1);

    messageLevel(0);

    //Coefficient matrix, row bounds

    LEMON_glp(add_rows)(lp, LEMON_glp(get_num_rows)(glp.lp));

    LEMON_glp(add_cols)(lp, LEMON_glp(get_num_cols)(glp.lp));

    int len;

    int ind[1+LEMON_glp(get_num_cols)(glp.lp)];

    Value val[1+LEMON_glp(get_num_cols)(glp.lp)];

    for (int i=1;i<=LEMON_glp(get_num_rows)(glp.lp);++i)

      {

	len=LEMON_glp(get_mat_row)(glp.lp,i,ind,val);

	LEMON_glp(set_mat_row)(lp, i,len,ind,val);

	LEMON_glp(set_row_bnds)(lp,i,

				LEMON_glp(get_row_type)(glp.lp,i),

				LEMON_glp(get_row_lb)(glp.lp,i),

				LEMON_glp(get_row_ub)(glp.lp,i));

      }

    //Objective function, coloumn bounds

    LEMON_glp(set_obj_dir)(lp, LEMON_glp(get_obj_dir)(glp.lp));

    //Objectif function's constant term treated separately

    LEMON_glp(set_obj_coef)(lp,0,LEMON_glp(get_obj_coef)(glp.lp,0));

    for (int i=1;i<=LEMON_glp(get_num_cols)(glp.lp);++i)

      {

	LEMON_glp(set_obj_coef)(lp,i,

				LEMON_glp(get_obj_coef)(glp.lp,i));

	LEMON_glp(set_col_bnds)(lp,i,

				LEMON_glp(get_col_type)(glp.lp,i),

				LEMON_glp(get_col_lb)(glp.lp,i),

				LEMON_glp(get_col_ub)(glp.lp,i));

      }

  }

  LpGlpk::~LpGlpk() {

    LEMON_glp(delete_prob)(lp);

  }

  int LpGlpk::_addCol() { 

    int i=LEMON_glp(add_cols)(lp, 1);

    LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), 0.0, 0.0);

    solved = false;

    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=LEMON_glp(add_rows)(lp, 1);

    solved = false;

    return i;

  }

  void LpGlpk::_eraseCol(int i) {

    int ca[2];

    ca[1]=i;

    LEMON_glp(del_cols)(lp, 1, ca);

    solved = false;

  }

  void LpGlpk::_eraseRow(int i) {

    int ra[2];

    ra[1]=i;

    LEMON_glp(del_rows)(lp, 1, ra);

    solved = false;

  }

  void LpGlpk::_getColName(int c, std::string & name) const

  {

    const char *n = LEMON_glp(get_col_name)(lp,c);

    name = n?n:"";

  }

  void LpGlpk::_setColName(int c, const std::string & name)

  {

    LEMON_glp(set_col_name)(lp,c,const_cast<char*>(name.c_str()));

  }

  int LpGlpk::_colByName(const std::string& name) const

  {

    int k = LEMON_glp(find_col)(lp, const_cast<char*>(name.c_str()));

    return k > 0 ? k : -1; 

  }

  void LpGlpk::_setRowCoeffs(int i, ConstRowIterator b, ConstRowIterator e) 

  {

    std::vector<int> indices;

    std::vector<Value> values;

    indices.push_back(0);

    values.push_back(0);

    for(ConstRowIterator it=b; it!=e; ++it) {

      indices.push_back(it->first);

      values.push_back(it->second);

    }

    LEMON_glp(set_mat_row)(lp, i, values.size() - 1, 

				&indices[0], &values[0]);

    solved = false;

  }

  void LpGlpk::_getRowCoeffs(int ix, RowIterator b) const

  {

    int length = LEMON_glp(get_mat_row)(lp, ix, 0, 0);

    std::vector<int> indices(length + 1);

    std::vector<Value> values(length + 1);

    LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);

    for (int i = 1; i <= length; ++i) {

      *b = std::make_pair(indices[i], values[i]);

      ++b;

    }

  }

  void LpGlpk::_setColCoeffs(int ix, ConstColIterator b, ConstColIterator e) {

    std::vector<int> indices;

    std::vector<Value> values;

    indices.push_back(0);

    values.push_back(0);

    for(ConstColIterator it=b; it!=e; ++it) {

      indices.push_back(it->first);

      values.push_back(it->second);

    }

    LEMON_glp(set_mat_col)(lp, ix, values.size() - 1, 

				&indices[0], &values[0]);

    solved = false;

  }

  void LpGlpk::_getColCoeffs(int ix, ColIterator b) const

  {

    int length = LEMON_glp(get_mat_col)(lp, ix, 0, 0);

    std::vector<int> indices(length + 1);

    std::vector<Value> values(length + 1);

    LEMON_glp(get_mat_col)(lp, ix, &indices[0], &values[0]);

    for (int i = 1; i <= length; ++i) {

      *b = std::make_pair(indices[i], values[i]);

      ++b;

    }

  }

  void LpGlpk::_setCoeff(int ix, int jx, Value value) 

  {

    if (LEMON_glp(get_num_cols)(lp) < LEMON_glp(get_num_rows)(lp)) {

      int length=LEMON_glp(get_mat_row)(lp, ix, 0, 0);

      std::vector<int> indices(length + 2);

      std::vector<Value> values(length + 2);

      LEMON_glp(get_mat_row)(lp, ix, &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]==jx){

          found=true;

          values[i]=value;

          break;

        }

      }

      if (!found){

        ++length;

        indices[length]=jx;

        values[length]=value;

      }

      LEMON_glp(set_mat_row)(lp, ix, length, &indices[0], &values[0]);

    } else {

      int length=LEMON_glp(get_mat_col)(lp, jx, 0, 0);

      std::vector<int> indices(length + 2);

      std::vector<Value> values(length + 2);

      LEMON_glp(get_mat_col)(lp, jx, &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]==jx){

          found=true;

          values[i]=value;

          break;

        }

      }

      if (!found){

        ++length;

        indices[length]=ix;

        values[length]=value;

      }

      LEMON_glp(set_mat_col)(lp, jx, length, &indices[0], &values[0]);

    }

    solved = false;

  }

  LpGlpk::Value LpGlpk::_getCoeff(int ix, int jx) const

  {

    int length=LEMON_glp(get_mat_row)(lp, ix, 0, 0);

    std::vector<int> indices(length + 1);

    std::vector<Value> values(length + 1);

    LEMON_glp(get_mat_row)(lp, ix, &indices[0], &values[0]);

    //The following code does not suppose that the elements of the

    //array indices are sorted

    for (int i = 1; i <= length; ++i) {

      if (indices[i]==jx){

	return values[i];

      }

    }

    return 0;

  }

  void LpGlpk::_setColLowerBound(int i, Value lo)

  {

    if (lo==INF) {

      //FIXME error

    }

    int b=LEMON_glp(get_col_type)(lp, i);

    double up=LEMON_glp(get_col_ub)(lp, i);	

    if (lo==-INF) {

      switch (b) {

      case LEMON_GLP(FR):

      case LEMON_GLP(LO):

	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), lo, up);

	break;

      case LEMON_GLP(UP):

	break;

      case LEMON_GLP(DB):

      case LEMON_GLP(FX):

	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);

	break;

      default: ;

	//FIXME error

      }

    } else {

      switch (b) {

      case LEMON_GLP(FR):

      case LEMON_GLP(LO):

	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(LO), lo, up);

	break;

      case LEMON_GLP(UP):	  

      case LEMON_GLP(DB):

      case LEMON_GLP(FX):

	if (lo==up) 

	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FX), lo, up);

	else 

	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(DB), lo, up);

	break;

      default: ;

	//FIXME error

      }

    }

    solved = false;

  }

  LpGlpk::Value LpGlpk::_getColLowerBound(int i) const

  {

    int b=LEMON_glp(get_col_type)(lp, i);

      switch (b) {

      case LEMON_GLP(LO):

      case LEMON_GLP(DB):

      case LEMON_GLP(FX):

	return LEMON_glp(get_col_lb)(lp, i);	

      default: ;

	return -INF;

      }

  }

  void LpGlpk::_setColUpperBound(int i, Value up)

  {

    if (up==-INF) {

      //FIXME error

    }

    int b=LEMON_glp(get_col_type)(lp, i);

    double lo=LEMON_glp(get_col_lb)(lp, i);

    if (up==INF) {

      switch (b) {

      case LEMON_GLP(FR):

      case LEMON_GLP(LO):

	break;

      case LEMON_GLP(UP):

	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FR), lo, up);

	break;

      case LEMON_GLP(DB):

      case LEMON_GLP(FX):

	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(LO), lo, up);

	break;

      default: ;

	//FIXME error

      }

    } else {

      switch (b) {

      case LEMON_GLP(FR):

	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);

	break;

      case LEMON_GLP(UP):

	LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(UP), lo, up);

	break;

      case LEMON_GLP(LO):

      case LEMON_GLP(DB):

      case LEMON_GLP(FX):

	if (lo==up) 

	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(FX), lo, up);

	else 

	  LEMON_glp(set_col_bnds)(lp, i, LEMON_GLP(DB), lo, up);

	break;

      default: ;

	//FIXME error

      }

    }

    solved = false;

  }

  LpGlpk::Value LpGlpk::_getColUpperBound(int i) const

  {

    int b=LEMON_glp(get_col_type)(lp, i);

      switch (b) {

      case LEMON_GLP(UP):

      case LEMON_GLP(DB):

      case LEMON_GLP(FX):

	return LEMON_glp(get_col_ub)(lp, i);	

      default: ;

	return INF;

      }

  }

  void LpGlpk::_setRowBounds(int i, Value lb, Value ub)

  {

    //Bad parameter

    if (lb==INF || ub==-INF) {

      //FIXME error

    }

    if (lb == -INF){

      if (ub == INF){

	LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(FR), lb, ub);

      }

      else{

	LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(UP), lb, ub);

      }

    }

    else{

      if (ub==INF){

	LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(LO), lb, ub);

      }

      else{

	if (lb == ub){

	  LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(FX), lb, ub);

	}

	else{

	  LEMON_glp(set_row_bnds)(lp, i, LEMON_GLP(DB), lb, ub);

	}

      }

    }

    solved = false;

  }

  void LpGlpk::_getRowBounds(int i, Value &lb, Value &ub) const

  {

    int b=LEMON_glp(get_row_type)(lp, i);

    switch (b) {

    case LEMON_GLP(FR):

    case LEMON_GLP(UP):

      lb = -INF;

	break;

    default: 

      lb=LEMON_glp(get_row_lb)(lp, i);

    }

    switch (b) {

    case LEMON_GLP(FR):

    case LEMON_GLP(LO):

      ub = INF;

	break;

    default: 

      ub=LEMON_glp(get_row_ub)(lp, i);

    }

  }

  void LpGlpk::_setObjCoeff(int i, Value obj_coef)

  {

    //i=0 means the constant term (shift)

    LEMON_glp(set_obj_coef)(lp, i, obj_coef);

    solved = false;

  }

  LpGlpk::Value LpGlpk::_getObjCoeff(int i) const {

    //i=0 means the constant term (shift)

    return LEMON_glp(get_obj_coef)(lp, i);

  }

  void LpGlpk::_clearObj()

  {

    for (int i=0;i<=LEMON_glp(get_num_cols)(lp);++i){

      LEMON_glp(set_obj_coef)(lp, i, 0);

    }

    solved = false;

  }

  LpGlpk::SolveExitStatus LpGlpk::_solve()

  {

    // A way to check the problem to be solved

    //LEMON_glp(write_cpxlp(lp,"naittvan.cpx");    

    LEMON_lpx(std_basis)(lp);

    int i =  LEMON_lpx(simplex)(lp);

    switch (i) {

    case LEMON_LPX(E_OK): 

      solved = true;

      return SOLVED;

    default:

      return UNSOLVED;

    }

  }

  LpGlpk::Value LpGlpk::_getPrimal(int i) const

  {

    return LEMON_glp(get_col_prim)(lp,i);

  }

  LpGlpk::Value LpGlpk::_getDual(int i) const

  {

    return LEMON_glp(get_row_dual)(lp,i);

  }

  LpGlpk::Value LpGlpk::_getPrimalValue() const

  {

    return LEMON_glp(get_obj_val)(lp);

  }

  bool LpGlpk::_isBasicCol(int i) const

  {

    return (LEMON_glp(get_col_stat)(lp, i)==LEMON_GLP(BS));

  }

  LpGlpk::SolutionStatus LpGlpk::_getPrimalStatus() const

  {

    if (!solved) return UNDEFINED;

    int stat=  LEMON_lpx(get_status)(lp);

    switch (stat) {

    case LEMON_LPX(UNDEF)://Undefined (no solve has been run yet)

      return UNDEFINED;

    case LEMON_LPX(NOFEAS)://There is no feasible solution (primal, I guess)

    case LEMON_LPX(INFEAS)://Infeasible 

      return INFEASIBLE;

    case LEMON_LPX(UNBND)://Unbounded

      return INFINITE;

    case LEMON_LPX(FEAS)://Feasible

      return FEASIBLE;

    case LEMON_LPX(OPT)://Feasible

      return OPTIMAL;

    default:

      return UNDEFINED; //to avoid gcc warning

      //FIXME error

    }

  }

  LpGlpk::SolutionStatus LpGlpk::_getDualStatus() const

  {

    if (!solved) return UNDEFINED;

    switch (LEMON_lpx(get_dual_stat)(lp)) {

    case LEMON_LPX(D_UNDEF)://Undefined (no solve has been run yet)

      return UNDEFINED;

    case LEMON_LPX(D_NOFEAS)://There is no dual feasible solution 

//    case LEMON_LPX(D_INFEAS://Infeasible 

      return INFEASIBLE;

    case LEMON_LPX(D_FEAS)://Feasible    

      switch (LEMON_lpx(get_status)(lp)) {

      case LEMON_LPX(NOFEAS):

	return INFINITE;

      case LEMON_LPX(OPT):

	return OPTIMAL;

      default:

	return FEASIBLE;

      }

    default:

      return UNDEFINED; //to avoid gcc warning

      //FIXME error

    }

  }

  LpGlpk::ProblemTypes LpGlpk::_getProblemType() const

  {

    if (!solved) return UNKNOWN;

      //int stat=  LEMON_glp(get_status(lp);

    int statp=  LEMON_lpx(get_prim_stat)(lp);

    int statd=  LEMON_lpx(get_dual_stat)(lp);

    if (statp==LEMON_LPX(P_FEAS) && statd==LEMON_LPX(D_FEAS))

	return PRIMAL_DUAL_FEASIBLE;

    if (statp==LEMON_LPX(P_FEAS) && statd==LEMON_LPX(D_NOFEAS))

	return PRIMAL_FEASIBLE_DUAL_INFEASIBLE;

    if (statp==LEMON_LPX(P_NOFEAS) && statd==LEMON_LPX(D_FEAS))

	return PRIMAL_INFEASIBLE_DUAL_FEASIBLE;

    if (statp==LEMON_LPX(P_NOFEAS) && statd==LEMON_LPX(D_NOFEAS))

	return PRIMAL_DUAL_INFEASIBLE;

    //In all other cases

    return UNKNOWN;

  }

  void LpGlpk::_setMax()

  {

    solved = false;

    LEMON_glp(set_obj_dir)(lp, LEMON_GLP(MAX));

  }

  void LpGlpk::_setMin()

  {

    solved = false;

    LEMON_glp(set_obj_dir)(lp, LEMON_GLP(MIN));

  }

  bool LpGlpk::_isMax() const

  {

    return (LEMON_glp(get_obj_dir)(lp)==LEMON_GLP(MAX));

  }

  void LpGlpk::messageLevel(int m)

  {

    LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_MSGLEV), m);

  }

  void LpGlpk::presolver(bool b)

  {

    LEMON_lpx(set_int_parm)(lp, LEMON_LPX(K_PRESOL), b);

  }

} //END OF NAMESPACE LEMON