Changeset 2363:2aabce558574 in lemon-0.x for lemon

Timestamp:

02/15/07 15:22:08 (17 years ago)

Author:

Balazs Dezso

Branch:

default

Phase:

public

Convert:

svn:c9d7d8f5-90d6-0310-b91f-818b3a526b0e/lemon/trunk@3174

Message:

Changes on the LP interface

_FixId => LpId?

• handling of not common ids soplex

LpGlpk? row and col erase bug fix

• calling lpx_std_basis before simplex

LpSoplex?

• added getter functions
• better m4 file
• integration to the tests
• better handling of unsolved lps

Location:

lemon

Files:

• Makefile.am (modified) (2 diffs)
• bits/lp_id.h (added)
• lp_base.h (modified) (11 diffs)
• lp_cplex.cc (modified) (1 diff)
• lp_glpk.cc (modified) (4 diffs)
• lp_skeleton.h (modified) (1 diff)
• lp_soplex.cc (modified) (13 diffs)
• lp_soplex.h (modified) (3 diffs)
• lp_utils.h (modified) (3 diffs)

lemon/Makefile.am

@@ -16 +16 @@
         lemon/random.cc
 
-lemon_libemon_la_CXXFLAGS = $(GLPK_CFLAGS) $(CPLEX_CFLAGS) $(SOPLEX_CFLAGS)
-lemon_libemon_la_LDFLAGS = $(GLPK_LIBS) $(CPLEX_LIBS) $(SOPLEX_CFLAGS)
+lemon_libemon_la_CXXFLAGS = $(GLPK_CFLAGS) $(CPLEX_CFLAGS) $(SOPLEX_CXXFLAGS)
+lemon_libemon_la_LDFLAGS = $(GLPK_LIBS) $(CPLEX_LIBS) $(SOPLEX_LIBS)
 
 if HAVE_GLPK
@@ -124 +124 @@
         lemon/bits/item_reader.h \
         lemon/bits/item_writer.h \
+        lemon/bits/lp_id.h \
         lemon/bits/map_extender.h \
         lemon/bits/mingw32_time.h \

lemon/lp_base.h

@@ -28 +28 @@
 #include<cmath>
 
-#include<lemon/bits/utility.h>
 #include<lemon/error.h>
 #include<lemon/bits/invalid.h>
+#include<lemon/bits/utility.h>
+#include<lemon/bits/lp_id.h>
 
 ///\file
