diff --git a/lemon/lp_glpk.cc b/lemon/lp_glpk.cc
new file mode 100644
--- /dev/null
+++ b/lemon/lp_glpk.cc
@@ -0,0 +1,644 @@
+/* -*- mode: C++; indent-tabs-mode: nil; -*-
+ *
+ * This file is a part of LEMON, a generic C++ optimization library.
+ *
+ * Copyright (C) 2003-2008
+ * 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>
+
+extern "C" {
+#include <glpk.h>
+}
+
+#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);
+    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);
+    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;
+    std::vector<int> ind(1+LEMON_glp(get_num_cols)(glp.lp));
+    std::vector<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.begin(),&*val.begin());
+        LEMON_glp(set_mat_row)(lp, i,len,&*ind.begin(),&*val.begin());
+        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));
+      }
+    rows = glp.rows;
+    cols = glp.cols;
+  }
+
+  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]==ix){
+          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