lemon/lp_base.h
author deba
Mon, 27 Feb 2006 10:36:01 +0000
changeset 1986 9b56cca61e2e
parent 1908 e225719bde6b
child 1993 2115143eceea
permissions -rw-r--r--
An additional simplier interface for static size graphs.
Node operator()(int) for getting node by index
int index(Node node) for getting index by node
     1 /* -*- C++ -*-
     2  *
     3  * This file is a part of LEMON, a generic C++ optimization library
     4  *
     5  * Copyright (C) 2003-2006
     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<vector>
    23 #include<map>
    24 #include<limits>
    25 #include<cmath>
    26 
    27 #include<lemon/utility.h>
    28 #include<lemon/error.h>
    29 #include<lemon/invalid.h>
    30 
    31 ///\file
    32 ///\brief The interface of the LP solver interface.
    33 ///\ingroup gen_opt_group
    34 namespace lemon {
    35   
    36   ///Internal data structure to convert floating id's to fix one's
    37     
    38   ///\todo This might be implemented to be also usable in other places.
    39   class _FixId 
    40   {
    41   protected:
    42     std::vector<int> index;
    43     std::vector<int> cross;
    44     int first_free;
    45   public:
    46     _FixId() : first_free(-1) {};
    47     ///Convert a floating id to a fix one
    48 
    49     ///\param n is a floating id
    50     ///\return the corresponding fix id
    51     int fixId(int n) const {return cross[n];}
    52     ///Convert a fix id to a floating one
    53 
    54     ///\param n is a fix id
    55     ///\return the corresponding floating id
    56     int floatingId(int n) const { return index[n];}
    57     ///Add a new floating id.
    58 
    59     ///\param n is a floating id
    60     ///\return the fix id of the new value
    61     ///\todo Multiple additions should also be handled.
    62     int insert(int n)
    63     {
    64       if(n>=int(cross.size())) {
    65 	cross.resize(n+1);
    66 	if(first_free==-1) {
    67 	  cross[n]=index.size();
    68 	  index.push_back(n);
    69 	}
    70 	else {
    71 	  cross[n]=first_free;
    72 	  int next=index[first_free];
    73 	  index[first_free]=n;
    74 	  first_free=next;
    75 	}
    76 	return cross[n];
    77       }
    78       ///\todo Create an own exception type.
    79       else throw LogicError(); //floatingId-s must form a continuous range;
    80     }
    81     ///Remove a fix id.
    82 
    83     ///\param n is a fix id
    84     ///
    85     void erase(int n) 
    86     {
    87       int fl=index[n];
    88       index[n]=first_free;
    89       first_free=n;
    90       for(int i=fl+1;i<int(cross.size());++i) {
    91 	cross[i-1]=cross[i];
    92 	index[cross[i]]--;
    93       }
    94       cross.pop_back();
    95     }
    96     ///An upper bound on the largest fix id.
    97 
    98     ///\todo Do we need this?
    99     ///
   100     std::size_t maxFixId() { return cross.size()-1; }
   101   
   102   };
   103     
   104   ///Common base class for LP solvers
   105   
   106   ///\todo Much more docs
   107   ///\ingroup gen_opt_group
   108   class LpSolverBase {
   109 
   110   public:
   111 
   112     ///Possible outcomes of an LP solving procedure
   113     enum SolveExitStatus {
   114       ///This means that the problem has been successfully solved: either
   115       ///an optimal solution has been found or infeasibility/unboundedness
   116       ///has been proved.
   117       SOLVED = 0,
   118       ///Any other case (including the case when some user specified limit has been exceeded)
   119       UNSOLVED = 1
   120     };
   121       
   122       ///\e
   123     enum SolutionStatus {
   124       ///Feasible solution has'n been found (but may exist).
   125 
   126       ///\todo NOTFOUND might be a better name.
   127       ///
   128       UNDEFINED = 0,
   129       ///The problem has no feasible solution
   130       INFEASIBLE = 1,
   131       ///Feasible solution found
   132       FEASIBLE = 2,
   133       ///Optimal solution exists and found
   134       OPTIMAL = 3,
   135       ///The cost function is unbounded
   136 
   137       ///\todo Give a feasible solution and an infinite ray (and the
   138       ///corresponding bases)
   139       INFINITE = 4
   140     };
   141 
   142     ///\e The type of the investigated LP problem
   143     enum ProblemTypes {
   144       ///Primal-dual feasible
   145       PRIMAL_DUAL_FEASIBLE = 0,
   146       ///Primal feasible dual infeasible
   147       PRIMAL_FEASIBLE_DUAL_INFEASIBLE = 1,
   148       ///Primal infeasible dual feasible
   149       PRIMAL_INFEASIBLE_DUAL_FEASIBLE = 2,
   150       ///Primal-dual infeasible
   151       PRIMAL_DUAL_INFEASIBLE = 3,
   152       ///Could not determine so far
   153       UNKNOWN = 4
   154     };
   155 
   156     ///The floating point type used by the solver
   157     typedef double Value;
   158     ///The infinity constant
   159     static const Value INF;
   160     ///The not a number constant
   161     static const Value NaN;
   162     
   163     ///Refer to a column of the LP.
   164 
   165     ///This type is used to refer to a column of the LP.
   166     ///
   167     ///Its value remains valid and correct even after the addition or erase of
   168     ///other columns.
   169     ///
   170     ///\todo Document what can one do with a Col (INVALID, comparing,
   171     ///it is similar to Node/Edge)
   172     class Col {
   173     protected:
   174       int id;
   175       friend class LpSolverBase;
   176     public:
   177       typedef Value ExprValue;
   178       typedef True LpSolverCol;
   179       Col() {}
   180       Col(const Invalid&) : id(-1) {}
   181       bool operator< (Col c) const  {return id< c.id;}
   182       bool operator> (Col c) const  {return id> c.id;}
   183       bool operator==(Col c) const  {return id==c.id;}
   184       bool operator!=(Col c) const  {return id!=c.id;}
   185     };
   186 
   187     ///Refer to a row of the LP.
   188 
   189     ///This type is used to refer to a row of the LP.
   190     ///
   191     ///Its value remains valid and correct even after the addition or erase of
   192     ///other rows.
   193     ///
   194     ///\todo Document what can one do with a Row (INVALID, comparing,
   195     ///it is similar to Node/Edge)
   196     class Row {
   197     protected:
   198       int id;
   199       friend class LpSolverBase;
   200     public:
   201       typedef Value ExprValue;
   202       typedef True LpSolverRow;
   203       Row() {}
   204       Row(const Invalid&) : id(-1) {}
   205 
   206       bool operator< (Row c) const  {return id< c.id;}
   207       bool operator> (Row c) const  {return id> c.id;}
   208       bool operator==(Row c) const  {return id==c.id;}
   209       bool operator!=(Row c) const  {return id!=c.id;} 
   210    };
   211     
   212     ///Linear expression of variables and a constant component
   213     
   214     ///This data structure strores a linear expression of the variables
   215     ///(\ref Col "Col"s) and also has a constant component.
   216     ///
   217     ///There are several ways to access and modify the contents of this
   218     ///container.
   219     ///- Its it fully compatible with \c std::map<Col,double>, so for expamle
   220     ///if \c e is an Expr and \c v and \c w are of type \ref Col, then you can
   221     ///read and modify the coefficients like
   222     ///these.
   223     ///\code
   224     ///e[v]=5;
   225     ///e[v]+=12;
   226     ///e.erase(v);
   227     ///\endcode
   228     ///or you can also iterate through its elements.
   229     ///\code
   230     ///double s=0;
   231     ///for(LpSolverBase::Expr::iterator i=e.begin();i!=e.end();++i)
   232     ///  s+=i->second;
   233     ///\endcode
   234     ///(This code computes the sum of all coefficients).
   235     ///- Numbers (<tt>double</tt>'s)
   236     ///and variables (\ref Col "Col"s) directly convert to an
   237     ///\ref Expr and the usual linear operations are defined, so  
   238     ///\code
   239     ///v+w
   240     ///2*v-3.12*(v-w/2)+2
   241     ///v*2.1+(3*v+(v*12+w+6)*3)/2
   242     ///\endcode
   243     ///are valid \ref Expr "Expr"essions.
   244     ///The usual assignment operations are also defined.
   245     ///\code
   246     ///e=v+w;
   247     ///e+=2*v-3.12*(v-w/2)+2;
   248     ///e*=3.4;
   249     ///e/=5;
   250     ///\endcode
   251     ///- The constant member can be set and read by \ref constComp()
   252     ///\code
   253     ///e.constComp()=12;
   254     ///double c=e.constComp();
   255     ///\endcode
   256     ///
   257     ///\note \ref clear() not only sets all coefficients to 0 but also
   258     ///clears the constant components.
   259     ///
   260     ///\sa Constr
   261     ///
   262     class Expr : public std::map<Col,Value>
   263     {
   264     public:
   265       typedef LpSolverBase::Col Key; 
   266       typedef LpSolverBase::Value Value;
   267       
   268     protected:
   269       typedef std::map<Col,Value> Base;
   270       
   271       Value const_comp;
   272   public:
   273       typedef True IsLinExpression;
   274       ///\e
   275       Expr() : Base(), const_comp(0) { }
   276       ///\e
   277       Expr(const Key &v) : const_comp(0) {
   278 	Base::insert(std::make_pair(v, 1));
   279       }
   280       ///\e
   281       Expr(const Value &v) : const_comp(v) {}
   282       ///\e
   283       void set(const Key &v,const Value &c) {
   284 	Base::insert(std::make_pair(v, c));
   285       }
   286       ///\e
   287       Value &constComp() { return const_comp; }
   288       ///\e
   289       const Value &constComp() const { return const_comp; }
   290       
   291       ///Removes the components with zero coefficient.
   292       void simplify() {
   293 	for (Base::iterator i=Base::begin(); i!=Base::end();) {
   294 	  Base::iterator j=i;
   295 	  ++j;
   296 	  if ((*i).second==0) Base::erase(i);
   297 	  j=i;
   298 	}
   299       }
   300 
   301       ///Removes the coefficients closer to zero than \c tolerance.
   302       void simplify(double &tolerance) {
   303 	for (Base::iterator i=Base::begin(); i!=Base::end();) {
   304 	  Base::iterator j=i;
   305 	  ++j;
   306 	  if (std::fabs((*i).second)<tolerance) Base::erase(i);
   307 	  j=i;
   308 	}
   309       }
   310 
   311       ///Sets all coefficients and the constant component to 0.
   312       void clear() {
   313 	Base::clear();
   314 	const_comp=0;
   315       }
   316 
   317       ///\e
   318       Expr &operator+=(const Expr &e) {
   319 	for (Base::const_iterator j=e.begin(); j!=e.end(); ++j)
   320 	  (*this)[j->first]+=j->second;
   321 	const_comp+=e.const_comp;
   322 	return *this;
   323       }
   324       ///\e
   325       Expr &operator-=(const Expr &e) {
   326 	for (Base::const_iterator j=e.begin(); j!=e.end(); ++j)
   327 	  (*this)[j->first]-=j->second;
   328 	const_comp-=e.const_comp;
   329 	return *this;
   330       }
   331       ///\e
   332       Expr &operator*=(const Value &c) {
   333 	for (Base::iterator j=Base::begin(); j!=Base::end(); ++j)
   334 	  j->second*=c;
   335 	const_comp*=c;
   336 	return *this;
   337       }
   338       ///\e
   339       Expr &operator/=(const Value &c) {
   340 	for (Base::iterator j=Base::begin(); j!=Base::end(); ++j)
   341 	  j->second/=c;
   342 	const_comp/=c;
   343 	return *this;
   344       }
   345     };
   346     
   347     ///Linear constraint
   348 
   349     ///This data stucture represents a linear constraint in the LP.
   350     ///Basically it is a linear expression with a lower or an upper bound
   351     ///(or both). These parts of the constraint can be obtained by the member
   352     ///functions \ref expr(), \ref lowerBound() and \ref upperBound(),
   353     ///respectively.
   354     ///There are two ways to construct a constraint.
   355     ///- You can set the linear expression and the bounds directly
   356     ///  by the functions above.
   357     ///- The operators <tt>\<=</tt>, <tt>==</tt> and  <tt>\>=</tt>
   358     ///  are defined between expressions, or even between constraints whenever
   359     ///  it makes sense. Therefore if \c e and \c f are linear expressions and
   360     ///  \c s and \c t are numbers, then the followings are valid expressions
   361     ///  and thus they can be used directly e.g. in \ref addRow() whenever
   362     ///  it makes sense.
   363     ///\code
   364     ///  e<=s
   365     ///  e<=f
   366     ///  e==f
   367     ///  s<=e<=t
   368     ///  e>=t
   369     ///\endcode
   370     ///\warning The validity of a constraint is checked only at run time, so
   371     ///e.g. \ref addRow(<tt>x[1]\<=x[2]<=5</tt>) will compile, but will throw a
   372     ///\ref LogicError exception.
   373     class Constr
   374     {
   375     public:
   376       typedef LpSolverBase::Expr Expr;
   377       typedef Expr::Key Key;
   378       typedef Expr::Value Value;
   379       
   380 //       static const Value INF;
   381 //       static const Value NaN;
   382 
   383     protected:
   384       Expr _expr;
   385       Value _lb,_ub;
   386     public:
   387       ///\e
   388       Constr() : _expr(), _lb(NaN), _ub(NaN) {}
   389       ///\e
   390       Constr(Value lb,const Expr &e,Value ub) :
   391 	_expr(e), _lb(lb), _ub(ub) {}
   392       ///\e
   393       Constr(const Expr &e,Value ub) : 
   394 	_expr(e), _lb(NaN), _ub(ub) {}
   395       ///\e
   396       Constr(Value lb,const Expr &e) :
   397 	_expr(e), _lb(lb), _ub(NaN) {}
   398       ///\e
   399       Constr(const Expr &e) : 
   400 	_expr(e), _lb(NaN), _ub(NaN) {}
   401       ///\e
   402       void clear() 
   403       {
   404 	_expr.clear();
   405 	_lb=_ub=NaN;
   406       }
   407 
   408       ///Reference to the linear expression 
   409       Expr &expr() { return _expr; }
   410       ///Cont reference to the linear expression 
   411       const Expr &expr() const { return _expr; }
   412       ///Reference to the lower bound.
   413 
   414       ///\return
   415       ///- \ref INF "INF": the constraint is lower unbounded.
   416       ///- \ref NaN "NaN": lower bound has not been set.
   417       ///- finite number: the lower bound
   418       Value &lowerBound() { return _lb; }
   419       ///The const version of \ref lowerBound()
   420       const Value &lowerBound() const { return _lb; }
   421       ///Reference to the upper bound.
   422 
   423       ///\return
   424       ///- \ref INF "INF": the constraint is upper unbounded.
   425       ///- \ref NaN "NaN": upper bound has not been set.
   426       ///- finite number: the upper bound
   427       Value &upperBound() { return _ub; }
   428       ///The const version of \ref upperBound()
   429       const Value &upperBound() const { return _ub; }
   430       ///Is the constraint lower bounded?
   431       bool lowerBounded() const { 
   432 	using namespace std;
   433 	return finite(_lb);
   434       }
   435       ///Is the constraint upper bounded?
   436       bool upperBounded() const {
   437 	using namespace std;
   438 	return finite(_ub);
   439       }
   440     };
   441     
   442     ///Linear expression of rows
   443     
   444     ///This data structure represents a column of the matrix,
   445     ///thas is it strores a linear expression of the dual variables
   446     ///(\ref Row "Row"s).
   447     ///
   448     ///There are several ways to access and modify the contents of this
   449     ///container.
   450     ///- Its it fully compatible with \c std::map<Row,double>, so for expamle
   451     ///if \c e is an DualExpr and \c v
   452     ///and \c w are of type \ref Row, then you can
   453     ///read and modify the coefficients like
   454     ///these.
   455     ///\code
   456     ///e[v]=5;
   457     ///e[v]+=12;
   458     ///e.erase(v);
   459     ///\endcode
   460     ///or you can also iterate through its elements.
   461     ///\code
   462     ///double s=0;
   463     ///for(LpSolverBase::DualExpr::iterator i=e.begin();i!=e.end();++i)
   464     ///  s+=i->second;
   465     ///\endcode
   466     ///(This code computes the sum of all coefficients).
   467     ///- Numbers (<tt>double</tt>'s)
   468     ///and variables (\ref Row "Row"s) directly convert to an
   469     ///\ref DualExpr and the usual linear operations are defined, so
   470     ///\code
   471     ///v+w
   472     ///2*v-3.12*(v-w/2)
   473     ///v*2.1+(3*v+(v*12+w)*3)/2
   474     ///\endcode
   475     ///are valid \ref DualExpr "DualExpr"essions.
   476     ///The usual assignment operations are also defined.
   477     ///\code
   478     ///e=v+w;
   479     ///e+=2*v-3.12*(v-w/2);
   480     ///e*=3.4;
   481     ///e/=5;
   482     ///\endcode
   483     ///
   484     ///\sa Expr
   485     ///
   486     class DualExpr : public std::map<Row,Value>
   487     {
   488     public:
   489       typedef LpSolverBase::Row Key; 
   490       typedef LpSolverBase::Value Value;
   491       
   492     protected:
   493       typedef std::map<Row,Value> Base;
   494       
   495     public:
   496       typedef True IsLinExpression;
   497       ///\e
   498       DualExpr() : Base() { }
   499       ///\e
   500       DualExpr(const Key &v) {
   501 	Base::insert(std::make_pair(v, 1));
   502       }
   503       ///\e
   504       void set(const Key &v,const Value &c) {
   505 	Base::insert(std::make_pair(v, c));
   506       }
   507       
   508       ///Removes the components with zero coefficient.
   509       void simplify() {
   510 	for (Base::iterator i=Base::begin(); i!=Base::end();) {
   511 	  Base::iterator j=i;
   512 	  ++j;
   513 	  if ((*i).second==0) Base::erase(i);
   514 	  j=i;
   515 	}
   516       }
   517 
   518       ///Removes the coefficients closer to zero than \c tolerance.
   519       void simplify(double &tolerance) {
   520 	for (Base::iterator i=Base::begin(); i!=Base::end();) {
   521 	  Base::iterator j=i;
   522 	  ++j;
   523 	  if (std::fabs((*i).second)<tolerance) Base::erase(i);
   524 	  j=i;
   525 	}
   526       }
   527 
   528 
   529       ///Sets all coefficients to 0.
   530       void clear() {
   531 	Base::clear();
   532       }
   533 
   534       ///\e
   535       DualExpr &operator+=(const DualExpr &e) {
   536 	for (Base::const_iterator j=e.begin(); j!=e.end(); ++j)
   537 	  (*this)[j->first]+=j->second;
   538 	return *this;
   539       }
   540       ///\e
   541       DualExpr &operator-=(const DualExpr &e) {
   542 	for (Base::const_iterator j=e.begin(); j!=e.end(); ++j)
   543 	  (*this)[j->first]-=j->second;
   544 	return *this;
   545       }
   546       ///\e
   547       DualExpr &operator*=(const Value &c) {
   548 	for (Base::iterator j=Base::begin(); j!=Base::end(); ++j)
   549 	  j->second*=c;
   550 	return *this;
   551       }
   552       ///\e
   553       DualExpr &operator/=(const Value &c) {
   554 	for (Base::iterator j=Base::begin(); j!=Base::end(); ++j)
   555 	  j->second/=c;
   556 	return *this;
   557       }
   558     };
   559     
   560 
   561   protected:
   562     _FixId rows;
   563     _FixId cols;
   564 
   565     //Abstract virtual functions
   566     virtual LpSolverBase &_newLp() = 0;
   567     virtual LpSolverBase &_copyLp(){
   568       ///\todo This should be implemented here, too,  when we have problem retrieving routines. It can be overriden.
   569 
   570       //Starting:
   571       LpSolverBase & newlp(_newLp());
   572       return newlp;
   573       //return *(LpSolverBase*)0;
   574     };
   575 
   576     virtual int _addCol() = 0;
   577     virtual int _addRow() = 0;
   578     virtual void _eraseCol(int col) = 0;
   579     virtual void _eraseRow(int row) = 0;
   580     virtual void _getColName(int col,       std::string & name) = 0;
   581     virtual void _setColName(int col, const std::string & name) = 0;
   582     virtual void _setRowCoeffs(int i, 
   583 			       int length,
   584                                int  const * indices, 
   585                                Value  const * values ) = 0;
   586     virtual void _setColCoeffs(int i, 
   587 			       int length,
   588                                int  const * indices, 
   589                                Value  const * values ) = 0;
   590     virtual void _setCoeff(int row, int col, Value value) = 0;
   591     virtual void _setColLowerBound(int i, Value value) = 0;
   592     virtual void _setColUpperBound(int i, Value value) = 0;
   593 //     virtual void _setRowLowerBound(int i, Value value) = 0;
   594 //     virtual void _setRowUpperBound(int i, Value value) = 0;
   595     virtual void _setRowBounds(int i, Value lower, Value upper) = 0;
   596     virtual void _setObjCoeff(int i, Value obj_coef) = 0;
   597     virtual void _clearObj()=0;
   598 //     virtual void _setObj(int length,
   599 //                          int  const * indices, 
   600 //                          Value  const * values ) = 0;
   601     virtual SolveExitStatus _solve() = 0;
   602     virtual Value _getPrimal(int i) = 0;
   603     virtual Value _getDual(int i) = 0;
   604     virtual Value _getPrimalValue() = 0;
   605     virtual bool _isBasicCol(int i) = 0;
   606     virtual SolutionStatus _getPrimalStatus() = 0;
   607     virtual SolutionStatus _getDualStatus() = 0;
   608     ///\todo This could be implemented here, too, using _getPrimalStatus() and
   609     ///_getDualStatus()
   610     virtual ProblemTypes _getProblemType() = 0;
   611 
   612     virtual void _setMax() = 0;
   613     virtual void _setMin() = 0;
   614     
   615     //Own protected stuff
   616     
   617     //Constant component of the objective function
   618     Value obj_const_comp;
   619     
   620 
   621 
   622     
   623   public:
   624 
   625     ///\e
   626     LpSolverBase() : obj_const_comp(0) {}
   627 
   628     ///\e
   629     virtual ~LpSolverBase() {}
   630 
   631     ///Creates a new LP problem
   632     LpSolverBase &newLp() {return _newLp();}
   633     ///Makes a copy of the LP problem
   634     LpSolverBase &copyLp() {return _copyLp();}
   635     
   636     ///\name Build up and modify the LP
   637 
   638     ///@{
   639 
   640     ///Add a new empty column (i.e a new variable) to the LP
   641     Col addCol() { Col c; c.id=cols.insert(_addCol()); return c;}
   642 
   643     ///\brief Adds several new columns
   644     ///(i.e a variables) at once
   645     ///
   646     ///This magic function takes a container as its argument
   647     ///and fills its elements
   648     ///with new columns (i.e. variables)
   649     ///\param t can be
   650     ///- a standard STL compatible iterable container with
   651     ///\ref Col as its \c values_type
   652     ///like
   653     ///\code
   654     ///std::vector<LpSolverBase::Col>
   655     ///std::list<LpSolverBase::Col>
   656     ///\endcode
   657     ///- a standard STL compatible iterable container with
   658     ///\ref Col as its \c mapped_type
   659     ///like
   660     ///\code
   661     ///std::map<AnyType,LpSolverBase::Col>
   662     ///\endcode
   663     ///- an iterable lemon \ref concept::WriteMap "write map" like 
   664     ///\code
   665     ///ListGraph::NodeMap<LpSolverBase::Col>
   666     ///ListGraph::EdgeMap<LpSolverBase::Col>
   667     ///\endcode
   668     ///\return The number of the created column.
   669 #ifdef DOXYGEN
   670     template<class T>
   671     int addColSet(T &t) { return 0;} 
   672 #else
   673     template<class T>
   674     typename enable_if<typename T::value_type::LpSolverCol,int>::type
   675     addColSet(T &t,dummy<0> = 0) {
   676       int s=0;
   677       for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addCol();s++;}
   678       return s;
   679     }
   680     template<class T>
   681     typename enable_if<typename T::value_type::second_type::LpSolverCol,
   682 		       int>::type
   683     addColSet(T &t,dummy<1> = 1) { 
   684       int s=0;
   685       for(typename T::iterator i=t.begin();i!=t.end();++i) {
   686 	i->second=addCol();
   687 	s++;
   688       }
   689       return s;
   690     }
   691     template<class T>
   692     typename enable_if<typename T::MapIt::Value::LpSolverCol,
   693 		       int>::type
   694     addColSet(T &t,dummy<2> = 2) { 
   695       int s=0;
   696       for(typename T::MapIt i(t); i!=INVALID; ++i)
   697 	{
   698 	  i.set(addCol());
   699 	  s++;
   700 	}
   701       return s;
   702     }
   703 #endif
   704 
   705     ///Set a column (i.e a dual constraint) of the LP
   706 
   707     ///\param c is the column to be modified
   708     ///\param e is a dual linear expression (see \ref DualExpr)
   709     ///a better one.
   710     void col(Col c,const DualExpr &e) {
   711       std::vector<int> indices;
   712       std::vector<Value> values;
   713       indices.push_back(0);
   714       values.push_back(0);
   715       for(DualExpr::const_iterator i=e.begin(); i!=e.end(); ++i)
   716 	if((*i).second!=0) {
   717 	  indices.push_back(rows.floatingId((*i).first.id));
   718 	  values.push_back((*i).second);
   719 	}
   720       _setColCoeffs(cols.floatingId(c.id),indices.size()-1,
   721 		    &indices[0],&values[0]);
   722     }
   723 
   724     ///Add a new column to the LP
   725 
   726     ///\param e is a dual linear expression (see \ref DualExpr)
   727     ///\param obj is the corresponding component of the objective
   728     ///function. It is 0 by default.
   729     ///\return The created column.
   730     Col addCol(const DualExpr &e, Value obj=0) {
   731       Col c=addCol();
   732       col(c,e);
   733       objCoeff(c,obj);
   734       return c;
   735     }
   736 
   737     ///Add a new empty row (i.e a new constraint) to the LP
   738 
   739     ///This function adds a new empty row (i.e a new constraint) to the LP.
   740     ///\return The created row
   741     Row addRow() { Row r; r.id=rows.insert(_addRow()); return r;}
   742 
   743     ///\brief Add several new rows
   744     ///(i.e a constraints) at once
   745     ///
   746     ///This magic function takes a container as its argument
   747     ///and fills its elements
   748     ///with new row (i.e. variables)
   749     ///\param t can be
   750     ///- a standard STL compatible iterable container with
   751     ///\ref Row as its \c values_type
   752     ///like
   753     ///\code
   754     ///std::vector<LpSolverBase::Row>
   755     ///std::list<LpSolverBase::Row>
   756     ///\endcode
   757     ///- a standard STL compatible iterable container with
   758     ///\ref Row as its \c mapped_type
   759     ///like
   760     ///\code
   761     ///std::map<AnyType,LpSolverBase::Row>
   762     ///\endcode
   763     ///- an iterable lemon \ref concept::WriteMap "write map" like 
   764     ///\code
   765     ///ListGraph::NodeMap<LpSolverBase::Row>
   766     ///ListGraph::EdgeMap<LpSolverBase::Row>
   767     ///\endcode
   768     ///\return The number of rows created.
   769 #ifdef DOXYGEN
   770     template<class T>
   771     int addRowSet(T &t) { return 0;} 
   772 #else
   773     template<class T>
   774     typename enable_if<typename T::value_type::LpSolverRow,int>::type
   775     addRowSet(T &t,dummy<0> = 0) {
   776       int s=0;
   777       for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addRow();s++;}
   778       return s;
   779     }
   780     template<class T>
   781     typename enable_if<typename T::value_type::second_type::LpSolverRow,
   782 		       int>::type
   783     addRowSet(T &t,dummy<1> = 1) { 
   784       int s=0;
   785       for(typename T::iterator i=t.begin();i!=t.end();++i) {
   786 	i->second=addRow();
   787 	s++;
   788       }
   789       return s;
   790     }
   791     template<class T>
   792     typename enable_if<typename T::MapIt::Value::LpSolverRow,
   793 		       int>::type
   794     addRowSet(T &t,dummy<2> = 2) { 
   795       int s=0;
   796       for(typename T::MapIt i(t); i!=INVALID; ++i)
   797 	{
   798 	  i.set(addRow());
   799 	  s++;
   800 	}
   801       return s;
   802     }
   803 #endif
   804 
   805     ///Set a row (i.e a constraint) of the LP
   806 
   807     ///\param r is the row to be modified
   808     ///\param l is lower bound (-\ref INF means no bound)
   809     ///\param e is a linear expression (see \ref Expr)
   810     ///\param u is the upper bound (\ref INF means no bound)
   811     ///\bug This is a temportary function. The interface will change to
   812     ///a better one.
   813     ///\todo Option to control whether a constraint with a single variable is
   814     ///added or not.
   815     void row(Row r, Value l,const Expr &e, Value u) {
   816       std::vector<int> indices;
   817       std::vector<Value> values;
   818       indices.push_back(0);
   819       values.push_back(0);
   820       for(Expr::const_iterator i=e.begin(); i!=e.end(); ++i)
   821 	if((*i).second!=0) { ///\bug EPSILON would be necessary here!!!
   822 	  indices.push_back(cols.floatingId((*i).first.id));
   823 	  values.push_back((*i).second);
   824 	}
   825       _setRowCoeffs(rows.floatingId(r.id),indices.size()-1,
   826 		    &indices[0],&values[0]);
   827 //       _setRowLowerBound(rows.floatingId(r.id),l-e.constComp());
   828 //       _setRowUpperBound(rows.floatingId(r.id),u-e.constComp());
   829        _setRowBounds(rows.floatingId(r.id),l-e.constComp(),u-e.constComp());
   830     }
   831 
   832     ///Set a row (i.e a constraint) of the LP
   833 
   834     ///\param r is the row to be modified
   835     ///\param c is a linear expression (see \ref Constr)
   836     void row(Row r, const Constr &c) {
   837       row(r,
   838 	     c.lowerBounded()?c.lowerBound():-INF,
   839 	     c.expr(),
   840 	     c.upperBounded()?c.upperBound():INF);
   841     }
   842 
   843     ///Add a new row (i.e a new constraint) to the LP
   844 
   845     ///\param l is the lower bound (-\ref INF means no bound)
   846     ///\param e is a linear expression (see \ref Expr)
   847     ///\param u is the upper bound (\ref INF means no bound)
   848     ///\return The created row.
   849     ///\bug This is a temportary function. The interface will change to
   850     ///a better one.
   851     Row addRow(Value l,const Expr &e, Value u) {
   852       Row r=addRow();
   853       row(r,l,e,u);
   854       return r;
   855     }
   856 
   857     ///Add a new row (i.e a new constraint) to the LP
   858 
   859     ///\param c is a linear expression (see \ref Constr)
   860     ///\return The created row.
   861     Row addRow(const Constr &c) {
   862       Row r=addRow();
   863       row(r,c);
   864       return r;
   865     }
   866     ///Erase a coloumn (i.e a variable) from the LP
   867 
   868     ///\param c is the coloumn to be deleted
   869     ///\todo Please check this
   870     void eraseCol(Col c) {
   871       _eraseCol(cols.floatingId(c.id));
   872       cols.erase(c.id);
   873     }
   874     ///Erase a  row (i.e a constraint) from the LP
   875 
   876     ///\param r is the row to be deleted
   877     ///\todo Please check this
   878     void eraseRow(Row r) {
   879       _eraseRow(rows.floatingId(r.id));
   880       rows.erase(r.id);
   881     }
   882 
   883     /// Get the name of a column
   884     
   885     ///\param c is the coresponding coloumn 
   886     ///\return The name of the colunm
   887     std::string ColName(Col c){
   888       std::string name;
   889       _getColName(cols.floatingId(c.id), name);
   890       return name;
   891     }
   892     
   893     /// Set the name of a column
   894     
   895     ///\param c is the coresponding coloumn 
   896     ///\param name The name to be given
   897     void ColName(Col c, const std::string & name){
   898       _setColName(cols.floatingId(c.id), name);
   899     }
   900     
   901     /// Set an element of the coefficient matrix of the LP
   902 
   903     ///\param r is the row of the element to be modified
   904     ///\param c is the coloumn of the element to be modified
   905     ///\param val is the new value of the coefficient
   906 
   907     void Coeff(Row r, Col c, Value val){
   908       _setCoeff(rows.floatingId(r.id),cols.floatingId(c.id), val);
   909     }
   910 
   911     /// Set the lower bound of a column (i.e a variable)
   912 
   913     /// The lower bound of a variable (column) has to be given by an 
   914     /// extended number of type Value, i.e. a finite number of type 
   915     /// Value or -\ref INF.
   916     void colLowerBound(Col c, Value value) {
   917       _setColLowerBound(cols.floatingId(c.id),value);
   918     }
   919     
   920     ///\brief Set the lower bound of  several columns
   921     ///(i.e a variables) at once
   922     ///
   923     ///This magic function takes a container as its argument
   924     ///and applies the function on all of its elements.
   925     /// The lower bound of a variable (column) has to be given by an 
   926     /// extended number of type Value, i.e. a finite number of type 
   927     /// Value or -\ref INF.
   928 #ifdef DOXYGEN
   929     template<class T>
   930     void colLowerBound(T &t, Value value) { return 0;} 
   931 #else
   932     template<class T>
   933     typename enable_if<typename T::value_type::LpSolverCol,void>::type
   934     colLowerBound(T &t, Value value,dummy<0> = 0) {
   935       for(typename T::iterator i=t.begin();i!=t.end();++i) {
   936 	colLowerBound(*i, value);
   937       }
   938     }
   939     template<class T>
   940     typename enable_if<typename T::value_type::second_type::LpSolverCol,
   941 		       void>::type
   942     colLowerBound(T &t, Value value,dummy<1> = 1) { 
   943       for(typename T::iterator i=t.begin();i!=t.end();++i) {
   944 	colLowerBound(i->second, value);
   945       }
   946     }
   947     template<class T>
   948     typename enable_if<typename T::MapIt::Value::LpSolverCol,
   949 		       void>::type
   950     colLowerBound(T &t, Value value,dummy<2> = 2) { 
   951       for(typename T::MapIt i(t); i!=INVALID; ++i){
   952 	colLowerBound(*i, value);
   953       }
   954     }
   955 #endif
   956     
   957     /// Set the upper bound of a column (i.e a variable)
   958 
   959     /// The upper bound of a variable (column) has to be given by an 
   960     /// extended number of type Value, i.e. a finite number of type 
   961     /// Value or \ref INF.
   962     void colUpperBound(Col c, Value value) {
   963       _setColUpperBound(cols.floatingId(c.id),value);
   964     };
   965 
   966     ///\brief Set the lower bound of  several columns
   967     ///(i.e a variables) at once
   968     ///
   969     ///This magic function takes a container as its argument
   970     ///and applies the function on all of its elements.
   971     /// The upper bound of a variable (column) has to be given by an 
   972     /// extended number of type Value, i.e. a finite number of type 
   973     /// Value or \ref INF.
   974 #ifdef DOXYGEN
   975     template<class T>
   976     void colUpperBound(T &t, Value value) { return 0;} 
   977 #else
   978     template<class T>
   979     typename enable_if<typename T::value_type::LpSolverCol,void>::type
   980     colUpperBound(T &t, Value value,dummy<0> = 0) {
   981       for(typename T::iterator i=t.begin();i!=t.end();++i) {
   982 	colUpperBound(*i, value);
   983       }
   984     }
   985     template<class T>
   986     typename enable_if<typename T::value_type::second_type::LpSolverCol,
   987 		       void>::type
   988     colUpperBound(T &t, Value value,dummy<1> = 1) { 
   989       for(typename T::iterator i=t.begin();i!=t.end();++i) {
   990 	colUpperBound(i->second, value);
   991       }
   992     }
   993     template<class T>
   994     typename enable_if<typename T::MapIt::Value::LpSolverCol,
   995 		       void>::type
   996     colUpperBound(T &t, Value value,dummy<2> = 2) { 
   997       for(typename T::MapIt i(t); i!=INVALID; ++i){
   998 	colUpperBound(*i, value);
   999       }
  1000     }
  1001 #endif
  1002 
  1003     /// Set the lower and the upper bounds of a column (i.e a variable)
  1004 
  1005     /// The lower and the upper bounds of
  1006     /// a variable (column) have to be given by an 
  1007     /// extended number of type Value, i.e. a finite number of type 
  1008     /// Value, -\ref INF or \ref INF.
  1009     void colBounds(Col c, Value lower, Value upper) {
  1010       _setColLowerBound(cols.floatingId(c.id),lower);
  1011       _setColUpperBound(cols.floatingId(c.id),upper);
  1012     }
  1013     
  1014     ///\brief Set the lower and the upper bound of several columns
  1015     ///(i.e a variables) at once
  1016     ///
  1017     ///This magic function takes a container as its argument
  1018     ///and applies the function on all of its elements.
  1019     /// The lower and the upper bounds of
  1020     /// a variable (column) have to be given by an 
  1021     /// extended number of type Value, i.e. a finite number of type 
  1022     /// Value, -\ref INF or \ref INF.
  1023 #ifdef DOXYGEN
  1024     template<class T>
  1025     void colBounds(T &t, Value lower, Value upper) { return 0;} 
  1026 #else
  1027     template<class T>
  1028     typename enable_if<typename T::value_type::LpSolverCol,void>::type
  1029     colBounds(T &t, Value lower, Value upper,dummy<0> = 0) {
  1030       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1031 	colBounds(*i, lower, upper);
  1032       }
  1033     }
  1034     template<class T>
  1035     typename enable_if<typename T::value_type::second_type::LpSolverCol,
  1036 		       void>::type
  1037     colBounds(T &t, Value lower, Value upper,dummy<1> = 1) { 
  1038       for(typename T::iterator i=t.begin();i!=t.end();++i) {
  1039 	colBounds(i->second, lower, upper);
  1040       }
  1041     }
  1042     template<class T>
  1043     typename enable_if<typename T::MapIt::Value::LpSolverCol,
  1044 		       void>::type
  1045     colBounds(T &t, Value lower, Value upper,dummy<2> = 2) { 
  1046       for(typename T::MapIt i(t); i!=INVALID; ++i){
  1047 	colBounds(*i, lower, upper);
  1048       }
  1049     }
  1050 #endif
  1051     
  1052 //     /// Set the lower bound of a row (i.e a constraint)
  1053 
  1054 //     /// The lower bound of a linear expression (row) has to be given by an 
  1055 //     /// extended number of type Value, i.e. a finite number of type 
  1056 //     /// Value or -\ref INF.
  1057 //     void rowLowerBound(Row r, Value value) {
  1058 //       _setRowLowerBound(rows.floatingId(r.id),value);
  1059 //     };
  1060 //     /// Set the upper bound of a row (i.e a constraint)
  1061 
  1062 //     /// The upper bound of a linear expression (row) has to be given by an 
  1063 //     /// extended number of type Value, i.e. a finite number of type 
  1064 //     /// Value or \ref INF.
  1065 //     void rowUpperBound(Row r, Value value) {
  1066 //       _setRowUpperBound(rows.floatingId(r.id),value);
  1067 //     };
  1068 
  1069     /// Set the lower and the upper bounds of a row (i.e a constraint)
  1070 
  1071     /// The lower and the upper bounds of
  1072     /// a constraint (row) have to be given by an 
  1073     /// extended number of type Value, i.e. a finite number of type 
  1074     /// Value, -\ref INF or \ref INF.
  1075     void rowBounds(Row c, Value lower, Value upper) {
  1076       _setRowBounds(rows.floatingId(c.id),lower, upper);
  1077       // _setRowUpperBound(rows.floatingId(c.id),upper);
  1078     }
  1079     
  1080     ///Set an element of the objective function
  1081     void objCoeff(Col c, Value v) {_setObjCoeff(cols.floatingId(c.id),v); };
  1082     ///Set the objective function
  1083     
  1084     ///\param e is a linear expression of type \ref Expr.
  1085     ///\bug Is should be called obj()
  1086     void setObj(Expr e) {
  1087       _clearObj();
  1088       for (Expr::iterator i=e.begin(); i!=e.end(); ++i)
  1089 	objCoeff((*i).first,(*i).second);
  1090       obj_const_comp=e.constComp();
  1091     }
  1092 
  1093     ///Maximize
  1094     void max() { _setMax(); }
  1095     ///Minimize
  1096     void min() { _setMin(); }
  1097 
  1098     
  1099     ///@}
  1100 
  1101 
  1102     ///\name Solve the LP
  1103 
  1104     ///@{
  1105 
  1106     ///\e Solve the LP problem at hand
  1107     ///
  1108     ///\return The result of the optimization procedure. Possible values and their meanings can be found in the documentation of \ref SolveExitStatus.
  1109     ///
  1110     ///\todo Which method is used to solve the problem
  1111     SolveExitStatus solve() { return _solve(); }
  1112     
  1113     ///@}
  1114     
  1115     ///\name Obtain the solution
  1116 
  1117     ///@{
  1118 
  1119     /// The status of the primal problem (the original LP problem)
  1120     SolutionStatus primalStatus() {
  1121       return _getPrimalStatus();
  1122     }
  1123 
  1124     /// The status of the dual (of the original LP) problem 
  1125     SolutionStatus dualStatus() {
  1126       return _getDualStatus();
  1127     }
  1128 
  1129     ///The type of the original LP problem
  1130     ProblemTypes problemType() {
  1131       return _getProblemType();
  1132     }
  1133 
  1134     ///\e
  1135     Value primal(Col c) { return _getPrimal(cols.floatingId(c.id)); }
  1136 
  1137     ///\e
  1138     Value dual(Row r) { return _getDual(rows.floatingId(r.id)); }
  1139 
  1140     ///\e
  1141     bool isBasicCol(Col c) { return _isBasicCol(cols.floatingId(c.id)); }
  1142 
  1143     ///\e
  1144 
  1145     ///\return
  1146     ///- \ref INF or -\ref INF means either infeasibility or unboundedness
  1147     /// of the primal problem, depending on whether we minimize or maximize.
  1148     ///- \ref NaN if no primal solution is found.
  1149     ///- The (finite) objective value if an optimal solution is found.
  1150     Value primalValue() { return _getPrimalValue()+obj_const_comp;}
  1151     ///@}
  1152     
  1153   };  
  1154 
  1155   ///\e
  1156   
  1157   ///\relates LpSolverBase::Expr
  1158   ///
  1159   inline LpSolverBase::Expr operator+(const LpSolverBase::Expr &a,
  1160 				      const LpSolverBase::Expr &b) 
  1161   {
  1162     LpSolverBase::Expr tmp(a);
  1163     tmp+=b;
  1164     return tmp;
  1165   }
  1166   ///\e
  1167   
  1168   ///\relates LpSolverBase::Expr
  1169   ///
  1170   inline LpSolverBase::Expr operator-(const LpSolverBase::Expr &a,
  1171 				      const LpSolverBase::Expr &b) 
  1172   {
  1173     LpSolverBase::Expr tmp(a);
  1174     tmp-=b;
  1175     return tmp;
  1176   }
  1177   ///\e
  1178   
  1179   ///\relates LpSolverBase::Expr
  1180   ///
  1181   inline LpSolverBase::Expr operator*(const LpSolverBase::Expr &a,
  1182 				      const LpSolverBase::Value &b) 
  1183   {
  1184     LpSolverBase::Expr tmp(a);
  1185     tmp*=b;
  1186     return tmp;
  1187   }
  1188   
  1189   ///\e
  1190   
  1191   ///\relates LpSolverBase::Expr
  1192   ///
  1193   inline LpSolverBase::Expr operator*(const LpSolverBase::Value &a,
  1194 				      const LpSolverBase::Expr &b) 
  1195   {
  1196     LpSolverBase::Expr tmp(b);
  1197     tmp*=a;
  1198     return tmp;
  1199   }
  1200   ///\e
  1201   
  1202   ///\relates LpSolverBase::Expr
  1203   ///
  1204   inline LpSolverBase::Expr operator/(const LpSolverBase::Expr &a,
  1205 				      const LpSolverBase::Value &b) 
  1206   {
  1207     LpSolverBase::Expr tmp(a);
  1208     tmp/=b;
  1209     return tmp;
  1210   }
  1211   
  1212   ///\e
  1213   
  1214   ///\relates LpSolverBase::Constr
  1215   ///
  1216   inline LpSolverBase::Constr operator<=(const LpSolverBase::Expr &e,
  1217 					 const LpSolverBase::Expr &f) 
  1218   {
  1219     return LpSolverBase::Constr(-LpSolverBase::INF,e-f,0);
  1220   }
  1221 
  1222   ///\e
  1223   
  1224   ///\relates LpSolverBase::Constr
  1225   ///
  1226   inline LpSolverBase::Constr operator<=(const LpSolverBase::Value &e,
  1227 					 const LpSolverBase::Expr &f) 
  1228   {
  1229     return LpSolverBase::Constr(e,f);
  1230   }
  1231 
  1232   ///\e
  1233   
  1234   ///\relates LpSolverBase::Constr
  1235   ///
  1236   inline LpSolverBase::Constr operator<=(const LpSolverBase::Expr &e,
  1237 					 const LpSolverBase::Value &f) 
  1238   {
  1239     return LpSolverBase::Constr(e,f);
  1240   }
  1241 
  1242   ///\e
  1243   
  1244   ///\relates LpSolverBase::Constr
  1245   ///
  1246   inline LpSolverBase::Constr operator>=(const LpSolverBase::Expr &e,
  1247 					 const LpSolverBase::Expr &f) 
  1248   {
  1249     return LpSolverBase::Constr(-LpSolverBase::INF,f-e,0);
  1250   }
  1251 
  1252 
  1253   ///\e
  1254   
  1255   ///\relates LpSolverBase::Constr
  1256   ///
  1257   inline LpSolverBase::Constr operator>=(const LpSolverBase::Value &e,
  1258 					 const LpSolverBase::Expr &f) 
  1259   {
  1260     return LpSolverBase::Constr(f,e);
  1261   }
  1262 
  1263 
  1264   ///\e
  1265   
  1266   ///\relates LpSolverBase::Constr
  1267   ///
  1268   inline LpSolverBase::Constr operator>=(const LpSolverBase::Expr &e,
  1269 					 const LpSolverBase::Value &f) 
  1270   {
  1271     return LpSolverBase::Constr(f,e);
  1272   }
  1273 
  1274   ///\e
  1275   
  1276   ///\relates LpSolverBase::Constr
  1277   ///
  1278   inline LpSolverBase::Constr operator==(const LpSolverBase::Expr &e,
  1279 					 const LpSolverBase::Expr &f) 
  1280   {
  1281     return LpSolverBase::Constr(0,e-f,0);
  1282   }
  1283 
  1284   ///\e
  1285   
  1286   ///\relates LpSolverBase::Constr
  1287   ///
  1288   inline LpSolverBase::Constr operator<=(const LpSolverBase::Value &n,
  1289 					 const LpSolverBase::Constr&c) 
  1290   {
  1291     LpSolverBase::Constr tmp(c);
  1292     ///\todo Create an own exception type.
  1293     if(!isnan(tmp.lowerBound())) throw LogicError();
  1294     else tmp.lowerBound()=n;
  1295     return tmp;
  1296   }
  1297   ///\e
  1298   
  1299   ///\relates LpSolverBase::Constr
  1300   ///
  1301   inline LpSolverBase::Constr operator<=(const LpSolverBase::Constr& c,
  1302 					 const LpSolverBase::Value &n)
  1303   {
  1304     LpSolverBase::Constr tmp(c);
  1305     ///\todo Create an own exception type.
  1306     if(!isnan(tmp.upperBound())) throw LogicError();
  1307     else tmp.upperBound()=n;
  1308     return tmp;
  1309   }
  1310 
  1311   ///\e
  1312   
  1313   ///\relates LpSolverBase::Constr
  1314   ///
  1315   inline LpSolverBase::Constr operator>=(const LpSolverBase::Value &n,
  1316 					 const LpSolverBase::Constr&c) 
  1317   {
  1318     LpSolverBase::Constr tmp(c);
  1319     ///\todo Create an own exception type.
  1320     if(!isnan(tmp.upperBound())) throw LogicError();
  1321     else tmp.upperBound()=n;
  1322     return tmp;
  1323   }
  1324   ///\e
  1325   
  1326   ///\relates LpSolverBase::Constr
  1327   ///
  1328   inline LpSolverBase::Constr operator>=(const LpSolverBase::Constr& c,
  1329 					 const LpSolverBase::Value &n)
  1330   {
  1331     LpSolverBase::Constr tmp(c);
  1332     ///\todo Create an own exception type.
  1333     if(!isnan(tmp.lowerBound())) throw LogicError();
  1334     else tmp.lowerBound()=n;
  1335     return tmp;
  1336   }
  1337 
  1338   ///\e
  1339   
  1340   ///\relates LpSolverBase::DualExpr
  1341   ///
  1342   inline LpSolverBase::DualExpr operator+(const LpSolverBase::DualExpr &a,
  1343 				      const LpSolverBase::DualExpr &b) 
  1344   {
  1345     LpSolverBase::DualExpr tmp(a);
  1346     tmp+=b;
  1347     return tmp;
  1348   }
  1349   ///\e
  1350   
  1351   ///\relates LpSolverBase::DualExpr
  1352   ///
  1353   inline LpSolverBase::DualExpr operator-(const LpSolverBase::DualExpr &a,
  1354 				      const LpSolverBase::DualExpr &b) 
  1355   {
  1356     LpSolverBase::DualExpr tmp(a);
  1357     tmp-=b;
  1358     return tmp;
  1359   }
  1360   ///\e
  1361   
  1362   ///\relates LpSolverBase::DualExpr
  1363   ///
  1364   inline LpSolverBase::DualExpr operator*(const LpSolverBase::DualExpr &a,
  1365 				      const LpSolverBase::Value &b) 
  1366   {
  1367     LpSolverBase::DualExpr tmp(a);
  1368     tmp*=b;
  1369     return tmp;
  1370   }
  1371   
  1372   ///\e
  1373   
  1374   ///\relates LpSolverBase::DualExpr
  1375   ///
  1376   inline LpSolverBase::DualExpr operator*(const LpSolverBase::Value &a,
  1377 				      const LpSolverBase::DualExpr &b) 
  1378   {
  1379     LpSolverBase::DualExpr tmp(b);
  1380     tmp*=a;
  1381     return tmp;
  1382   }
  1383   ///\e
  1384   
  1385   ///\relates LpSolverBase::DualExpr
  1386   ///
  1387   inline LpSolverBase::DualExpr operator/(const LpSolverBase::DualExpr &a,
  1388 				      const LpSolverBase::Value &b) 
  1389   {
  1390     LpSolverBase::DualExpr tmp(a);
  1391     tmp/=b;
  1392     return tmp;
  1393   }
  1394   
  1395 
  1396 } //namespace lemon
  1397 
  1398 #endif //LEMON_LP_BASE_H