@@ -37 +38 @@
 namespace lemon {
 
-
-  ///Internal data structure to convert floating id's to fix one's
-    
-  ///\todo This might be implemented to be also usable in other places.
-  class _FixId 
-  {
-  protected:
-    int _first_index;
-    int first_free;
-  public:
-    std::vector<int> index;
-    std::vector<int> cross;
-    _FixId() : _first_index(-1), first_free(-1) {};
-    ///Convert a floating id to a fix one
-
-    ///\param n is a floating id
-    ///\return the corresponding fix id
-    int fixId(int n) const {return cross[n];}
-    ///Convert a fix id to a floating one
-
-    ///\param n is a fix id
-    ///\return the corresponding floating id
-    int floatingId(int n) const { return index[n];}
-    ///Add a new floating id.
-
-    ///\param n is a floating id
-    ///\return the fix id of the new value
-    ///\todo Multiple additions should also be handled.
-    int insert(int n)
-    {
-      if(cross.empty()) _first_index=n;
-      if(n>=int(cross.size())) {
-        cross.resize(n+1);
-        if(first_free==-1) {
-          cross[n]=index.size();
-          index.push_back(n);
-        }
-        else {
-          cross[n]=first_free;
-          int next=index[first_free];
-          index[first_free]=n;
-          first_free=next;
-        }
-        return cross[n];
-      }
-      else {
-        ///\todo Create an own exception type.
-        throw LogicError(); //floatingId-s must form a continuous range;
-      }
-    }
-    ///Remove a fix id.
-
-    ///\param n is a fix id
-    ///
-    void erase(int n) 
-    {
-      int fl=index[n];
-      index[n]=first_free;
-      first_free=n;
-      for(int i=fl+1;i<int(cross.size());++i) {
-        cross[i-1]=cross[i];
-        index[cross[i]]--;
-      }
-      cross.pop_back();
-    }
-    ///An upper bound on the largest fix id.
-
-    ///\todo Do we need this?
-    ///
-    std::size_t maxFixId() { return cross.size()-1; }
-  
-    ///Returns the first (smallest) inserted index
-
-    ///Returns the first (smallest) inserted index
-    ///or -1 if no index has been inserted before.
-    int firstIndex() {return _first_index;}
-  };
-
   ///Common base class for LP solvers
   
@@ -122 +45 @@
 
   protected:
-    _FixId rows;
-    _FixId cols;
-
+
+    _lp_bits::LpId rows;
+    _lp_bits::LpId cols;
+    
   public:
 
@@ -216 +140 @@
       ColIt(LpSolverBase &lp) : _lp(&lp)
       {
-        id = _lp->cols.cross.empty()?-1:
-          _lp->cols.fixId(_lp->cols.firstIndex());
+        _lp->cols.firstFix(id);
       }
       ColIt(const Invalid&) : Col(INVALID) {}
       ColIt &operator++() 
       {
-        int fid = _lp->cols.floatingId(id)+1;
-        id = unsigned(fid)<_lp->cols.cross.size() ? _lp->cols.fixId(fid) : -1;
+        _lp->cols.nextFix(id);
         return *this;
       }
@@ -747 +669 @@
     virtual void _setColUpperBound(int i, Value value) = 0;
     virtual Value _getColUpperBound(int i) = 0;
-//     virtual void _setRowLowerBound(int i, Value value) = 0;
-//     virtual void _setRowUpperBound(int i, Value value) = 0;
     virtual void _setRowBounds(int i, Value lower, Value upper) = 0;
     virtual void _getRowBounds(int i, Value &lower, Value &upper)=0;
@@ -796 +716 @@
 
     ///Add a new empty column (i.e a new variable) to the LP
-    Col addCol() { Col c; c.id=cols.insert(_addCol()); return c;}
+    Col addCol() { Col c; _addCol(); c.id = cols.addId(); return c;}
 
     ///\brief Adds several new columns
@@ -888 +808 @@
     ///This function adds a new empty row (i.e a new constraint) to the LP.
     ///\return The created row
-    Row addRow() { Row r; r.id=rows.insert(_addRow()); return r;}
+    Row addRow() { Row r; _addRow(); r.id = rows.addId(); return r;}
 
     ///\brief Add several new rows
@@ -966 +886 @@
       _setRowCoeffs(_lpId(r), LpRowIterator(e.begin(), *this),
                     LpRowIterator(e.end(), *this));
-//       _setRowLowerBound(_lpId(r),l-e.constComp());
-//       _setRowUpperBound(_lpId(r),u-e.constComp());
        _setRowBounds(_lpId(r),l-e.constComp(),u-e.constComp());
     }
@@ -1009 +927 @@
     void eraseCol(Col c) {
       _eraseCol(_lpId(c));
-      cols.erase(c.id);
+      cols.eraseId(c.id);
     }
     ///Erase a  row (i.e a constraint) from the LP
@@ -1017 +935 @@
     void eraseRow(Row r) {
       _eraseRow(_lpId(r));
-      rows.erase(r.id);
+      rows.eraseId(r.id);
     }
 
@@ -1217 +1135 @@
 #endif
     
-//     /// Set the lower bound of a row (i.e a constraint)
-
-//     /// The lower bound of a linear expression (row) has to be given by an 
-//     /// extended number of type Value, i.e. a finite number of type 
-//     /// Value or -\ref INF.
-//     void rowLowerBound(Row r, Value value) {
-//       _setRowLowerBound(_lpId(r),value);
-//     };
-//     /// Set the upper bound of a row (i.e a constraint)
-
-//     /// The upper bound of a linear expression (row) has to be given by an 
-//     /// extended number of type Value, i.e. a finite number of type 
-//     /// Value or \ref INF.
-//     void rowUpperBound(Row r, Value value) {
-//       _setRowUpperBound(_lpId(r),value);
-//     };
 
     /// Set the lower and the upper bounds of a row (i.e a constraint)
 
