lemon/lp_base.h
author Alpar Juttner <alpar@cs.elte.hu>
Fri, 13 Feb 2009 14:56:27 +0000
changeset 486 bc0a130bbb46
parent 470 81627fa1b007
child 487 06787db0ef5f
permissions -rw-r--r--
Fix in lp_base.cc targeting AIX/xlC compatibility
     1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
     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 #ifndef LEMON_LP_BASE_H
    20 #define LEMON_LP_BASE_H
    21 
    22 #include<iostream>
    23 #include<vector>
    24 #include<map>
    25 #include<limits>
    26 #include<lemon/math.h>
    27 
    28 #include<lemon/error.h>
    29 #include<lemon/assert.h>
    30 
    31 #include<lemon/core.h>
    32 #include<lemon/bits/solver_bits.h>
    33 
    34 ///\file
    35 ///\brief The interface of the LP solver interface.
    36 ///\ingroup lp_group
    37 namespace lemon {
    38 
    39   ///Common base class for LP and MIP solvers
    40 
    41   ///Usually this class is not used directly, please use one of the concrete
    42   ///implementations of the solver interface.
    43   ///\ingroup lp_group
    44   class LpBase {
    45 
    46   protected:
    47 
    48     _solver_bits::VarIndex rows;
    49     _solver_bits::VarIndex cols;
    50 
    51   public:
    52 
    53     ///Possible outcomes of an LP solving procedure
    54     enum SolveExitStatus {
    55       ///This means that the problem has been successfully solved: either
    56       ///an optimal solution has been found or infeasibility/unboundedness
    57       ///has been proved.
    58       SOLVED = 0,
    59       ///Any other case (including the case when some user specified
    60       ///limit has been exceeded)
    61       UNSOLVED = 1
    62     };
    63 
    64     ///Direction of the optimization
    65     enum Sense {
    66       /// Minimization
    67       MIN,
    68       /// Maximization
    69       MAX
    70     };
    71 
    72     ///The floating point type used by the solver
    73     typedef double Value;
    74     ///The infinity constant
    75     static const Value INF;
    76     ///The not a number constant
    77     static const Value NaN;
    78 
    79     friend class Col;
    80     friend class ColIt;
    81     friend class Row;
    82     friend class RowIt;
    83 
    84     ///Refer to a column of the LP.
    85 
    86     ///This type is used to refer to a column of the LP.
    87     ///
    88     ///Its value remains valid and correct even after the addition or erase of
    89     ///other columns.
    90     ///
    91     ///\note This class is similar to other Item types in LEMON, like
    92     ///Node and Arc types in digraph.
    93     class Col {
    94       friend class LpBase;
    95     protected:
    96       int _id;
    97       explicit Col(int id) : _id(id) {}
    98     public:
    99       typedef Value ExprValue;
   100       typedef True LpCol;
   101       /// Default constructor
   102       
   103       /// \warning The default constructor sets the Col to an
   104       /// undefined value.
   105       Col() {}
   106       /// Invalid constructor \& conversion.
   107       
   108       /// This constructor initializes the Col to be invalid.
   109       /// \sa Invalid for more details.      
   110       Col(const Invalid&) : _id(-1) {}
   111       /// Equality operator
   112 
   113       /// Two \ref Col "Col"s are equal if and only if they point to
   114       /// the same LP column or both are invalid.
   115       bool operator==(Col c) const  {return _id == c._id;}
   116       /// Inequality operator
   117 
   118       /// \sa operator==(Col c)
   119       ///
   120       bool operator!=(Col c) const  {return _id != c._id;}
   121       /// Artificial ordering operator.
   122 
   123       /// To allow the use of this object in std::map or similar
   124       /// associative container we require this.
   125       ///
   126       /// \note This operator only have to define some strict ordering of
   127       /// the items; this order has nothing to do with the iteration
   128       /// ordering of the items.
   129       bool operator<(Col c) const  {return _id < c._id;}
   130     };
   131 
   132     ///Iterator for iterate over the columns of an LP problem
   133 
   134     /// Its usage is quite simple, for example you can count the number
   135     /// of columns in an LP \c lp:
   136     ///\code
   137     /// int count=0;
   138     /// for (LpBase::ColIt c(lp); c!=INVALID; ++c) ++count;
   139     ///\endcode
   140     class ColIt : public Col {
   141       const LpBase *_solver;
   142     public:
   143       /// Default constructor
   144       
   145       /// \warning The default constructor sets the iterator
   146       /// to an undefined value.
   147       ColIt() {}
   148       /// Sets the iterator to the first Col
   149       
   150       /// Sets the iterator to the first Col.
   151       ///
   152       ColIt(const LpBase &solver) : _solver(&solver)
   153       {
   154         _solver->cols.firstItem(_id);
   155       }
   156       /// Invalid constructor \& conversion
   157       
   158       /// Initialize the iterator to be invalid.
   159       /// \sa Invalid for more details.
   160       ColIt(const Invalid&) : Col(INVALID) {}
   161       /// Next column
   162       
   163       /// Assign the iterator to the next column.
   164       ///
   165       ColIt &operator++()
   166       {
   167         _solver->cols.nextItem(_id);
   168         return *this;
   169       }
   170     };
   171 
   172     /// \brief Returns the ID of the column.
   173     static int id(const Col& col) { return col._id; }
   174     /// \brief Returns the column with the given ID.
   175     ///
   176     /// \pre The argument should be a valid column ID in the LP problem.
   177     static Col colFromId(int id) { return Col(id); }
   178 
   179     ///Refer to a row of the LP.
   180 
   181     ///This type is used to refer to a row of the LP.
   182     ///
   183     ///Its value remains valid and correct even after the addition or erase of
   184     ///other rows.
   185     ///
   186     ///\note This class is similar to other Item types in LEMON, like
   187     ///Node and Arc types in digraph.
   188     class Row {
   189       friend class LpBase;
   190     protected:
   191       int _id;
   192       explicit Row(int id) : _id(id) {}
   193     public:
   194       typedef Value ExprValue;
   195       typedef True LpRow;
   196       /// Default constructor
   197       
   198       /// \warning The default constructor sets the Row to an
   199       /// undefined value.
   200       Row() {}
   201       /// Invalid constructor \& conversion.
   202       
   203       /// This constructor initializes the Row to be invalid.
   204       /// \sa Invalid for more details.      
   205       Row(const Invalid&) : _id(-1) {}
   206       /// Equality operator
   207 
   208       /// Two \ref Row "Row"s are equal if and only if they point to
   209       /// the same LP row or both are invalid.
   210       bool operator==(Row r) const  {return _id == r._id;}
   211       /// Inequality operator
   212       
   213       /// \sa operator==(Row r)
   214       ///
   215       bool operator!=(Row r) const  {return _id != r._id;}
   216       /// Artificial ordering operator.
   217 
   218       /// To allow the use of this object in std::map or similar
   219       /// associative container we require this.
   220       ///
   221       /// \note This operator only have to define some strict ordering of
   222       /// the items; this order has nothing to do with the iteration
   223       /// ordering of the items.
   224       bool operator<(Row r) const  {return _id < r._id;}
   225     };
   226 
   227     ///Iterator for iterate over the rows of an LP problem
   228 
   229     /// Its usage is quite simple, for example you can count the number
   230     /// of rows in an LP \c lp:
   231     ///\code
   232     /// int count=0;
   233     /// for (LpBase::RowIt c(lp); c!=INVALID; ++c) ++count;
   234     ///\endcode
   235     class RowIt : public Row {
   236       const LpBase *_solver;
   237     public:
   238       /// Default constructor
   239       
   240       /// \warning The default constructor sets the iterator
   241       /// to an undefined value.
   242       RowIt() {}
   243       /// Sets the iterator to the first Row
   244       
   245       /// Sets the iterator to the first Row.
   246       ///
   247       RowIt(const LpBase &solver) : _solver(&solver)
   248       {
   249         _solver->rows.firstItem(_id);
   250       }
   251       /// Invalid constructor \& conversion
   252       
   253       /// Initialize the iterator to be invalid.
   254       /// \sa Invalid for more details.
   255       RowIt(const Invalid&) : Row(INVALID) {}
   256       /// Next row
   257       
   258       /// Assign the iterator to the next row.
   259       ///
   260       RowIt &operator++()
   261       {
   262         _solver->rows.nextItem(_id);
   263         return *this;
   264       }
   265     };
   266 
   267     /// \brief Returns the ID of the row.
   268     static int id(const Row& row) { return row._id; }
   269     /// \brief Returns the row with the given ID.
   270     ///
   271     /// \pre The argument should be a valid row ID in the LP problem.
   272     static Row rowFromId(int id) { return Row(id); }
   273 
   274   public:
   275 
   276     ///Linear expression of variables and a constant component
   277 
   278     ///This data structure stores a linear expression of the variables
   279     ///(\ref Col "Col"s) and also has a constant component.
   280     ///
   281     ///There are several ways to access and modify the contents of this
   282     ///container.
   283     ///\code
   284     ///e[v]=5;
   285     ///e[v]+=12;
   286     ///e.erase(v);
   287     ///\endcode
   288     ///or you can also iterate through its elements.
   289     ///\code
   290     ///double s=0;
   291     ///for(LpBase::Expr::ConstCoeffIt i(e);i!=INVALID;++i)
   292     ///  s+=*i * primal(i);
   293     ///\endcode
   294     ///(This code computes the primal value of the expression).
   295     ///- Numbers (<tt>double</tt>'s)
   296     ///and variables (\ref Col "Col"s) directly convert to an
   297     ///\ref Expr and the usual linear operations are defined, so
   298     ///\code
   299     ///v+w
   300     ///2*v-3.12*(v-w/2)+2
   301     ///v*2.1+(3*v+(v*12+w+6)*3)/2
   302     ///\endcode
   303     ///are valid expressions.
   304     ///The usual assignment operations are also defined.
   305     ///\code
   306     ///e=v+w;
   307     ///e+=2*v-3.12*(v-w/2)+2;
   308     ///e*=3.4;
   309     ///e/=5;
   310     ///\endcode
   311     ///- The constant member can be set and read by dereference
   312     ///  operator (unary *)
   313     ///
   314     ///\code
   315     ///*e=12;
   316     ///double c=*e;
   317     ///\endcode
   318     ///
   319     ///\sa Constr
   320     class Expr {
   321       friend class LpBase;
   322     public:
   323       /// The key type of the expression
   324       typedef LpBase::Col Key;
   325       /// The value type of the expression
   326       typedef LpBase::Value Value;
   327 
   328     protected:
   329       Value const_comp;
   330       std::map<int, Value> comps;
   331 
   332     public:
   333       typedef True SolverExpr;
   334       /// Default constructor
   335       
   336       /// Construct an empty expression, the coefficients and
   337       /// the constant component are initialized to zero.
   338       Expr() : const_comp(0) {}
   339       /// Construct an expression from a column
   340 
   341       /// Construct an expression, which has a term with \c c variable
   342       /// and 1.0 coefficient.
   343       Expr(const Col &c) : const_comp(0) {
   344         typedef std::map<int, Value>::value_type pair_type;
   345         comps.insert(pair_type(id(c), 1));
   346       }
   347       /// Construct an expression from a constant
   348 
   349       /// Construct an expression, which's constant component is \c v.
   350       ///
   351       Expr(const Value &v) : const_comp(v) {}
   352       /// Returns the coefficient of the column
   353       Value operator[](const Col& c) const {
   354         std::map<int, Value>::const_iterator it=comps.find(id(c));
   355         if (it != comps.end()) {
   356           return it->second;
   357         } else {
   358           return 0;
   359         }
   360       }
   361       /// Returns the coefficient of the column
   362       Value& operator[](const Col& c) {
   363         return comps[id(c)];
   364       }
   365       /// Sets the coefficient of the column
   366       void set(const Col &c, const Value &v) {
   367         if (v != 0.0) {
   368           typedef std::map<int, Value>::value_type pair_type;
   369           comps.insert(pair_type(id(c), v));
   370         } else {
   371           comps.erase(id(c));
   372         }
   373       }
   374       /// Returns the constant component of the expression
   375       Value& operator*() { return const_comp; }
   376       /// Returns the constant component of the expression
   377       const Value& operator*() const { return const_comp; }
   378       /// \brief Removes the coefficients which's absolute value does
   379       /// not exceed \c epsilon. It also sets to zero the constant
   380       /// component, if it does not exceed epsilon in absolute value.
   381       void simplify(Value epsilon = 0.0) {
   382         std::map<int, Value>::iterator it=comps.begin();
   383         while (it != comps.end()) {
   384           std::map<int, Value>::iterator jt=it;
   385           ++jt;
   386           if (std::fabs((*it).second) <= epsilon) comps.erase(it);
   387           it=jt;
   388         }
   389         if (std::fabs(const_comp) <= epsilon) const_comp = 0;
   390       }
   391 
   392       void simplify(Value epsilon = 0.0) const {
   393         const_cast<Expr*>(this)->simplify(epsilon);
   394       }
   395 
   396       ///Sets all coefficients and the constant component to 0.
   397       void clear() {
   398         comps.clear();
   399         const_comp=0;
   400       }
   401 
   402       ///Compound assignment
   403       Expr &operator+=(const Expr &e) {
   404         for (std::map<int, Value>::const_iterator it=e.comps.begin();
   405              it!=e.comps.end(); ++it)
   406           comps[it->first]+=it->second;
   407         const_comp+=e.const_comp;
   408         return *this;
   409       }
   410       ///Compound assignment
   411       Expr &operator-=(const Expr &e) {
   412         for (std::map<int, Value>::const_iterator it=e.comps.begin();
   413              it!=e.comps.end(); ++it)
   414           comps[it->first]-=it->second;
   415         const_comp-=e.const_comp;
   416         return *this;
   417       }
   418       ///Multiply with a constant
   419       Expr &operator*=(const Value &v) {
   420         for (std::map<int, Value>::iterator it=comps.begin();
   421              it!=comps.end(); ++it)
   422           it->second*=v;
   423         const_comp*=v;
   424         return *this;
   425       }
   426       ///Division with a constant
   427       Expr &operator/=(const Value &c) {
   428         for (std::map<int, Value>::iterator it=comps.begin();
   429              it!=comps.end(); ++it)
   430           it->second/=c;
   431         const_comp/=c;
   432         return *this;
   433       }
   434 
   435       ///Iterator over the expression
   436       
   437       ///The iterator iterates over the terms of the expression. 
   438       /// 
   439       ///\code
   440       ///double s=0;
   441       ///for(LpBase::Expr::CoeffIt i(e);i!=INVALID;++i)
   442       ///  s+= *i * primal(i);
   443       ///\endcode
   444       class CoeffIt {
   445       private:
   446 
   447         std::map<int, Value>::iterator _it, _end;
   448 
   449       public:
   450 
   451         /// Sets the iterator to the first term
   452         
   453         /// Sets the iterator to the first term of the expression.
   454         ///
   455         CoeffIt(Expr& e)
   456           : _it(e.comps.begin()), _end(e.comps.end()){}
   457 
   458         /// Convert the iterator to the column of the term
   459         operator Col() const {
   460           return colFromId(_it->first);
   461         }
   462 
   463         /// Returns the coefficient of the term
   464         Value& operator*() { return _it->second; }
   465 
   466         /// Returns the coefficient of the term
   467         const Value& operator*() const { return _it->second; }
   468         /// Next term
   469         
   470         /// Assign the iterator to the next term.
   471         ///
   472         CoeffIt& operator++() { ++_it; return *this; }
   473 
   474         /// Equality operator
   475         bool operator==(Invalid) const { return _it == _end; }
   476         /// Inequality operator
   477         bool operator!=(Invalid) const { return _it != _end; }
   478       };
   479 
   480       /// Const iterator over the expression
   481       
   482       ///The iterator iterates over the terms of the expression. 
   483       /// 
   484       ///\code
   485       ///double s=0;
   486       ///for(LpBase::Expr::ConstCoeffIt i(e);i!=INVALID;++i)
   487       ///  s+=*i * primal(i);
   488       ///\endcode
   489       class ConstCoeffIt {
   490       private:
   491 
   492         std::map<int, Value>::const_iterator _it, _end;
   493 
   494       public:
   495 
   496         /// Sets the iterator to the first term
   497         
   498         /// Sets the iterator to the first term of the expression.
   499         ///
   500         ConstCoeffIt(const Expr& e)
   501           : _it(e.comps.begin()), _end(e.comps.end()){}
   502 
   503         /// Convert the iterator to the column of the term
   504         operator Col() const {
   505           return colFromId(_it->first);
   506         }
   507 
   508         /// Returns the coefficient of the term
   509         const Value& operator*() const { return _it->second; }
   510 
   511         /// Next term
   512         
   513         /// Assign the iterator to the next term.
   514         ///
   515         ConstCoeffIt& operator++() { ++_it; return *this; }
   516 
   517         /// Equality operator
   518         bool operator==(Invalid) const { return _it == _end; }
   519         /// Inequality operator
   520         bool operator!=(Invalid) const { return _it != _end; }
   521       };
   522 
   523     };
   524 
   525     ///Linear constraint
   526 
   527     ///This data stucture represents a linear constraint in the LP.
   528     ///Basically it is a linear expression with a lower or an upper bound
   529     ///(or both). These parts of the constraint can be obtained by the member
   530     ///functions \ref expr(), \ref lowerBound() and \ref upperBound(),
   531     ///respectively.
   532     ///There are two ways to construct a constraint.
   533     ///- You can set the linear expression and the bounds directly
   534     ///  by the functions above.
   535     ///- The operators <tt>\<=</tt>, <tt>==</tt> and  <tt>\>=</tt>
   536     ///  are defined between expressions, or even between constraints whenever
   537     ///  it makes sense. Therefore if \c e and \c f are linear expressions and
   538     ///  \c s and \c t are numbers, then the followings are valid expressions
   539     ///  and thus they can be used directly e.g. in \ref addRow() whenever
   540     ///  it makes sense.
   541     ///\code
   542     ///  e<=s
   543     ///  e<=f
   544     ///  e==f
   545     ///  s<=e<=t
   546     ///  e>=t
   547     ///\endcode
   548     ///\warning The validity of a constraint is checked only at run
   549     ///time, so e.g. \ref addRow(<tt>x[1]\<=x[2]<=5</tt>) will
   550     ///compile, but will fail an assertion.
   551     class Constr
   552     {
   553     public:
   554       typedef LpBase::Expr Expr;
   555       typedef Expr::Key Key;
   556       typedef Expr::Value Value;
   557 
   558     protected:
   559       Expr _expr;
   560       Value _lb,_ub;
   561     public:
   562       ///\e
   563       Constr() : _expr(), _lb(NaN), _ub(NaN) {}
   564       ///\e
   565       Constr(Value lb, const Expr &e, Value ub) :
   566         _expr(e), _lb(lb), _ub(ub) {}
   567       Constr(const Expr &e) :
   568         _expr(e), _lb(NaN), _ub(NaN) {}
   569       ///\e
   570       void clear()
   571       {
   572         _expr.clear();
   573         _lb=_ub=NaN;
   574       }
   575 
   576       ///Reference to the linear expression
   577       Expr &expr() { return _expr; }
   578       ///Cont reference to the linear expression
   579       const Expr &expr() const { return _expr; }
   580       ///Reference to the lower bound.
   581 
   582       ///\return
   583       ///- \ref INF "INF": the constraint is lower unbounded.
   584       ///- \ref NaN "NaN": lower bound has not been set.
   585       ///- finite number: the lower bound
   586       Value &lowerBound() { return _lb; }
   587       ///The const version of \ref lowerBound()
   588       const Value &lowerBound() const { return _lb; }
   589       ///Reference to the upper bound.
   590 
   591       ///\return
   592       ///- \ref INF "INF": the constraint is upper unbounded.
   593       ///- \ref NaN "NaN": upper bound has not been set.
   594       ///- finite number: the upper bound
   595       Value &upperBound() { return _ub; }
   596       ///The const version of \ref upperBound()
   597       const Value &upperBound() const { return _ub; }
   598       ///Is the constraint lower bounded?
   599       bool lowerBounded() const {
   600         return _lb != -INF && !isnan(_lb);
   601       }
   602       ///Is the constraint upper bounded?
   603       bool upperBounded() const {
   604         return _ub != INF && !isnan(_ub);
   605       }
   606 
   607     };
   608 
   609     ///Linear expression of rows
   610 
   611     ///This data structure represents a column of the matrix,
   612     ///thas is it strores a linear expression of the dual variables
   613     ///(\ref Row "Row"s).
   614     ///
   615     ///There are several ways to access and modify the contents of this
   616     ///container.
   617     ///\code
   618     ///e[v]=5;
   619     ///e[v]+=12;
   620     ///e.erase(v);
   621     ///\endcode
   622     ///or you can also iterate through its elements.
   623     ///\code
   624     ///double s=0;
   625     ///for(LpBase::DualExpr::ConstCoeffIt i(e);i!=INVALID;++i)
   626     ///  s+=*i;
   627     ///\endcode
   628     ///(This code computes the sum of all coefficients).
   629     ///- Numbers (<tt>double</tt>'s)
   630     ///and variables (\ref Row "Row"s) directly convert to an
   631     ///\ref DualExpr and the usual linear operations are defined, so
   632     ///\code
   633     ///v+w
   634     ///2*v-3.12*(v-w/2)
   635     ///v*2.1+(3*v+(v*12+w)*3)/2
   636     ///\endcode
   637     ///are valid \ref DualExpr dual expressions.
   638     ///The usual assignment operations are also defined.
   639     ///\code
   640     ///e=v+w;
   641     ///e+=2*v-3.12*(v-w/2);
   642     ///e*=3.4;
   643     ///e/=5;
   644     ///\endcode
   645     ///
   646     ///\sa Expr
   647     class DualExpr {
   648       friend class LpBase;
   649     public:
   650       /// The key type of the expression
   651       typedef LpBase::Row Key;
   652       /// The value type of the expression
   653       typedef LpBase::Value Value;
   654 
   655     protected:
   656       std::map<int, Value> comps;
   657 
   658     public:
   659       typedef True SolverExpr;
   660       /// Default constructor
   661       
   662       /// Construct an empty expression, the coefficients are
   663       /// initialized to zero.
   664       DualExpr() {}
   665       /// Construct an expression from a row
   666 
   667       /// Construct an expression, which has a term with \c r dual
   668       /// variable and 1.0 coefficient.
   669       DualExpr(const Row &r) {
   670         typedef std::map<int, Value>::value_type pair_type;
   671         comps.insert(pair_type(id(r), 1));
   672       }
   673       /// Returns the coefficient of the row
   674       Value operator[](const Row& r) const {
   675         std::map<int, Value>::const_iterator it = comps.find(id(r));
   676         if (it != comps.end()) {
   677           return it->second;
   678         } else {
   679           return 0;
   680         }
   681       }
   682       /// Returns the coefficient of the row
   683       Value& operator[](const Row& r) {
   684         return comps[id(r)];
   685       }
   686       /// Sets the coefficient of the row
   687       void set(const Row &r, const Value &v) {
   688         if (v != 0.0) {
   689           typedef std::map<int, Value>::value_type pair_type;
   690           comps.insert(pair_type(id(r), v));
   691         } else {
   692           comps.erase(id(r));
   693         }
   694       }
   695       /// \brief Removes the coefficients which's absolute value does
   696       /// not exceed \c epsilon. 
   697       void simplify(Value epsilon = 0.0) {
   698         std::map<int, Value>::iterator it=comps.begin();
   699         while (it != comps.end()) {
   700           std::map<int, Value>::iterator jt=it;
   701           ++jt;
   702           if (std::fabs((*it).second) <= epsilon) comps.erase(it);
   703           it=jt;
   704         }
   705       }
   706 
   707       void simplify(Value epsilon = 0.0) const {
   708         const_cast<DualExpr*>(this)->simplify(epsilon);
   709       }
   710 
   711       ///Sets all coefficients to 0.
   712       void clear() {
   713         comps.clear();
   714       }
   715       ///Compound assignment
   716       DualExpr &operator+=(const DualExpr &e) {
   717         for (std::map<int, Value>::const_iterator it=e.comps.begin();
   718              it!=e.comps.end(); ++it)
   719           comps[it->first]+=it->second;
   720         return *this;
   721       }
   722       ///Compound assignment
   723       DualExpr &operator-=(const DualExpr &e) {
   724         for (std::map<int, Value>::const_iterator it=e.comps.begin();
   725              it!=e.comps.end(); ++it)
   726           comps[it->first]-=it->second;
   727         return *this;
   728       }
   729       ///Multiply with a constant
   730       DualExpr &operator*=(const Value &v) {
   731         for (std::map<int, Value>::iterator it=comps.begin();
   732              it!=comps.end(); ++it)
   733           it->second*=v;
   734         return *this;
   735       }
   736       ///Division with a constant
   737       DualExpr &operator/=(const Value &v) {
   738         for (std::map<int, Value>::iterator it=comps.begin();
   739              it!=comps.end(); ++it)
   740           it->second/=v;
   741         return *this;
   742       }
   743 
   744       ///Iterator over the expression
   745       
   746       ///The iterator iterates over the terms of the expression. 
   747       /// 
   748       ///\code
   749       ///double s=0;
   750       ///for(LpBase::DualExpr::CoeffIt i(e);i!=INVALID;++i)
   751       ///  s+= *i * dual(i);
   752       ///\endcode
   753       class CoeffIt {
   754       private:
   755 
   756         std::map<int, Value>::iterator _it, _end;
   757 
   758       public:
   759 
   760         /// Sets the iterator to the first term
   761         
   762         /// Sets the iterator to the first term of the expression.
   763         ///
   764         CoeffIt(DualExpr& e)
   765           : _it(e.comps.begin()), _end(e.comps.end()){}
   766 
   767         /// Convert the iterator to the row of the term
   768         operator Row() const {
   769           return rowFromId(_it->first);
   770         }
   771 
   772         /// Returns the coefficient of the term
   773         Value& operator*() { return _it->second; }
   774 
   775         /// Returns the coefficient of the term
   776         const Value& operator*() const { return _it->second; }
   777 
   778         /// Next term
   779         
   780         /// Assign the iterator to the next term.
   781         ///
   782         CoeffIt& operator++() { ++_it; return *this; }
   783 
   784         /// Equality operator
   785         bool operator==(Invalid) const { return _it == _end; }
   786         /// Inequality operator
   787         bool operator!=(Invalid) const { return _it != _end; }
   788       };
   789 
   790       ///Iterator over the expression
   791       
   792       ///The iterator iterates over the terms of the expression. 
   793       /// 
   794       ///\code
   795       ///double s=0;
   796       ///for(LpBase::DualExpr::ConstCoeffIt i(e);i!=INVALID;++i)
   797       ///  s+= *i * dual(i);
   798       ///\endcode
   799       class ConstCoeffIt {
   800       private:
   801 
   802         std::map<int, Value>::const_iterator _it, _end;
   803 
   804       public:
   805 
   806         /// Sets the iterator to the first term
   807         
   808         /// Sets the iterator to the first term of the expression.
   809         ///
   810         ConstCoeffIt(const DualExpr& e)
   811           : _it(e.comps.begin()), _end(e.comps.end()){}
   812 
   813         /// Convert the iterator to the row of the term
   814         operator Row() const {
   815           return rowFromId(_it->first);
   816         }
   817 
   818         /// Returns the coefficient of the term
   819         const Value& operator*() const { return _it->second; }
   820 
   821         /// Next term
   822         
   823         /// Assign the iterator to the next term.
   824         ///
   825         ConstCoeffIt& operator++() { ++_it; return *this; }
   826 
   827         /// Equality operator
   828         bool operator==(Invalid) const { return _it == _end; }
   829         /// Inequality operator
   830         bool operator!=(Invalid) const { return _it != _end; }
   831       };
   832     };
   833 
   834 
   835   protected:
   836 
   837     class InsertIterator {
   838     private:
   839 
   840       std::map<int, Value>& _host;
   841       const _solver_bits::VarIndex& _index;
   842 
   843     public:
   844 
   845       typedef std::output_iterator_tag iterator_category;
   846       typedef void difference_type;
   847       typedef void value_type;
   848       typedef void reference;
   849       typedef void pointer;
   850 
   851       InsertIterator(std::map<int, Value>& host,
   852                    const _solver_bits::VarIndex& index)
   853         : _host(host), _index(index) {}
   854 
   855       InsertIterator& operator=(const std::pair<int, Value>& value) {
   856         typedef std::map<int, Value>::value_type pair_type;
   857         _host.insert(pair_type(_index[value.first], value.second));
   858         return *this;
   859       }
   860 
   861       InsertIterator& operator*() { return *this; }
   862       InsertIterator& operator++() { return *this; }
   863       InsertIterator operator++(int) { return *this; }
   864 
   865     };
   866 
   867     class ExprIterator {
   868     private:
   869       std::map<int, Value>::const_iterator _host_it;
   870       const _solver_bits::VarIndex& _index;
   871     public:
   872 
   873       typedef std::bidirectional_iterator_tag iterator_category;
   874       typedef std::ptrdiff_t difference_type;
   875       typedef const std::pair<int, Value> value_type;
   876       typedef value_type reference;
   877 
   878       class pointer {
   879       public:
   880         pointer(value_type& _value) : value(_value) {}
   881         value_type* operator->() { return &value; }
   882       private:
   883         value_type value;
   884       };
   885 
   886       ExprIterator(const std::map<int, Value>::const_iterator& host_it,
   887                    const _solver_bits::VarIndex& index)
   888         : _host_it(host_it), _index(index) {}
   889 
   890       reference operator*() {
   891         return std::make_pair(_index(_host_it->first), _host_it->second);
   892       }
   893 
   894       pointer operator->() {
   895         return pointer(operator*());
   896       }
   897 
   898       ExprIterator& operator++() { ++_host_it; return *this; }
   899       ExprIterator operator++(int) {
   900         ExprIterator tmp(*this); ++_host_it; return tmp;
   901       }
   902 
   903       ExprIterator& operator--() { --_host_it; return *this; }
   904       ExprIterator operator--(int) {
   905         ExprIterator tmp(*this); --_host_it; return tmp;
   906       }
   907 
   908       bool operator==(const ExprIterator& it) const {
   909         return _host_it == it._host_it;
   910       }
   911 
   912       bool operator!=(const ExprIterator& it) const {
   913         return _host_it != it._host_it;
   914       }
   915 
   916     };
   917 
   918   protected:
   919 
   920     //Abstract virtual functions
   921     virtual LpBase* _newSolver() const = 0;
   922     virtual LpBase* _cloneSolver() const = 0;
   923 
   924     virtual int _addColId(int col) { return cols.addIndex(col); }
   925     virtual int _addRowId(int row) { return rows.addIndex(row); }
   926 
   927     virtual void _eraseColId(int col) { cols.eraseIndex(col); }
   928     virtual void _eraseRowId(int row) { rows.eraseIndex(row); }
   929 
   930     virtual int _addCol() = 0;
   931     virtual int _addRow() = 0;
   932 
   933     virtual void _eraseCol(int col) = 0;
   934     virtual void _eraseRow(int row) = 0;
   935 
   936     virtual void _getColName(int col, std::string& name) const = 0;
   937     virtual void _setColName(int col, const std::string& name) = 0;
   938     virtual int _colByName(const std::string& name) const = 0;
   939 
   940     virtual void _getRowName(int row, std::string& name) const = 0;
   941     virtual void _setRowName(int row, const std::string& name) = 0;
   942     virtual int _rowByName(const std::string& name) const = 0;
   943 
   944     virtual void _setRowCoeffs(int i, ExprIterator b, ExprIterator e) = 0;
   945     virtual void _getRowCoeffs(int i, InsertIterator b) const = 0;
   946 
   947     virtual void _setColCoeffs(int i, ExprIterator b, ExprIterator e) = 0;
   948     virtual void _getColCoeffs(int i, InsertIterator b) const = 0;
   949 
   950     virtual void _setCoeff(int row, int col, Value value) = 0;
   951     virtual Value _getCoeff(int row, int col) const = 0;
   952 
   953     virtual void _setColLowerBound(int i, Value value) = 0;
   954     virtual Value _getColLowerBound(int i) const = 0;
   955 
   956     virtual void _setColUpperBound(int i, Value value) = 0;
   957     virtual Value _getColUpperBound(int i) const = 0;
   958 
   959     virtual void _setRowLowerBound(int i, Value value) = 0;
   960     virtual Value _getRowLowerBound(int i) const = 0;
   961 
   962     virtual void _setRowUpperBound(int i, Value value) = 0;
   963     virtual Value _getRowUpperBound(int i) const = 0;
   964 
   965     virtual void _setObjCoeffs(ExprIterator b, ExprIterator e) = 0;
   966     virtual void _getObjCoeffs(InsertIterator b) const = 0;
   967 
   968     virtual void _setObjCoeff(int i, Value obj_coef) = 0;
   969     virtual Value _getObjCoeff(int i) const = 0;
   970 
   971     virtual void _setSense(Sense) = 0;
   972     virtual Sense _getSense() const = 0;
   973 
   974     virtual void _clear() = 0;
   975 
   976     virtual const char* _solverName() const = 0;
   977 
   978     //Own protected stuff
   979 
   980     //Constant component of the objective function
   981     Value obj_const_comp;
   982 
   983     LpBase() : rows(), cols(), obj_const_comp(0) {}
   984 
   985   public:
   986 
   987     /// Virtual destructor
   988     virtual ~LpBase() {}
   989 
   990     ///Creates a new LP problem
   991     LpBase* newSolver() {return _newSolver();}
   992     ///Makes a copy of the LP problem
   993     LpBase* cloneSolver() {return _cloneSolver();}
   994 
   995     ///Gives back the name of the solver.
   996     const char* solverName() const {return _solverName();}
   997 
   998     ///\name Build up and modify the LP
   999 
  1000     ///@{
  1001 
  1002     ///Add a new empty column (i.e a new variable) to the LP
  1003     Col addCol() { Col c; c._id = _addColId(_addCol()); return c;}
  1004 
  1005     ///\brief Adds several new columns (i.e variables) at once
  1006     ///
  1007     ///This magic function takes a container as its argument and fills
  1008     ///its elements with new columns (i.e. variables)
  1009     ///\param t can be
  1010     ///- a standard STL compatible iterable container with
  1011     ///\ref Col as its \c values_type like
  1012     ///\code
  1013     ///std::vector<LpBase::Col>
  1014     ///std::list<LpBase::Col>
  1015     ///\endcode
  1016     ///- a standard STL compatible iterable container with
  1017     ///\ref Col as its \c mapped_type like
  1018     ///\code
  1019     ///std::map<AnyType,LpBase::Col>
  1020     ///\endcode
  1021     ///- an iterable lemon \ref concepts::WriteMap "write map" like
  1022     ///\code
  1023     ///ListGraph::NodeMap<LpBase::Col>
  1024     ///ListGraph::ArcMap<LpBase::Col>
  1025     ///\endcode
  1026     ///\return The number of the created column.
  1027 #ifdef DOXYGEN
  1028     template<class T>
  1029     int addColSet(T &t) { return 0;}
  1030 #else
  1031     template<class T>
  1032     typename enable_if<typename T::value_type::LpCol,int>::type
  1033     addColSet(T &t,dummy<0> = 0) {
  1034       int s=0;
  1035       for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addCol();s++;}
  1036       return s;
  1037     }
  1038     template<class T>
  1039     typename enable_if<typename T::value_type::second_type::LpCol,
  1040                        int>::type
  1041     addColSet(T &t,dummy<1> = 1) {
  1042       int s=0;
  1043       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1044         i->second=addCol();
  1045         s++;
  1046       }
  1047       return s;
  1048     }
  1049     template<class T>
  1050     typename enable_if<typename T::MapIt::Value::LpCol,
  1051                        int>::type
  1052     addColSet(T &t,dummy<2> = 2) {
  1053       int s=0;
  1054       for(typename T::MapIt i(t); i!=INVALID; ++i)
  1055         {
  1056           i.set(addCol());
  1057           s++;
  1058         }
  1059       return s;
  1060     }
  1061 #endif
  1062 
  1063     ///Set a column (i.e a dual constraint) of the LP
  1064 
  1065     ///\param c is the column to be modified
  1066     ///\param e is a dual linear expression (see \ref DualExpr)
  1067     ///a better one.
  1068     void col(Col c, const DualExpr &e) {
  1069       e.simplify();
  1070       _setColCoeffs(cols(id(c)), ExprIterator(e.comps.begin(), rows),
  1071                     ExprIterator(e.comps.end(), rows));
  1072     }
  1073 
  1074     ///Get a column (i.e a dual constraint) of the LP
  1075 
  1076     ///\param c is the column to get
  1077     ///\return the dual expression associated to the column
  1078     DualExpr col(Col c) const {
  1079       DualExpr e;
  1080       _getColCoeffs(cols(id(c)), InsertIterator(e.comps, rows));
  1081       return e;
  1082     }
  1083 
  1084     ///Add a new column to the LP
  1085 
  1086     ///\param e is a dual linear expression (see \ref DualExpr)
  1087     ///\param o is the corresponding component of the objective
  1088     ///function. It is 0 by default.
  1089     ///\return The created column.
  1090     Col addCol(const DualExpr &e, Value o = 0) {
  1091       Col c=addCol();
  1092       col(c,e);
  1093       objCoeff(c,o);
  1094       return c;
  1095     }
  1096 
  1097     ///Add a new empty row (i.e a new constraint) to the LP
  1098 
  1099     ///This function adds a new empty row (i.e a new constraint) to the LP.
  1100     ///\return The created row
  1101     Row addRow() { Row r; r._id = _addRowId(_addRow()); return r;}
  1102 
  1103     ///\brief Add several new rows (i.e constraints) at once
  1104     ///
  1105     ///This magic function takes a container as its argument and fills
  1106     ///its elements with new row (i.e. variables)
  1107     ///\param t can be
  1108     ///- a standard STL compatible iterable container with
  1109     ///\ref Row as its \c values_type like
  1110     ///\code
  1111     ///std::vector<LpBase::Row>
  1112     ///std::list<LpBase::Row>
  1113     ///\endcode
  1114     ///- a standard STL compatible iterable container with
  1115     ///\ref Row as its \c mapped_type like
  1116     ///\code
  1117     ///std::map<AnyType,LpBase::Row>
  1118     ///\endcode
  1119     ///- an iterable lemon \ref concepts::WriteMap "write map" like
  1120     ///\code
  1121     ///ListGraph::NodeMap<LpBase::Row>
  1122     ///ListGraph::ArcMap<LpBase::Row>
  1123     ///\endcode
  1124     ///\return The number of rows created.
  1125 #ifdef DOXYGEN
  1126     template<class T>
  1127     int addRowSet(T &t) { return 0;}
  1128 #else
  1129     template<class T>
  1130     typename enable_if<typename T::value_type::LpRow,int>::type
  1131     addRowSet(T &t, dummy<0> = 0) {
  1132       int s=0;
  1133       for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addRow();s++;}
  1134       return s;
  1135     }
  1136     template<class T>
  1137     typename enable_if<typename T::value_type::second_type::LpRow, int>::type
  1138     addRowSet(T &t, dummy<1> = 1) {
  1139       int s=0;
  1140       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1141         i->second=addRow();
  1142         s++;
  1143       }
  1144       return s;
  1145     }
  1146     template<class T>
  1147     typename enable_if<typename T::MapIt::Value::LpRow, int>::type
  1148     addRowSet(T &t, dummy<2> = 2) {
  1149       int s=0;
  1150       for(typename T::MapIt i(t); i!=INVALID; ++i)
  1151         {
  1152           i.set(addRow());
  1153           s++;
  1154         }
  1155       return s;
  1156     }
  1157 #endif
  1158 
  1159     ///Set a row (i.e a constraint) of the LP
  1160 
  1161     ///\param r is the row to be modified
  1162     ///\param l is lower bound (-\ref INF means no bound)
  1163     ///\param e is a linear expression (see \ref Expr)
  1164     ///\param u is the upper bound (\ref INF means no bound)
  1165     void row(Row r, Value l, const Expr &e, Value u) {
  1166       e.simplify();
  1167       _setRowCoeffs(rows(id(r)), ExprIterator(e.comps.begin(), cols),
  1168                     ExprIterator(e.comps.end(), cols));
  1169       _setRowLowerBound(rows(id(r)),l - *e);
  1170       _setRowUpperBound(rows(id(r)),u - *e);
  1171     }
  1172 
  1173     ///Set a row (i.e a constraint) of the LP
  1174 
  1175     ///\param r is the row to be modified
  1176     ///\param c is a linear expression (see \ref Constr)
  1177     void row(Row r, const Constr &c) {
  1178       row(r, c.lowerBounded()?c.lowerBound():-INF,
  1179           c.expr(), c.upperBounded()?c.upperBound():INF);
  1180     }
  1181 
  1182 
  1183     ///Get a row (i.e a constraint) of the LP
  1184 
  1185     ///\param r is the row to get
  1186     ///\return the expression associated to the row
  1187     Expr row(Row r) const {
  1188       Expr e;
  1189       _getRowCoeffs(rows(id(r)), InsertIterator(e.comps, cols));
  1190       return e;
  1191     }
  1192 
  1193     ///Add a new row (i.e a new constraint) to the LP
  1194 
  1195     ///\param l is the lower bound (-\ref INF means no bound)
  1196     ///\param e is a linear expression (see \ref Expr)
  1197     ///\param u is the upper bound (\ref INF means no bound)
  1198     ///\return The created row.
  1199     Row addRow(Value l,const Expr &e, Value u) {
  1200       Row r=addRow();
  1201       row(r,l,e,u);
  1202       return r;
  1203     }
  1204 
  1205     ///Add a new row (i.e a new constraint) to the LP
  1206 
  1207     ///\param c is a linear expression (see \ref Constr)
  1208     ///\return The created row.
  1209     Row addRow(const Constr &c) {
  1210       Row r=addRow();
  1211       row(r,c);
  1212       return r;
  1213     }
  1214     ///Erase a column (i.e a variable) from the LP
  1215 
  1216     ///\param c is the column to be deleted
  1217     void erase(Col c) {
  1218       _eraseCol(cols(id(c)));
  1219       _eraseColId(cols(id(c)));
  1220     }
  1221     ///Erase a row (i.e a constraint) from the LP
  1222 
  1223     ///\param r is the row to be deleted
  1224     void erase(Row r) {
  1225       _eraseRow(rows(id(r)));
  1226       _eraseRowId(rows(id(r)));
  1227     }
  1228 
  1229     /// Get the name of a column
  1230 
  1231     ///\param c is the coresponding column
  1232     ///\return The name of the colunm
  1233     std::string colName(Col c) const {
  1234       std::string name;
  1235       _getColName(cols(id(c)), name);
  1236       return name;
  1237     }
  1238 
  1239     /// Set the name of a column
  1240 
  1241     ///\param c is the coresponding column
  1242     ///\param name The name to be given
  1243     void colName(Col c, const std::string& name) {
  1244       _setColName(cols(id(c)), name);
  1245     }
  1246 
  1247     /// Get the column by its name
  1248 
  1249     ///\param name The name of the column
  1250     ///\return the proper column or \c INVALID
  1251     Col colByName(const std::string& name) const {
  1252       int k = _colByName(name);
  1253       return k != -1 ? Col(cols[k]) : Col(INVALID);
  1254     }
  1255 
  1256     /// Get the name of a row
  1257 
  1258     ///\param r is the coresponding row
  1259     ///\return The name of the row
  1260     std::string rowName(Row r) const {
  1261       std::string name;
  1262       _getRowName(rows(id(r)), name);
  1263       return name;
  1264     }
  1265 
  1266     /// Set the name of a row
  1267 
  1268     ///\param r is the coresponding row
  1269     ///\param name The name to be given
  1270     void rowName(Row r, const std::string& name) {
  1271       _setRowName(rows(id(r)), name);
  1272     }
  1273 
  1274     /// Get the row by its name
  1275 
  1276     ///\param name The name of the row
  1277     ///\return the proper row or \c INVALID
  1278     Row rowByName(const std::string& name) const {
  1279       int k = _rowByName(name);
  1280       return k != -1 ? Row(rows[k]) : Row(INVALID);
  1281     }
  1282 
  1283     /// Set an element of the coefficient matrix of the LP
  1284 
  1285     ///\param r is the row of the element to be modified
  1286     ///\param c is the column of the element to be modified
  1287     ///\param val is the new value of the coefficient
  1288     void coeff(Row r, Col c, Value val) {
  1289       _setCoeff(rows(id(r)),cols(id(c)), val);
  1290     }
  1291 
  1292     /// Get an element of the coefficient matrix of the LP
  1293 
  1294     ///\param r is the row of the element
  1295     ///\param c is the column of the element
  1296     ///\return the corresponding coefficient
  1297     Value coeff(Row r, Col c) const {
  1298       return _getCoeff(rows(id(r)),cols(id(c)));
  1299     }
  1300 
  1301     /// Set the lower bound of a column (i.e a variable)
  1302 
  1303     /// The lower bound of a variable (column) has to be given by an
  1304     /// extended number of type Value, i.e. a finite number of type
  1305     /// Value or -\ref INF.
  1306     void colLowerBound(Col c, Value value) {
  1307       _setColLowerBound(cols(id(c)),value);
  1308     }
  1309 
  1310     /// Get the lower bound of a column (i.e a variable)
  1311 
  1312     /// This function returns the lower bound for column (variable) \c c
  1313     /// (this might be -\ref INF as well).
  1314     ///\return The lower bound for column \c c
  1315     Value colLowerBound(Col c) const {
  1316       return _getColLowerBound(cols(id(c)));
  1317     }
  1318 
  1319     ///\brief Set the lower bound of  several columns
  1320     ///(i.e variables) at once
  1321     ///
  1322     ///This magic function takes a container as its argument
  1323     ///and applies the function on all of its elements.
  1324     ///The lower bound of a variable (column) has to be given by an
  1325     ///extended number of type Value, i.e. a finite number of type
  1326     ///Value or -\ref INF.
  1327 #ifdef DOXYGEN
  1328     template<class T>
  1329     void colLowerBound(T &t, Value value) { return 0;}
  1330 #else
  1331     template<class T>
  1332     typename enable_if<typename T::value_type::LpCol,void>::type
  1333     colLowerBound(T &t, Value value,dummy<0> = 0) {
  1334       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1335         colLowerBound(*i, value);
  1336       }
  1337     }
  1338     template<class T>
  1339     typename enable_if<typename T::value_type::second_type::LpCol,
  1340                        void>::type
  1341     colLowerBound(T &t, Value value,dummy<1> = 1) {
  1342       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1343         colLowerBound(i->second, value);
  1344       }
  1345     }
  1346     template<class T>
  1347     typename enable_if<typename T::MapIt::Value::LpCol,
  1348                        void>::type
  1349     colLowerBound(T &t, Value value,dummy<2> = 2) {
  1350       for(typename T::MapIt i(t); i!=INVALID; ++i){
  1351         colLowerBound(*i, value);
  1352       }
  1353     }
  1354 #endif
  1355 
  1356     /// Set the upper bound of a column (i.e a variable)
  1357 
  1358     /// The upper bound of a variable (column) has to be given by an
  1359     /// extended number of type Value, i.e. a finite number of type
  1360     /// Value or \ref INF.
  1361     void colUpperBound(Col c, Value value) {
  1362       _setColUpperBound(cols(id(c)),value);
  1363     };
  1364 
  1365     /// Get the upper bound of a column (i.e a variable)
  1366 
  1367     /// This function returns the upper bound for column (variable) \c c
  1368     /// (this might be \ref INF as well).
  1369     /// \return The upper bound for column \c c
  1370     Value colUpperBound(Col c) const {
  1371       return _getColUpperBound(cols(id(c)));
  1372     }
  1373 
  1374     ///\brief Set the upper bound of  several columns
  1375     ///(i.e variables) at once
  1376     ///
  1377     ///This magic function takes a container as its argument
  1378     ///and applies the function on all of its elements.
  1379     ///The upper bound of a variable (column) has to be given by an
  1380     ///extended number of type Value, i.e. a finite number of type
  1381     ///Value or \ref INF.
  1382 #ifdef DOXYGEN
  1383     template<class T>
  1384     void colUpperBound(T &t, Value value) { return 0;}
  1385 #else
  1386     template<class T>
  1387     typename enable_if<typename T::value_type::LpCol,void>::type
  1388     colUpperBound(T &t, Value value,dummy<0> = 0) {
  1389       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1390         colUpperBound(*i, value);
  1391       }
  1392     }
  1393     template<class T>
  1394     typename enable_if<typename T::value_type::second_type::LpCol,
  1395                        void>::type
  1396     colUpperBound(T &t, Value value,dummy<1> = 1) {
  1397       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1398         colUpperBound(i->second, value);
  1399       }
  1400     }
  1401     template<class T>
  1402     typename enable_if<typename T::MapIt::Value::LpCol,
  1403                        void>::type
  1404     colUpperBound(T &t, Value value,dummy<2> = 2) {
  1405       for(typename T::MapIt i(t); i!=INVALID; ++i){
  1406         colUpperBound(*i, value);
  1407       }
  1408     }
  1409 #endif
  1410 
  1411     /// Set the lower and the upper bounds of a column (i.e a variable)
  1412 
  1413     /// The lower and the upper bounds of
  1414     /// a variable (column) have to be given by an
  1415     /// extended number of type Value, i.e. a finite number of type
  1416     /// Value, -\ref INF or \ref INF.
  1417     void colBounds(Col c, Value lower, Value upper) {
  1418       _setColLowerBound(cols(id(c)),lower);
  1419       _setColUpperBound(cols(id(c)),upper);
  1420     }
  1421 
  1422     ///\brief Set the lower and the upper bound of several columns
  1423     ///(i.e variables) at once
  1424     ///
  1425     ///This magic function takes a container as its argument
  1426     ///and applies the function on all of its elements.
  1427     /// The lower and the upper bounds of
  1428     /// a variable (column) have to be given by an
  1429     /// extended number of type Value, i.e. a finite number of type
  1430     /// Value, -\ref INF or \ref INF.
  1431 #ifdef DOXYGEN
  1432     template<class T>
  1433     void colBounds(T &t, Value lower, Value upper) { return 0;}
  1434 #else
  1435     template<class T>
  1436     typename enable_if<typename T::value_type::LpCol,void>::type
  1437     colBounds(T &t, Value lower, Value upper,dummy<0> = 0) {
  1438       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1439         colBounds(*i, lower, upper);
  1440       }
  1441     }
  1442     template<class T>
  1443     typename enable_if<typename T::value_type::second_type::LpCol, void>::type
  1444     colBounds(T &t, Value lower, Value upper,dummy<1> = 1) {
  1445       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1446         colBounds(i->second, lower, upper);
  1447       }
  1448     }
  1449     template<class T>
  1450     typename enable_if<typename T::MapIt::Value::LpCol, void>::type
  1451     colBounds(T &t, Value lower, Value upper,dummy<2> = 2) {
  1452       for(typename T::MapIt i(t); i!=INVALID; ++i){
  1453         colBounds(*i, lower, upper);
  1454       }
  1455     }
  1456 #endif
  1457 
  1458     /// Set the lower bound of a row (i.e a constraint)
  1459 
  1460     /// The lower bound of a constraint (row) has to be given by an
  1461     /// extended number of type Value, i.e. a finite number of type
  1462     /// Value or -\ref INF.
  1463     void rowLowerBound(Row r, Value value) {
  1464       _setRowLowerBound(rows(id(r)),value);
  1465     }
  1466 
  1467     /// Get the lower bound of a row (i.e a constraint)
  1468 
  1469     /// This function returns the lower bound for row (constraint) \c c
  1470     /// (this might be -\ref INF as well).
  1471     ///\return The lower bound for row \c r
  1472     Value rowLowerBound(Row r) const {
  1473       return _getRowLowerBound(rows(id(r)));
  1474     }
  1475 
  1476     /// Set the upper bound of a row (i.e a constraint)
  1477 
  1478     /// The upper bound of a constraint (row) has to be given by an
  1479     /// extended number of type Value, i.e. a finite number of type
  1480     /// Value or -\ref INF.
  1481     void rowUpperBound(Row r, Value value) {
  1482       _setRowUpperBound(rows(id(r)),value);
  1483     }
  1484 
  1485     /// Get the upper bound of a row (i.e a constraint)
  1486 
  1487     /// This function returns the upper bound for row (constraint) \c c
  1488     /// (this might be -\ref INF as well).
  1489     ///\return The upper bound for row \c r
  1490     Value rowUpperBound(Row r) const {
  1491       return _getRowUpperBound(rows(id(r)));
  1492     }
  1493 
  1494     ///Set an element of the objective function
  1495     void objCoeff(Col c, Value v) {_setObjCoeff(cols(id(c)),v); };
  1496 
  1497     ///Get an element of the objective function
  1498     Value objCoeff(Col c) const { return _getObjCoeff(cols(id(c))); };
  1499 
  1500     ///Set the objective function
  1501 
  1502     ///\param e is a linear expression of type \ref Expr.
  1503     ///
  1504     void obj(const Expr& e) {
  1505       _setObjCoeffs(ExprIterator(e.comps.begin(), cols),
  1506                     ExprIterator(e.comps.end(), cols));
  1507       obj_const_comp = *e;
  1508     }
  1509 
  1510     ///Get the objective function
  1511 
  1512     ///\return the objective function as a linear expression of type
  1513     ///Expr.
  1514     Expr obj() const {
  1515       Expr e;
  1516       _getObjCoeffs(InsertIterator(e.comps, cols));
  1517       *e = obj_const_comp;
  1518       return e;
  1519     }
  1520 
  1521 
  1522     ///Set the direction of optimization
  1523     void sense(Sense sense) { _setSense(sense); }
  1524 
  1525     ///Query the direction of the optimization
  1526     Sense sense() const {return _getSense(); }
  1527 
  1528     ///Set the sense to maximization
  1529     void max() { _setSense(MAX); }
  1530 
  1531     ///Set the sense to maximization
  1532     void min() { _setSense(MIN); }
  1533 
  1534     ///Clears the problem
  1535     void clear() { _clear(); }
  1536 
  1537     ///@}
  1538 
  1539   };
  1540 
  1541   /// Addition
  1542 
  1543   ///\relates LpBase::Expr
  1544   ///
  1545   inline LpBase::Expr operator+(const LpBase::Expr &a, const LpBase::Expr &b) {
  1546     LpBase::Expr tmp(a);
  1547     tmp+=b;
  1548     return tmp;
  1549   }
  1550   ///Substraction
  1551 
  1552   ///\relates LpBase::Expr
  1553   ///
  1554   inline LpBase::Expr operator-(const LpBase::Expr &a, const LpBase::Expr &b) {
  1555     LpBase::Expr tmp(a);
  1556     tmp-=b;
  1557     return tmp;
  1558   }
  1559   ///Multiply with constant
  1560 
  1561   ///\relates LpBase::Expr
  1562   ///
  1563   inline LpBase::Expr operator*(const LpBase::Expr &a, const LpBase::Value &b) {
  1564     LpBase::Expr tmp(a);
  1565     tmp*=b;
  1566     return tmp;
  1567   }
  1568 
  1569   ///Multiply with constant
  1570 
  1571   ///\relates LpBase::Expr
  1572   ///
  1573   inline LpBase::Expr operator*(const LpBase::Value &a, const LpBase::Expr &b) {
  1574     LpBase::Expr tmp(b);
  1575     tmp*=a;
  1576     return tmp;
  1577   }
  1578   ///Divide with constant
  1579 
  1580   ///\relates LpBase::Expr
  1581   ///
  1582   inline LpBase::Expr operator/(const LpBase::Expr &a, const LpBase::Value &b) {
  1583     LpBase::Expr tmp(a);
  1584     tmp/=b;
  1585     return tmp;
  1586   }
  1587 
  1588   ///Create constraint
  1589 
  1590   ///\relates LpBase::Constr
  1591   ///
  1592   inline LpBase::Constr operator<=(const LpBase::Expr &e,
  1593                                    const LpBase::Expr &f) {
  1594     return LpBase::Constr(0, f - e, LpBase::INF);
  1595   }
  1596 
  1597   ///Create constraint
  1598 
  1599   ///\relates LpBase::Constr
  1600   ///
  1601   inline LpBase::Constr operator<=(const LpBase::Value &e,
  1602                                    const LpBase::Expr &f) {
  1603     return LpBase::Constr(e, f, LpBase::NaN);
  1604   }
  1605 
  1606   ///Create constraint
  1607 
  1608   ///\relates LpBase::Constr
  1609   ///
  1610   inline LpBase::Constr operator<=(const LpBase::Expr &e,
  1611                                    const LpBase::Value &f) {
  1612     return LpBase::Constr(- LpBase::INF, e, f);
  1613   }
  1614 
  1615   ///Create constraint
  1616 
  1617   ///\relates LpBase::Constr
  1618   ///
  1619   inline LpBase::Constr operator>=(const LpBase::Expr &e,
  1620                                    const LpBase::Expr &f) {
  1621     return LpBase::Constr(0, e - f, LpBase::INF);
  1622   }
  1623 
  1624 
  1625   ///Create constraint
  1626 
  1627   ///\relates LpBase::Constr
  1628   ///
  1629   inline LpBase::Constr operator>=(const LpBase::Value &e,
  1630                                    const LpBase::Expr &f) {
  1631     return LpBase::Constr(LpBase::NaN, f, e);
  1632   }
  1633 
  1634 
  1635   ///Create constraint
  1636 
  1637   ///\relates LpBase::Constr
  1638   ///
  1639   inline LpBase::Constr operator>=(const LpBase::Expr &e,
  1640                                    const LpBase::Value &f) {
  1641     return LpBase::Constr(f, e, LpBase::INF);
  1642   }
  1643 
  1644   ///Create constraint
  1645 
  1646   ///\relates LpBase::Constr
  1647   ///
  1648   inline LpBase::Constr operator==(const LpBase::Expr &e,
  1649                                    const LpBase::Value &f) {
  1650     return LpBase::Constr(f, e, f);
  1651   }
  1652 
  1653   ///Create constraint
  1654 
  1655   ///\relates LpBase::Constr
  1656   ///
  1657   inline LpBase::Constr operator==(const LpBase::Expr &e,
  1658                                    const LpBase::Expr &f) {
  1659     return LpBase::Constr(0, f - e, 0);
  1660   }
  1661 
  1662   ///Create constraint
  1663 
  1664   ///\relates LpBase::Constr
  1665   ///
  1666   inline LpBase::Constr operator<=(const LpBase::Value &n,
  1667                                    const LpBase::Constr &c) {
  1668     LpBase::Constr tmp(c);
  1669     LEMON_ASSERT(isnan(tmp.lowerBound()), "Wrong LP constraint");
  1670     tmp.lowerBound()=n;
  1671     return tmp;
  1672   }
  1673   ///Create constraint
  1674 
  1675   ///\relates LpBase::Constr
  1676   ///
  1677   inline LpBase::Constr operator<=(const LpBase::Constr &c,
  1678                                    const LpBase::Value &n)
  1679   {
  1680     LpBase::Constr tmp(c);
  1681     LEMON_ASSERT(isnan(tmp.upperBound()), "Wrong LP constraint");
  1682     tmp.upperBound()=n;
  1683     return tmp;
  1684   }
  1685 
  1686   ///Create constraint
  1687 
  1688   ///\relates LpBase::Constr
  1689   ///
  1690   inline LpBase::Constr operator>=(const LpBase::Value &n,
  1691                                    const LpBase::Constr &c) {
  1692     LpBase::Constr tmp(c);
  1693     LEMON_ASSERT(isnan(tmp.upperBound()), "Wrong LP constraint");
  1694     tmp.upperBound()=n;
  1695     return tmp;
  1696   }
  1697   ///Create constraint
  1698 
  1699   ///\relates LpBase::Constr
  1700   ///
  1701   inline LpBase::Constr operator>=(const LpBase::Constr &c,
  1702                                    const LpBase::Value &n)
  1703   {
  1704     LpBase::Constr tmp(c);
  1705     LEMON_ASSERT(isnan(tmp.lowerBound()), "Wrong LP constraint");
  1706     tmp.lowerBound()=n;
  1707     return tmp;
  1708   }
  1709 
  1710   ///Addition
  1711 
  1712   ///\relates LpBase::DualExpr
  1713   ///
  1714   inline LpBase::DualExpr operator+(const LpBase::DualExpr &a,
  1715                                     const LpBase::DualExpr &b) {
  1716     LpBase::DualExpr tmp(a);
  1717     tmp+=b;
  1718     return tmp;
  1719   }
  1720   ///Substraction
  1721 
  1722   ///\relates LpBase::DualExpr
  1723   ///
  1724   inline LpBase::DualExpr operator-(const LpBase::DualExpr &a,
  1725                                     const LpBase::DualExpr &b) {
  1726     LpBase::DualExpr tmp(a);
  1727     tmp-=b;
  1728     return tmp;
  1729   }
  1730   ///Multiply with constant
  1731 
  1732   ///\relates LpBase::DualExpr
  1733   ///
  1734   inline LpBase::DualExpr operator*(const LpBase::DualExpr &a,
  1735                                     const LpBase::Value &b) {
  1736     LpBase::DualExpr tmp(a);
  1737     tmp*=b;
  1738     return tmp;
  1739   }
  1740 
  1741   ///Multiply with constant
  1742 
  1743   ///\relates LpBase::DualExpr
  1744   ///
  1745   inline LpBase::DualExpr operator*(const LpBase::Value &a,
  1746                                     const LpBase::DualExpr &b) {
  1747     LpBase::DualExpr tmp(b);
  1748     tmp*=a;
  1749     return tmp;
  1750   }
  1751   ///Divide with constant
  1752 
  1753   ///\relates LpBase::DualExpr
  1754   ///
  1755   inline LpBase::DualExpr operator/(const LpBase::DualExpr &a,
  1756                                     const LpBase::Value &b) {
  1757     LpBase::DualExpr tmp(a);
  1758     tmp/=b;
  1759     return tmp;
  1760   }
  1761 
  1762   /// \ingroup lp_group
  1763   ///
  1764   /// \brief Common base class for LP solvers
  1765   ///
  1766   /// This class is an abstract base class for LP solvers. This class
  1767   /// provides a full interface for set and modify an LP problem,
  1768   /// solve it and retrieve the solution. You can use one of the
  1769   /// descendants as a concrete implementation, or the \c Lp
  1770   /// default LP solver. However, if you would like to handle LP
  1771   /// solvers as reference or pointer in a generic way, you can use
  1772   /// this class directly.
  1773   class LpSolver : virtual public LpBase {
  1774   public:
  1775 
  1776     /// The problem types for primal and dual problems
  1777     enum ProblemType {
  1778       ///Feasible solution hasn't been found (but may exist).
  1779       UNDEFINED = 0,
  1780       ///The problem has no feasible solution
  1781       INFEASIBLE = 1,
  1782       ///Feasible solution found
  1783       FEASIBLE = 2,
  1784       ///Optimal solution exists and found
  1785       OPTIMAL = 3,
  1786       ///The cost function is unbounded
  1787       UNBOUNDED = 4
  1788     };
  1789 
  1790     ///The basis status of variables
  1791     enum VarStatus {
  1792       /// The variable is in the basis
  1793       BASIC, 
  1794       /// The variable is free, but not basic
  1795       FREE,
  1796       /// The variable has active lower bound 
  1797       LOWER,
  1798       /// The variable has active upper bound
  1799       UPPER,
  1800       /// The variable is non-basic and fixed
  1801       FIXED
  1802     };
  1803 
  1804   protected:
  1805 
  1806     virtual SolveExitStatus _solve() = 0;
  1807 
  1808     virtual Value _getPrimal(int i) const = 0;
  1809     virtual Value _getDual(int i) const = 0;
  1810 
  1811     virtual Value _getPrimalRay(int i) const = 0;
  1812     virtual Value _getDualRay(int i) const = 0;
  1813 
  1814     virtual Value _getPrimalValue() const = 0;
  1815 
  1816     virtual VarStatus _getColStatus(int i) const = 0;
  1817     virtual VarStatus _getRowStatus(int i) const = 0;
  1818 
  1819     virtual ProblemType _getPrimalType() const = 0;
  1820     virtual ProblemType _getDualType() const = 0;
  1821 
  1822   public:
  1823 
  1824     ///\name Solve the LP
  1825 
  1826     ///@{
  1827 
  1828     ///\e Solve the LP problem at hand
  1829     ///
  1830     ///\return The result of the optimization procedure. Possible
  1831     ///values and their meanings can be found in the documentation of
  1832     ///\ref SolveExitStatus.
  1833     SolveExitStatus solve() { return _solve(); }
  1834 
  1835     ///@}
  1836 
  1837     ///\name Obtain the solution
  1838 
  1839     ///@{
  1840 
  1841     /// The type of the primal problem
  1842     ProblemType primalType() const {
  1843       return _getPrimalType();
  1844     }
  1845 
  1846     /// The type of the dual problem
  1847     ProblemType dualType() const {
  1848       return _getDualType();
  1849     }
  1850 
  1851     /// Return the primal value of the column
  1852 
  1853     /// Return the primal value of the column.
  1854     /// \pre The problem is solved.
  1855     Value primal(Col c) const { return _getPrimal(cols(id(c))); }
  1856 
  1857     /// Return the primal value of the expression
  1858 
  1859     /// Return the primal value of the expression, i.e. the dot
  1860     /// product of the primal solution and the expression.
  1861     /// \pre The problem is solved.
  1862     Value primal(const Expr& e) const {
  1863       double res = *e;
  1864       for (Expr::ConstCoeffIt c(e); c != INVALID; ++c) {
  1865         res += *c * primal(c);
  1866       }
  1867       return res;
  1868     }
  1869     /// Returns a component of the primal ray
  1870     
  1871     /// The primal ray is solution of the modified primal problem,
  1872     /// where we change each finite bound to 0, and we looking for a
  1873     /// negative objective value in case of minimization, and positive
  1874     /// objective value for maximization. If there is such solution,
  1875     /// that proofs the unsolvability of the dual problem, and if a
  1876     /// feasible primal solution exists, then the unboundness of
  1877     /// primal problem.
  1878     ///
  1879     /// \pre The problem is solved and the dual problem is infeasible.
  1880     /// \note Some solvers does not provide primal ray calculation
  1881     /// functions.
  1882     Value primalRay(Col c) const { return _getPrimalRay(cols(id(c))); }
  1883 
  1884     /// Return the dual value of the row
  1885 
  1886     /// Return the dual value of the row.
  1887     /// \pre The problem is solved.
  1888     Value dual(Row r) const { return _getDual(rows(id(r))); }
  1889 
  1890     /// Return the dual value of the dual expression
  1891 
  1892     /// Return the dual value of the dual expression, i.e. the dot
  1893     /// product of the dual solution and the dual expression.
  1894     /// \pre The problem is solved.
  1895     Value dual(const DualExpr& e) const {
  1896       double res = 0.0;
  1897       for (DualExpr::ConstCoeffIt r(e); r != INVALID; ++r) {
  1898         res += *r * dual(r);
  1899       }
  1900       return res;
  1901     }
  1902 
  1903     /// Returns a component of the dual ray
  1904     
  1905     /// The dual ray is solution of the modified primal problem, where
  1906     /// we change each finite bound to 0 (i.e. the objective function
  1907     /// coefficients in the primal problem), and we looking for a
  1908     /// ositive objective value. If there is such solution, that
  1909     /// proofs the unsolvability of the primal problem, and if a
  1910     /// feasible dual solution exists, then the unboundness of
  1911     /// dual problem.
  1912     ///
  1913     /// \pre The problem is solved and the primal problem is infeasible.
  1914     /// \note Some solvers does not provide dual ray calculation
  1915     /// functions.
  1916     Value dualRay(Row r) const { return _getDualRay(rows(id(r))); }
  1917 
  1918     /// Return the basis status of the column
  1919 
  1920     /// \see VarStatus
  1921     VarStatus colStatus(Col c) const { return _getColStatus(cols(id(c))); }
  1922 
  1923     /// Return the basis status of the row
  1924 
  1925     /// \see VarStatus
  1926     VarStatus rowStatus(Row r) const { return _getRowStatus(rows(id(r))); }
  1927 
  1928     ///The value of the objective function
  1929 
  1930     ///\return
  1931     ///- \ref INF or -\ref INF means either infeasibility or unboundedness
  1932     /// of the primal problem, depending on whether we minimize or maximize.
  1933     ///- \ref NaN if no primal solution is found.
  1934     ///- The (finite) objective value if an optimal solution is found.
  1935     Value primal() const { return _getPrimalValue()+obj_const_comp;}
  1936     ///@}
  1937 
  1938     LpSolver* newSolver() {return _newSolver();}
  1939     LpSolver* cloneSolver() {return _cloneSolver();}
  1940 
  1941   protected:
  1942 
  1943     virtual LpSolver* _newSolver() const = 0;
  1944     virtual LpSolver* _cloneSolver() const = 0;
  1945   };
  1946 
  1947 
  1948   /// \ingroup lp_group
  1949   ///
  1950   /// \brief Common base class for MIP solvers
  1951   ///
  1952   /// This class is an abstract base class for MIP solvers. This class
  1953   /// provides a full interface for set and modify an MIP problem,
  1954   /// solve it and retrieve the solution. You can use one of the
  1955   /// descendants as a concrete implementation, or the \c Lp
  1956   /// default MIP solver. However, if you would like to handle MIP
  1957   /// solvers as reference or pointer in a generic way, you can use
  1958   /// this class directly.
  1959   class MipSolver : virtual public LpBase {
  1960   public:
  1961 
  1962     /// The problem types for MIP problems
  1963     enum ProblemType {
  1964       ///Feasible solution hasn't been found (but may exist).
  1965       UNDEFINED = 0,
  1966       ///The problem has no feasible solution
  1967       INFEASIBLE = 1,
  1968       ///Feasible solution found
  1969       FEASIBLE = 2,
  1970       ///Optimal solution exists and found
  1971       OPTIMAL = 3,
  1972       ///The cost function is unbounded
  1973       ///
  1974       ///The Mip or at least the relaxed problem is unbounded
  1975       UNBOUNDED = 4
  1976     };
  1977 
  1978     ///\name Solve the MIP
  1979 
  1980     ///@{
  1981 
  1982     /// Solve the MIP problem at hand
  1983     ///
  1984     ///\return The result of the optimization procedure. Possible
  1985     ///values and their meanings can be found in the documentation of
  1986     ///\ref SolveExitStatus.
  1987     SolveExitStatus solve() { return _solve(); }
  1988 
  1989     ///@}
  1990 
  1991     ///\name Setting column type
  1992     ///@{
  1993 
  1994     ///Possible variable (column) types (e.g. real, integer, binary etc.)
  1995     enum ColTypes {
  1996       ///Continuous variable (default)
  1997       REAL = 0,
  1998       ///Integer variable
  1999       INTEGER = 1
  2000     };
  2001 
  2002     ///Sets the type of the given column to the given type
  2003 
  2004     ///Sets the type of the given column to the given type.
  2005     ///
  2006     void colType(Col c, ColTypes col_type) {
  2007       _setColType(cols(id(c)),col_type);
  2008     }
  2009 
  2010     ///Gives back the type of the column.
  2011 
  2012     ///Gives back the type of the column.
  2013     ///
  2014     ColTypes colType(Col c) const {
  2015       return _getColType(cols(id(c)));
  2016     }
  2017     ///@}
  2018 
  2019     ///\name Obtain the solution
  2020 
  2021     ///@{
  2022 
  2023     /// The type of the MIP problem
  2024     ProblemType type() const {
  2025       return _getType();
  2026     }
  2027 
  2028     /// Return the value of the row in the solution
  2029 
  2030     ///  Return the value of the row in the solution.
  2031     /// \pre The problem is solved.
  2032     Value sol(Col c) const { return _getSol(cols(id(c))); }
  2033 
  2034     /// Return the value of the expression in the solution
  2035 
  2036     /// Return the value of the expression in the solution, i.e. the
  2037     /// dot product of the solution and the expression.
  2038     /// \pre The problem is solved.
  2039     Value sol(const Expr& e) const {
  2040       double res = *e;
  2041       for (Expr::ConstCoeffIt c(e); c != INVALID; ++c) {
  2042         res += *c * sol(c);
  2043       }
  2044       return res;
  2045     }
  2046     ///The value of the objective function
  2047     
  2048     ///\return
  2049     ///- \ref INF or -\ref INF means either infeasibility or unboundedness
  2050     /// of the problem, depending on whether we minimize or maximize.
  2051     ///- \ref NaN if no primal solution is found.
  2052     ///- The (finite) objective value if an optimal solution is found.
  2053     Value solValue() const { return _getSolValue()+obj_const_comp;}
  2054     ///@}
  2055 
  2056   protected:
  2057 
  2058     virtual SolveExitStatus _solve() = 0;
  2059     virtual ColTypes _getColType(int col) const = 0;
  2060     virtual void _setColType(int col, ColTypes col_type) = 0;
  2061     virtual ProblemType _getType() const = 0;
  2062     virtual Value _getSol(int i) const = 0;
  2063     virtual Value _getSolValue() const = 0;
  2064 
  2065   public:
  2066 
  2067     MipSolver* newSolver() {return _newSolver();}
  2068     MipSolver* cloneSolver() {return _cloneSolver();}
  2069 
  2070   protected:
  2071 
  2072     virtual MipSolver* _newSolver() const = 0;
  2073     virtual MipSolver* _cloneSolver() const = 0;
  2074   };
  2075 
  2076 
  2077 
  2078 } //namespace lemon
  2079 
  2080 #endif //LEMON_LP_BASE_H