-    /// The lower and the upper bound of
-    /// a constraint (row) have to be given by an 
-    /// extended number of type Value, i.e. a finite number of type 
-    /// Value, -\ref INF or \ref INF. There is no separate function for the 
-    /// lower and the upper bound because that would have been hard to implement 
-    /// for CPLEX.
+    /// The lower and the upper bound of a constraint (row) have to be
+    /// given by an extended number of type Value, i.e. a finite
+    /// number of type Value, -\ref INF or \ref INF. There is no
+    /// separate function for the lower and the upper bound because
+    /// that would have been hard to implement for CPLEX.
     void rowBounds(Row c, Value lower, Value upper) {
       _setRowBounds(_lpId(c),lower, upper);

lemon/lp_cplex.cc

@@ -25 +25 @@
   
   LpCplex::LpCplex() : LpSolverBase() {
-
     //    env = CPXopenCPLEXdevelop(&status);     
     env = CPXopenCPLEX(&status);     

lemon/lp_glpk.cc

@@ -25 +25 @@
 
 
-  LpGlpk::LpGlpk() : Parent(), 
-                     lp(lpx_create_prob()) {
-    ///\todo constrol function for this:
+  LpGlpk::LpGlpk() : Parent() {
+    rows = _lp_bits::LpId(1);
+    cols = _lp_bits::LpId(1);
+    lp = lpx_create_prob();
+    ///\todo control 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:
+  LpGlpk::LpGlpk(const LpGlpk &glp) : Parent() {
+    rows = _lp_bits::LpId(1);
+    cols = _lp_bits::LpId(1);
+    lp = lpx_create_prob();
+    ///\todo control function for this:
     lpx_set_int_parm(lp, LPX_K_DUAL, 1);
     messageLevel(0);
@@ -350 +354 @@
   }
   
-//   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)
   {
@@ -509 +427 @@
     }
   }
-//   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()
@@ -531 +433 @@
     //lpx_write_cpxlp(lp,"naittvan.cpx");    
 
+    lpx_std_basis(lp);
     int i =  lpx_simplex(lp);
+    
     switch (i) {
     case LPX_E_OK: 

lemon/lp_skeleton.h

@@ -31 +31 @@
     
   protected:
+
     ///\e
     virtual LpSolverBase &_newLp();

lemon/lp_soplex.cc

@@ -28 +28 @@
   
   LpSoplex::LpSoplex() : LpSolverBase() {
+    rows.setIdHandler(relocateIdHandler);
+    cols.setIdHandler(relocateIdHandler);
     soplex = new soplex::SoPlex;
+    solved = false;
   }
   
@@ -48 +51 @@
   int LpSoplex::_addCol() {
     soplex::LPCol col;
+    col.setLower(-soplex::infinity);
+    col.setUpper(soplex::infinity);
     soplex->addCol(col);
 
     colNames.push_back(std::string());
-    primal.push_back(0.0);
+    primal_value.push_back(0.0);
+    solved = false;
 
     return soplex->nCols() - 1;
@@ -58 +64 @@
   int LpSoplex::_addRow() {
     soplex::LPRow row;
+    row.setLhs(-soplex::infinity);
+    row.setRhs(soplex::infinity);
     soplex->addRow(row);
 
-    dual.push_back(0.0);
+    dual_value.push_back(0.0);
+    solved = false;
 
     return soplex->nRows() - 1;
@@ -68 +77 @@
   void LpSoplex::_eraseCol(int i) {
     soplex->removeCol(i);
-    primal[i] = primal.back();
-    primal.pop_back();
+    colNames[i] = colNames.back();
+    colNames.pop_back();
+    primal_value[i] = primal_value.back();
+    primal_value.pop_back();
+    solved = false;
   }
   
   void LpSoplex::_eraseRow(int i) {
     soplex->removeRow(i);
-    dual[i] = dual.back();
-    dual.pop_back();
+    dual_value[i] = dual_value.back();
+    dual_value.pop_back();
+    solved = false;
   }
   
@@ -94 +107 @@
       soplex->changeElement(i, it->first, it->second);
     }
+    solved = false;
   }
   
@@ -103 +117 @@
       soplex->changeElement(it->first, j, it->second);
     }
-  }
-  
-  void LpSoplex::_setCoeff(int row, int col, Value value) {
-    soplex->changeElement(row, col, value);
+    solved = false;
+  }
+  
+  void LpSoplex::_setCoeff(int i, int j, Value value) {
+    soplex->changeElement(i, j, value);
+    solved = false;
+  }
+
+  LpSoplex::Value LpSoplex::_getCoeff(int i, int j) {
+    return soplex->rowVector(i)[j];
   }
 
   void LpSoplex::_setColLowerBound(int i, Value value) {
-    soplex->changeLower(i, value); 
-  }
-  
+    soplex->changeLower(i, value != -INF ? value : -soplex::infinity); 
+    solved = false;
+  }
+  
+  LpSoplex::Value LpSoplex::_getColLowerBound(int i) {
+    double value = soplex->lower(i);
+    return value != -soplex::infinity ? value : -INF;
+  }
+
   void LpSoplex::_setColUpperBound(int i, Value value) {
-    soplex->changeUpper(i, value); 
+    soplex->changeUpper(i, value != INF ? value : soplex::infinity); 
+    solved = false;
+  }
+
+  LpSoplex::Value LpSoplex::_getColUpperBound(int i) {
+    double value = soplex->upper(i);
+    return value != soplex::infinity ? value : INF;
   }
 
   void LpSoplex::_setRowBounds(int i, Value lb, Value ub) {
-    soplex->changeRange(i, lb, ub);
+    soplex->changeRange(i, lb != -INF ? lb : -soplex::infinity,
+                        ub != INF ? ub : soplex::infinity);
+    solved = false;
+  }
+  void LpSoplex::_getRowBounds(int i, Value &lower, Value &upper) {
+    lower = soplex->lhs(i);
+    if (lower == -soplex::infinity) lower = -INF;
+    upper = soplex->rhs(i);
+    if (upper == -soplex::infinity) upper = INF;
   }
 
   void LpSoplex::_setObjCoeff(int i, Value obj_coef) {
     soplex->changeObj(i, obj_coef);
+    solved = false;
+  }
+
+  LpSoplex::Value LpSoplex::_getObjCoeff(int i) {
+    return soplex->obj(i);
   }
 
@@ -129 +174 @@
       soplex->changeObj(i, 0.0);
     }
+    solved = false;
   }
 
@@ -134 +180 @@
     soplex::SPxSolver::Status status = soplex->solve();
 
-    soplex::Vector pv(primal.size(), &primal[0]);
+    soplex::Vector pv(primal_value.size(), &primal_value[0]);
     soplex->getPrimal(pv);
 
-    soplex::Vector dv(dual.size(), &dual[0]);
+    soplex::Vector dv(dual_value.size(), &dual_value[0]);
     soplex->getDual(dv);
 
@@ -144 +190 @@
     case soplex::SPxSolver::INFEASIBLE:
     case soplex::SPxSolver::UNBOUNDED:
+      solved = true;
       return SOLVED;
     default:
@@ -151 +198 @@
 
   LpSoplex::Value LpSoplex::_getPrimal(int i) {
-    return primal[i];
+    return primal_value[i];
   }
 
   LpSoplex::Value LpSoplex::_getDual(int i) {
-    return dual[i];
+    return dual_value[i];
   }
   
@@ -167 +214 @@
 
   LpSoplex::SolutionStatus LpSoplex::_getPrimalStatus() {
+    if (!solved) return UNDEFINED;
     switch (soplex->status()) {
     case soplex::SPxSolver::OPTIMAL:
@@ -180 +228 @@
 
   LpSoplex::SolutionStatus LpSoplex::_getDualStatus() {
-    switch (0) {
-    case 0:
+    if (!solved) return UNDEFINED;
+    switch (soplex->status()) {
+    case soplex::SPxSolver::OPTIMAL:
+      return OPTIMAL;
+    case soplex::SPxSolver::UNBOUNDED:
+      return INFEASIBLE;
+    default:
       return UNDEFINED;
-      return OPTIMAL;
-      return INFEASIBLE;
-      return UNDEFINED;
     }
   }
 
   LpSoplex::ProblemTypes LpSoplex::_getProblemType() {
-    switch (0) {
-    case 0:
+    if (!solved) return UNKNOWN;
+    switch (soplex->status()) {
+    case soplex::SPxSolver::OPTIMAL:
       return PRIMAL_DUAL_FEASIBLE;
+    case soplex::SPxSolver::UNBOUNDED:
       return PRIMAL_FEASIBLE_DUAL_INFEASIBLE;
+    default:
       return UNKNOWN;
     }
@@ -200 +253 @@
   void LpSoplex::_setMax() {
     soplex->changeSense(soplex::SPxSolver::MAXIMIZE);
+    solved = false;
   }
   void LpSoplex::_setMin() {
     soplex->changeSense(soplex::SPxSolver::MINIMIZE);
-  }
+    solved = false;
+  }
+  bool LpSoplex::_isMax() {
+    return soplex->spxSense() == soplex::SPxSolver::MAXIMIZE;
+  }
+
   
 } //namespace lemon

lemon/lp_soplex.h

@@ -41 +41 @@
   protected:
 
+    _lp_bits::RelocateIdHandler relocateIdHandler;
+
     soplex::SoPlex* soplex;
+    bool solved;
 
     std::vector<std::string> colNames;
 
-    std::vector<Value> primal;
-    std::vector<Value> dual;
+    std::vector<Value> primal_value;
+    std::vector<Value> dual_value;
 
 
@@ -68 +71 @@
     virtual void _eraseCol(int i);
     virtual void _eraseRow(int i);
-    virtual void _getColName(int col,       std::string & name);
+    virtual void _getColName(int col, std::string & name);
     virtual void _setColName(int col, const std::string & name);
     virtual void _setRowCoeffs(int i, LpRowIterator b, LpRowIterator e);
     virtual void _setColCoeffs(int i, LpColIterator b, LpColIterator e);
     virtual void _setCoeff(int row, int col, Value value);
+    virtual Value _getCoeff(int row, int col);
     virtual void _setColLowerBound(int i, Value value);
+    virtual Value _getColLowerBound(int i);
     virtual void _setColUpperBound(int i, Value value);
+    virtual Value _getColUpperBound(int i);
     virtual void _setRowBounds(int i, Value lower, Value upper);
+    virtual void _getRowBounds(int i, Value &lower, Value &upper);
     virtual void _setObjCoeff(int i, Value obj_coef);
+    virtual Value _getObjCoeff(int i);
     virtual void _clearObj();
     
@@ -92 +100 @@
     virtual void _setMax();
     virtual void _setMin();
+    virtual bool _isMax();
 
   };

lemon/lp_utils.h

@@ -23 +23 @@
 
 #include <lemon/lemon_reader.h>
+#include <lemon/lemon_writer.h>
 
 namespace lemon {
@@ -320 +321 @@
     virtual void read(std::istream& is) {
       std::string line;
-      std::map<std::string, LpSolverBase::Col> vars;
       while (getline(is, line)) {
         std::istringstream ls(line);
@@ -360 +360 @@
   };
 
+
+//   class LpWriter : public LemonWriter::SectionWriter {
+//     typedef LemonWriter::SectionWriter Parent;
+//   public:
+
+
+//     /// \brief Constructor.
+//     ///
+//     /// Constructor for LpWriter. It creates the LpWriter and attach
+//     /// it into the given LemonWriter. The lp writer will add
+//     /// variables, constraints and objective function to the
+//     /// given lp solver.
+//     LpWriter(LemonWriter& _writer, LpSolverBase& _lp, 
+//              const std::string& _name = std::string())
+//       : Parent(_writer), lp(_lp), name(_name) {} 
+
+
+//     /// \brief Destructor.
+//     ///
+//     /// Destructor for NodeSetWriter.
+//     virtual ~LpWriter() {}
+
+//   private:
+//     LpWriter(const LpWriter&);
+//     void operator=(const LpWriter&);
+  
+//   protected:
+
+//     /// \brief Gives back true when the SectionWriter can process 
+//     /// the section with the given header line.
+//     ///
+//     /// It gives back the header line of the \c \@lp section.
+//     virtual std::string header() {
+//       std::ostringstream ls(line);
+//       ls << "@lp " << name;
+//       return ls.str();
+//     }
+
+//   private:
+
+//     std::ostream& writeConstraint(std::ostream& is, LpSolverBase::Constr& c) {
+//       char x;
+
+      
+//       LpSolverBase::Expr e1, e2;
+//       Relation op1;
+//       is >> std::ws;
+//       writexpression(is, e1);
+//       is >> std::ws;
+//       writeRelation(is, op1);
+//       is >> std::ws;
+//       writexpression(is, e2);
+//       if (!is.get(x)) {
+//         if (op1 == LE) {
+//           c = e1 <= e2;
+//         } else if (op1 == GE) {
+//           c = e1 >= e2;
+//         } else {
+//           c = e1 == e2;
+//         }
+//       } else {
+//         is.putback(x);
+//         LpSolverBase::Expr e3;
+//         Relation op2;
+//         writeRelation(is, op2);
+//         is >> std::ws;
+//         writexpression(is, e3);
+//         if (!e1.empty() || !e3.empty()) {
+//           throw DataFormatError("Wrong range format");
+//         }
+//         if (op2 != op1 || op1 == EQ) {
+//           throw DataFormatError("Wrong range format");
+//         }
+//         if (op1 == LE) {
+//           c = e1.constComp() <= e2 <= e3.constComp();
+//         } else {
+//           c = e1.constComp() >= e2 >= e3.constComp();
+//         }
+//       }
+//     }
+
+//     std::ostream& writexpression(std::ostream& is, LpSolverBase::Expr& e) {
+//       LpSolverBase::Col c;
+//       double d;
+//       char x;
+//       writelement(is, c, d);
+//       if (c != INVALID) {
+//         e += d * c;
+//       } else {
+//         e += d;
+//       }
+//       is >> std::ws;
+//       while (is.get(x) && (x == '+' || x == '-')) {
+//         is >> std::ws;
+//         writelement(is, c, d);
+//         if (c != INVALID) {
+//           e += (x == '+' ? d : -d) * c;
+//         } else {
+//           e += (x == '+' ? d : -d);
+//         }
+//         is >> std::ws;
+//       }
+//       if (!is) {
+//         is.clear();
+//       } else {
+//         is.putback(x);
+//       }
+//       return is;
+//     }
+
+//     std::ostream& writelement(std::ostream& is, 
+//                               LpSolverBase::Col& c, double& d) { 
+//       d = 1.0;
+//       c = INVALID;
+//       char x, y;
+//       if (!is.get(x)) throw DataFormatError("Cannot find lp element");
+//       if (x == '+' || x == '-') {
+//         is >> std::ws;
+//         d *= x == '-' ? -1 : 1;
+//         while (is.get(x) && (x == '+' || x == '-')) {
+//           d *= x == '-' ? -1 : 1;
+//           is >> std::ws;
+//         }
+//         if (!is) throw DataFormatError("Cannot find lp element");
+//       }
+//       if (numFirstChar(x)) {
+//         is.putback(x);
+//         double e;
+//         writeNum(is, e);
+//         d *= e;
+//       } else if (varFirstChar(x)) {
+//         is.putback(x);
+//         LpSolverBase::Col f;
+//         writeCol(is, f);
+//         c = f;
+//       } else {
+//         throw DataFormatError("Invalid expression format");          
+//       }
+//       is >> std::ws;
+//       while (is.get(y) && (y == '*' || y == '/')) {
+//         is >> std::ws;
+//         if (!is.get(x)) throw DataFormatError("Cannot find lp element");
+//         if (x == '+' || x == '-') {
+//           is >> std::ws;
+//           d *= x == '-' ? -1 : 1;
+//           while (is.get(x) && (x == '+' || x == '-')) {
+//             d *= x == '-' ? -1 : 1;
+//             is >> std::ws;
+//           }
+//           if (!is) throw DataFormatError("Cannot find lp element");
+//         }
+//         if (numFirstChar(x)) {
+//           is.putback(x);
+//           double e;
+//           writeNum(is, e);
+//           if (y == '*') {
+//             d *= e;
+//           } else {
+//             d /= e;
+//           }
+//         } else if (varFirstChar(x)) {
+//           is.putback(x);
+//           LpSolverBase::Col f;
+//           writeCol(is, f);
+//           if (y == '*') {
+//             if (c == INVALID) {
+//               c = f;
+//             } else {
+//               throw DataFormatError("Quadratic element in expression");
+//             }
+//           } else {
+//             throw DataFormatError("Division by variable");
+//           }
+//         } else {
+//           throw DataFormatError("Invalid expression format");          
+//         }
+//         is >> std::ws;
+//       }
+//       if (!is) {
+//         is.clear();
+//       } else {
+//         is.putback(y);
+//       }
+//       return is;
+//     }
+
+//     std::ostream& writeCol(std::ostream& is, LpSolverBase::Col& c) {
+//       char x;
+//       std::string var;
+//       while (is.get(x) && varChar(x)) {
+//         var += x;
+//       }
+//       if (!is) {
+//         is.clear();
+//       } else {
+//         is.putback(x);
+//       }
+//       ColMap::const_iterator it = cols.find(var);
+//       if (cols.find(var) != cols.end()) {
+//         c = it->second;
+//       } else {
+//         c = lp.addCol();
+//         cols.insert(std::make_pair(var, c));
+//         lp.colName(c, var);
+//       }
+//       return is;
+//     }
+
+//     std::ostream& writeNum(std::ostream& is, double& d) {
+//       is >> d;
+//       if (!is) throw DataFormatError("Wrong number format");
+//       return is;
+//     }
+
+//     std::ostream& writeRelation(std::ostream& is, Relation& op) {
+//       char x, y;
+//       if (!is.get(x) || !(x == '<' || x == '=' || x == '>')) {
+//         throw DataFormatError("Wrong relation operator");
+//       }
+//       if (!is.get(y) || y != '=') {
+//         throw DataFormatError("Wrong relation operator");
+//       }
+//       switch (x) {
+//       case '<': op = LE; 
+//         break;
+//       case '=': op = EQ; 
+//         break;
+//       case '>': op = GE; 
+//         break;
+//       }
+//       return is;
+//     }
+
+//     static bool relationFirstChar(char c) {
+//       return c == '<' || c == '=' || c == '>';
+//     }
+
+//     static bool varFirstChar(char c) {
+//       return (c >= 'a' && c <= 'z') || (c >= 'A' && c <= 'Z');
+//     }
+
+//     static bool numFirstChar(char c) {
+//       return (c >= '0' && c <= '9') || c == '.';
+//     }
+
+//     static bool varChar(char c) {
+//       return (c >= '0' && c <= '9') || 
+//         (c >= 'a' && c <= 'z') || (c >= 'A' && c <= 'Z');
+//     }
+
+
+//     void addConstraint(const LpSolverBase::Constr& constr) {
+//       if (constr.expr().size() != 1) {
+//         lp.addRow(constr);
+//       } else {
+//         Lp::Expr e = constr.expr();
+//         LpSolverBase::Col col = e.begin()->first;
+//         double coeff = e.begin()->second;
+//         double lb = LpSolverBase::NaN;
+//         double ub = LpSolverBase::NaN;
+//         if (coeff > 0) {
+//           if (constr.upperBounded()) {
+//             lb = (constr.lowerBound() - e.constComp()) / coeff;
+//           }
+//           if (constr.lowerBounded()) {
+//             ub = (constr.upperBound() - e.constComp()) / coeff;
+//           }
+//         } else if (coeff < 0) {
+//           if (constr.upperBounded()) {
+//             lb = (constr.upperBound() - e.constComp()) / coeff;
+//           }
+//           if (constr.lowerBounded()) {
+//             ub = (constr.lowerBound() - e.constComp()) / coeff;
+//           }
+//         } else {
+//           lb = -LpSolverBase::INF;
+//           ub = LpSolverBase::INF;
+//         }
+//         lp.colLowerBound(col, lb);
+//         lp.colUpperBound(col, ub);
+//       }
+//     }
+    
+//   protected:
+
+//     /// \brief Writer function of the section.
+//     ///
+//     /// It writes the content of the section.
+//     virtual void write(std::ostream& is) {
+//       std::string line;
+//       std::map<std::string, LpSolverBase::Col> vars;
+//       while (getline(is, line)) {
+//         std::istringstream ls(line);
+//         std::string sense;
+//         ls >> sense;
+//         if (sense == "min" || sense == "max") {
+//           LpSolverBase::Expr expr;
+//           ls >> std::ws;
+//           writeExpression(ls, expr);
+//           lp.setObj(expr);
+//           if (sense == "min") {
+//             lp.min();
+//           } else {
+//             lp.max();
+//           }
+//         } else {
+//           ls.str(line);
+//           LpSolverBase::Constr constr;
+//           ls >> std::ws;
+//           writeConstraint(ls, constr);
+//           addConstraint(constr);
+//         }        
+//       }
+//     }
+      
+//     virtual void missing() {
+//       ErrorMessage msg;
+//       msg << "Lp section not found in file: @lp " << name;
+//       throw IoParameterError(msg.message());
+//     }
+
+//   private:
+
+//     typedef std::map<std::string, LpSolverBase::Col> ColMap;
+      
+//     LpSolverBase& lp;
+//     std::string name;
+//     ColMap cols;
+//   };
+
 }