src/work/athos/lp/lp_base.h
author alpar
Fri, 01 Apr 2005 19:50:29 +0000
changeset 1293 8ede2a6b2594
parent 1291 16cde3e1aa9f
child 1294 2dec219d9ca2
permissions -rw-r--r--
- Make lp stuff compilable
- Some 'set's removed
     1 /* -*- C++ -*-
     2  * src/lemon/lp_base.h - Part of LEMON, a generic C++ optimization library
     3  *
     4  * Copyright (C) 2005 Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     5  * (Egervary Combinatorial Optimization Research Group, EGRES).
     6  *
     7  * Permission to use, modify and distribute this software is granted
     8  * provided that this copyright notice appears in all copies. For
     9  * precise terms see the accompanying LICENSE file.
    10  *
    11  * This software is provided "AS IS" with no warranty of any kind,
    12  * express or implied, and with no claim as to its suitability for any
    13  * purpose.
    14  *
    15  */
    16 
    17 #ifndef LEMON_LP_BASE_H
    18 #define LEMON_LP_BASE_H
    19 
    20 #include<vector>
    21 #include<map>
    22 #include<limits>
    23 #include<math.h>
    24 
    25 #include<lemon/utility.h>
    26 #include<lemon/error.h>
    27 #include<lemon/invalid.h>
    28 
    29 //#include"lin_expr.h"
    30 
    31 ///\file
    32 ///\brief The interface of the LP solver interface.
    33 namespace lemon {
    34   
    35   ///Internal data structure to convert floating id's to fix one's
    36     
    37   ///\todo This might be implemented to be also usable in other places.
    38   class _FixId 
    39   {
    40     std::vector<int> index;
    41     std::vector<int> cross;
    42     int first_free;
    43   public:
    44     _FixId() : first_free(-1) {};
    45     ///Convert a floating id to a fix one
    46 
    47     ///\param n is a floating id
    48     ///\return the corresponding fix id
    49     int fixId(int n) {return cross[n];}
    50     ///Convert a fix id to a floating one
    51 
    52     ///\param n is a fix id
    53     ///\return the corresponding floating id
    54     int floatingId(int n) { return index[n];}
    55     ///Add a new floating id.
    56 
    57     ///\param n is a floating id
    58     ///\return the fix id of the new value
    59     ///\todo Multiple additions should also be handled.
    60     int insert(int n)
    61     {
    62       if(n>=int(cross.size())) {
    63 	cross.resize(n+1);
    64 	if(first_free==-1) {
    65 	  cross[n]=index.size();
    66 	  index.push_back(n);
    67 	}
    68 	else {
    69 	  cross[n]=first_free;
    70 	  int next=index[first_free];
    71 	  index[first_free]=n;
    72 	  first_free=next;
    73 	}
    74 	return cross[n];
    75       }
    76       ///\todo Create an own exception type.
    77       else throw LogicError(); //floatingId-s must form a continuous range;
    78     }
    79     ///Remove a fix id.
    80 
    81     ///\param n is a fix id
    82     ///
    83     void erase(int n) 
    84     {
    85       int fl=index[n];
    86       index[n]=first_free;
    87       first_free=n;
    88       for(int i=fl+1;i<int(cross.size());++i) {
    89 	cross[i-1]=cross[i];
    90 	index[cross[i]]--;
    91       }
    92       cross.pop_back();
    93     }
    94     ///An upper bound on the largest fix id.
    95 
    96     ///\todo Do we need this?
    97     ///
    98     std::size_t maxFixId() { return cross.size()-1; }
    99   
   100   };
   101     
   102   ///Common base class for LP solvers
   103   class LpSolverBase {
   104     
   105   public:
   106 
   107     ///\e
   108     enum SolutionStatus {
   109       ///\e
   110       SOLVED = 0,
   111       ///\e
   112       UNSOLVED = 1
   113     };
   114       
   115     ///\e
   116     enum SolutionType {
   117       ///\e
   118       UNDEFINED = 0,
   119       ///\e
   120       INFEASIBLE = 1,
   121       ///\e
   122       FEASIBLE = 2,
   123       ///\e
   124       OPTIMAL = 3
   125     };
   126       
   127     ///The floating point type used by the solver
   128     typedef double Value;
   129     ///The infinity constant
   130     static const Value INF;
   131     ///The not a number constant
   132     static const Value NaN;
   133     
   134     ///Refer to a column of the LP.
   135 
   136     ///This type is used to refer to a column of the LP.
   137     ///
   138     ///Its value remains valid and correct even after the addition or erase of
   139     ///other columns.
   140     ///
   141     ///\todo Document what can one do with a Col (INVALID, comparing,
   142     ///it is similar to Node/Edge)
   143     class Col {
   144     protected:
   145       int id;
   146       friend class LpSolverBase;
   147     public:
   148       typedef Value ExprValue;
   149       typedef True LpSolverCol;
   150       Col() {}
   151       Col(const Invalid&) : id(-1) {}
   152       bool operator<(Col c) const  {return id<c.id;}
   153       bool operator==(Col c) const  {return id==c.id;}
   154       bool operator!=(Col c) const  {return id==c.id;}
   155     };
   156 
   157     ///Refer to a row of the LP.
   158 
   159     ///This type is used to refer to a row of the LP.
   160     ///
   161     ///Its value remains valid and correct even after the addition or erase of
   162     ///other rows.
   163     ///
   164     ///\todo Document what can one do with a Row (INVALID, comparing,
   165     ///it is similar to Node/Edge)
   166     class Row {
   167     protected:
   168       int id;
   169       friend class LpSolverBase;
   170     public:
   171       typedef Value ExprValue;
   172       typedef True LpSolverRow;
   173       Row() {}
   174       Row(const Invalid&) : id(-1) {}
   175       typedef True LpSolverRow;
   176       bool operator<(Row c) const  {return id<c.id;}
   177       bool operator==(Row c) const  {return id==c.id;}
   178       bool operator!=(Row c) const  {return id==c.id;} 
   179    };
   180     
   181     ///Linear expression of variables and a constant component
   182     
   183     ///This data structure strores a linear expression of the variables
   184     ///(\ref Col "Col"s) and also has a constant component.
   185     ///
   186     ///There are several ways to access and modify the contents of this
   187     ///container.
   188     ///- Its it fully compatible with \c std::map<Col,double>, so for expamle
   189     ///if \c e is an Expr and \c v and \c w are of type \ref Col then you can
   190     ///read and modify the coefficients like
   191     ///these.
   192     ///\code
   193     ///e[v]=5;
   194     ///e[v]+=12;
   195     ///e.erase(v);
   196     ///\endcode
   197     ///or you can also iterate through its elements.
   198     ///\code
   199     ///double s=0;
   200     ///for(LpSolverBase::Expr::iterator i=e.begin();i!=e.end();++i)
   201     ///  s+=i->second;
   202     ///\endcode
   203     ///(This code computes the sum of all coefficients).
   204     ///- Numbers (<tt>double</tt>'s)
   205     ///and variables (\ref Col "Col"s) directly convert to an
   206     ///\ref Expr and the usual linear operations are defined so  
   207     ///\code
   208     ///v+w
   209     ///2*v-3.12*(v-w/2)+2
   210     ///v*2.1+(3*v+(v*12+w+6)*3)/2
   211     ///\endcode
   212     ///are valid expressions. The usual assignment operations are also defined.
   213     ///\code
   214     ///e=v+w;
   215     ///e+=2*v-3.12*(v-w/2)+2;
   216     ///e*=3.4;
   217     ///e/=5;
   218     ///\endcode
   219     ///- The constant member can be set and read by \ref constComp()
   220     ///\code
   221     ///e.constComp()=12;
   222     ///double c=e.constComp();
   223     ///\endcode
   224     ///
   225     ///\note that \ref clear() not only sets all coefficients to 0 but also
   226     ///clears the constant components.
   227     class Expr : public std::map<Col,Value>
   228     {
   229     public:
   230       typedef LpSolverBase::Col Key; 
   231       typedef LpSolverBase::Value Value;
   232       
   233     protected:
   234       typedef std::map<Col,Value> Base;
   235       
   236       Value const_comp;
   237   public:
   238       typedef True IsLinExpression;
   239       ///\e
   240       Expr() : Base(), const_comp(0) { }
   241       ///\e
   242       Expr(const Key &v) : const_comp(0) {
   243 	Base::insert(std::make_pair(v, 1));
   244       }
   245       ///\e
   246       Expr(const Value &v) : const_comp(v) {}
   247       ///\e
   248       void set(const Key &v,const Value &c) {
   249 	Base::insert(std::make_pair(v, c));
   250       }
   251       ///\e
   252       Value &constComp() { return const_comp; }
   253       ///\e
   254       const Value &constComp() const { return const_comp; }
   255       
   256       ///Removes the components with zero coefficient.
   257       void simplify() {
   258 	for (Base::iterator i=Base::begin(); i!=Base::end();) {
   259 	  Base::iterator j=i;
   260 	  ++j;
   261 	  if ((*i).second==0) Base::erase(i);
   262 	  j=i;
   263 	}
   264       }
   265 
   266       ///Sets all coefficients and the constant component to 0.
   267       void clear() {
   268 	Base::clear();
   269 	const_comp=0;
   270       }
   271 
   272       ///\e
   273       Expr &operator+=(const Expr &e) {
   274 	for (Base::const_iterator j=e.begin(); j!=e.end(); ++j)
   275 	  (*this)[j->first]+=j->second;
   276 	///\todo it might be speeded up using "hints"
   277 	const_comp+=e.const_comp;
   278 	return *this;
   279       }
   280       ///\e
   281       Expr &operator-=(const Expr &e) {
   282 	for (Base::const_iterator j=e.begin(); j!=e.end(); ++j)
   283 	  (*this)[j->first]-=j->second;
   284 	const_comp-=e.const_comp;
   285 	return *this;
   286       }
   287       ///\e
   288       Expr &operator*=(const Value &c) {
   289 	for (Base::iterator j=Base::begin(); j!=Base::end(); ++j)
   290 	  j->second*=c;
   291 	const_comp*=c;
   292 	return *this;
   293       }
   294       ///\e
   295       Expr &operator/=(const Value &c) {
   296 	for (Base::iterator j=Base::begin(); j!=Base::end(); ++j)
   297 	  j->second/=c;
   298 	const_comp/=c;
   299 	return *this;
   300       }
   301     };
   302     
   303     ///Linear constraint
   304     //typedef LinConstr<Expr> Constr;
   305     class Constr
   306     {
   307     public:
   308       typedef LpSolverBase::Expr Expr;
   309       typedef Expr::Key Key;
   310       typedef Expr::Value Value;
   311       
   312       static const Value INF;
   313       static const Value NaN;
   314       //     static const Value INF=0;
   315       //     static const Value NaN=1;
   316       
   317     protected:
   318       Expr _expr;
   319       Value _lb,_ub;
   320     public:
   321       ///\e
   322       Constr() : _expr(), _lb(NaN), _ub(NaN) {}
   323       ///\e
   324       Constr(Value lb,const Expr &e,Value ub) :
   325 	_expr(e), _lb(lb), _ub(ub) {}
   326       ///\e
   327       Constr(const Expr &e,Value ub) : 
   328 	_expr(e), _lb(NaN), _ub(ub) {}
   329       ///\e
   330       Constr(Value lb,const Expr &e) :
   331 	_expr(e), _lb(lb), _ub(NaN) {}
   332       ///\e
   333       Constr(const Expr &e) : 
   334 	_expr(e), _lb(NaN), _ub(NaN) {}
   335       ///\e
   336       void clear() 
   337       {
   338 	_expr.clear();
   339 	_lb=_ub=NaN;
   340       }
   341       ///\e
   342       Expr &expr() { return _expr; }
   343       ///\e
   344       const Expr &expr() const { return _expr; }
   345       ///\e
   346       Value &lowerBound() { return _lb; }
   347       ///\e
   348       const Value &lowerBound() const { return _lb; }
   349       ///\e
   350       Value &upperBound() { return _ub; }
   351       ///\e
   352       const Value &upperBound() const { return _ub; }
   353       ///\e
   354       bool lowerBounded() const { return std::isfinite(_lb); }
   355       ///\e
   356       bool upperBounded() const { return std::isfinite(_ub); }
   357     };
   358     
   359 
   360   protected:
   361     _FixId rows;
   362     _FixId cols;
   363 
   364     /// \e
   365     virtual int _addCol() = 0;
   366     /// \e
   367     virtual int _addRow() = 0;
   368     /// \e
   369 
   370     /// \warning Arrays are indexed from 1 (datum at index 0 is ignored)
   371     ///
   372     virtual void _setRowCoeffs(int i, 
   373 			       int length,
   374                                int  const * indices, 
   375                                Value  const * values ) = 0;
   376     /// \e
   377 
   378     /// \warning Arrays are indexed from 1 (datum at index 0 is ignored)
   379     ///
   380     virtual void _setColCoeffs(int i, 
   381 			       int length,
   382                                int  const * indices, 
   383                                Value  const * values ) = 0;
   384     
   385     /// \e
   386 
   387     /// The lower bound of a variable (column) have to be given by an 
   388     /// extended number of type Value, i.e. a finite number of type 
   389     /// Value or -\ref INF.
   390     virtual void _setColLowerBound(int i, Value value) = 0;
   391     /// \e
   392 
   393     /// The upper bound of a variable (column) have to be given by an 
   394     /// extended number of type Value, i.e. a finite number of type 
   395     /// Value or \ref INF.
   396     virtual void _setColUpperBound(int i, Value value) = 0;
   397     /// \e
   398 
   399     /// The lower bound of a linear expression (row) have to be given by an 
   400     /// extended number of type Value, i.e. a finite number of type 
   401     /// Value or -\ref INF.
   402     virtual void _setRowLowerBound(int i, Value value) = 0;
   403     /// \e
   404 
   405     /// The upper bound of a linear expression (row) have to be given by an 
   406     /// extended number of type Value, i.e. a finite number of type 
   407     /// Value or \ref INF.
   408     virtual void _setRowUpperBound(int i, Value value) = 0;
   409 
   410     /// \e
   411     virtual void _setObjCoeff(int i, Value obj_coef) = 0;
   412 
   413     ///\e
   414     
   415     ///\bug Wrong interface
   416     ///
   417     virtual SolutionStatus _solve() = 0;
   418 
   419     ///\e
   420 
   421     ///\bug Wrong interface
   422     ///
   423     virtual Value _getPrimal(int i) = 0;
   424     ///\e
   425 
   426     ///\bug unimplemented!!!!
   427     void clearObj() {}
   428   public:
   429 
   430 
   431     ///\e
   432     virtual ~LpSolverBase() {}
   433 
   434     ///\name Building up and modification of the LP
   435 
   436     ///@{
   437 
   438     ///Add a new empty column (i.e a new variable) to the LP
   439     Col addCol() { Col c; c.id=cols.insert(_addCol()); return c;}
   440 
   441     ///\brief Fill the elements of a container with newly created columns
   442     ///(i.e a new variables)
   443     ///
   444     ///This magic function takes a container as its argument
   445     ///and fills its elements
   446     ///with new columns (i.e. variables)
   447     ///\param t can be
   448     ///- a standard STL compatible iterable container with
   449     ///\ref Col as its \c values_type
   450     ///like
   451     ///\code
   452     ///std::vector<LpSolverBase::Col>
   453     ///std::list<LpSolverBase::Col>
   454     ///\endcode
   455     ///- a standard STL compatible iterable container with
   456     ///\ref Col as its \c mapped_type
   457     ///like
   458     ///\code
   459     ///std::map<AnyType,LpSolverBase::Col>
   460     ///\endcode
   461     ///- an iterable lemon \ref concept::WriteMap "write map" like 
   462     ///\code
   463     ///ListGraph::NodeMap<LpSolverBase::Col>
   464     ///ListGraph::EdgeMap<LpSolverBase::Col>
   465     ///\endcode
   466     ///\return The number of the created column.
   467     ///\bug Iterable nodemap hasn't been implemented yet.
   468 #ifdef DOXYGEN
   469     template<class T>
   470     int addColSet(T &t) { return 0;} 
   471 #else
   472     template<class T>
   473     typename enable_if<typename T::value_type::LpSolverCol,int>::type
   474     addColSet(T &t,dummy<0> = 0) {
   475       int s=0;
   476       for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addCol();s++;}
   477       return s;
   478     }
   479     template<class T>
   480     typename enable_if<typename T::value_type::second_type::LpSolverCol,
   481 		       int>::type
   482     addColSet(T &t,dummy<1> = 1) { 
   483       int s=0;
   484       for(typename T::iterator i=t.begin();i!=t.end();++i) {
   485 	i->second=addCol();
   486 	s++;
   487       }
   488       return s;
   489     }
   490     template<class T>
   491     typename enable_if<typename T::ValueSet::value_type::LpSolverCol,
   492 		       int>::type
   493     addColSet(T &t,dummy<2> = 2) { 
   494       ///\bug <tt>return addColSet(t.valueSet());</tt> should also work.
   495       int s=0;
   496       for(typename T::ValueSet::iterator i=t.valueSet().begin();
   497 	  i!=t.valueSet().end();
   498 	  ++i)
   499 	{
   500 	  *i=addCol();
   501 	  s++;
   502 	}
   503       return s;
   504     }
   505 #endif
   506 
   507     ///Add a new empty row (i.e a new constaint) to the LP
   508 
   509     ///This function adds a new empty row (i.e a new constaint) to the LP.
   510     ///\return The created row
   511     Row addRow() { Row r; r.id=rows.insert(_addRow()); return r;}
   512 
   513     ///Set a row (i.e a constaint) of the LP
   514 
   515     ///\param r is the row to be modified
   516     ///\param l is lower bound (-\ref INF means no bound)
   517     ///\param e is a linear expression (see \ref Expr)
   518     ///\param u is the upper bound (\ref INF means no bound)
   519     ///\bug This is a temportary function. The interface will change to
   520     ///a better one.
   521     void setRow(Row r, Value l,const Expr &e, Value u) {
   522       std::vector<int> indices;
   523       std::vector<Value> values;
   524       indices.push_back(0);
   525       values.push_back(0);
   526       for(Expr::const_iterator i=e.begin(); i!=e.end(); ++i)
   527 	if((*i).second!=0) { ///\bug EPSILON would be necessary here!!!
   528 	  indices.push_back(cols.floatingId((*i).first.id));
   529 	  values.push_back((*i).second);
   530 	}
   531       _setRowCoeffs(rows.floatingId(r.id),indices.size()-1,
   532 		    &indices[0],&values[0]);
   533       _setRowLowerBound(rows.floatingId(r.id),l-e.constComp());
   534       _setRowUpperBound(rows.floatingId(r.id),u-e.constComp());
   535     }
   536 
   537     ///Set a row (i.e a constaint) of the LP
   538 
   539     ///\param r is the row to be modified
   540     ///\param c is a linear expression (see \ref Constr)
   541     void setRow(Row r, const Constr &c) {
   542       setRow(r,
   543 	     c.lowerBounded()?c.lowerBound():-INF,
   544 	     c.expr(),
   545 	     c.upperBounded()?c.upperBound():INF);
   546     }
   547 
   548     ///Add a new row (i.e a new constaint) to the LP
   549 
   550     ///\param l is the lower bound (-\ref INF means no bound)
   551     ///\param e is a linear expression (see \ref Expr)
   552     ///\param u is the upper bound (\ref INF means no bound)
   553     ///\return The created row.
   554     ///\bug This is a temportary function. The interface will change to
   555     ///a better one.
   556     Row addRow(Value l,const Expr &e, Value u) {
   557       Row r=addRow();
   558       setRow(r,l,e,u);
   559       return r;
   560     }
   561 
   562     ///Add a new row (i.e a new constaint) to the LP
   563 
   564     ///\param c is a linear expression (see \ref Constr)
   565     ///\return The created row.
   566     Row addRow(const Constr &c) {
   567       Row r=addRow();
   568       setRow(r,c);
   569       return r;
   570     }
   571 
   572     /// Set the lower bound of a column (i.e a variable)
   573 
   574     /// The upper bound of a variable (column) has to be given by an 
   575     /// extended number of type Value, i.e. a finite number of type 
   576     /// Value or -\ref INF.
   577     void colLowerBound(Col c, Value value) {
   578       _setColLowerBound(cols.floatingId(c.id),value);
   579     }
   580     /// Set the upper bound of a column (i.e a variable)
   581 
   582     /// The upper bound of a variable (column) has to be given by an 
   583     /// extended number of type Value, i.e. a finite number of type 
   584     /// Value or \ref INF.
   585     void colUpperBound(Col c, Value value) {
   586       _setColUpperBound(cols.floatingId(c.id),value);
   587     };
   588     /// Set the lower and the upper bounds of a column (i.e a variable)
   589 
   590     /// The lower and the upper bounds of
   591     /// a variable (column) have to be given by an 
   592     /// extended number of type Value, i.e. a finite number of type 
   593     /// Value, -\ref INF or \ref INF.
   594     void colBounds(Col c, Value lower, Value upper) {
   595       _setColLowerBound(cols.floatingId(c.id),lower);
   596       _setColUpperBound(cols.floatingId(c.id),upper);
   597     }
   598     
   599     /// Set the lower bound of a row (i.e a constraint)
   600 
   601     /// The lower bound of a linear expression (row) has to be given by an 
   602     /// extended number of type Value, i.e. a finite number of type 
   603     /// Value or -\ref INF.
   604     void rowLowerBound(Row r, Value value) {
   605       _setRowLowerBound(rows.floatingId(r.id),value);
   606     };
   607     /// Set the upper bound of a row (i.e a constraint)
   608 
   609     /// The upper bound of a linear expression (row) has to be given by an 
   610     /// extended number of type Value, i.e. a finite number of type 
   611     /// Value or \ref INF.
   612     void rowUpperBound(Row r, Value value) {
   613       _setRowUpperBound(rows.floatingId(r.id),value);
   614     };
   615     /// Set the lower and the upper bounds of a row (i.e a variable)
   616 
   617     /// The lower and the upper bounds of
   618     /// a constraint (row) have to be given by an 
   619     /// extended number of type Value, i.e. a finite number of type 
   620     /// Value, -\ref INF or \ref INF.
   621     void rowBounds(Row c, Value lower, Value upper) {
   622       _setRowLowerBound(rows.floatingId(c.id),lower);
   623       _setRowUpperBound(rows.floatingId(c.id),upper);
   624     }
   625     
   626     ///Set an element of the objective function
   627     void objCoeff(Col c, Value v) {_setObjCoeff(cols.floatingId(c.id),v); };
   628     ///Set the objective function
   629     
   630     ///\param e is a linear expression of type \ref Expr.
   631     ///\todo What to do with the constant component?
   632     void setObj(Expr e) {
   633       clearObj();
   634       for (Expr::iterator i=e.begin(); i!=e.end(); ++i)
   635 	objCoeff((*i).first,(*i).second);
   636     }
   637 
   638     ///@}
   639 
   640 
   641     ///\name Solving the LP
   642 
   643     ///@{
   644 
   645     ///\e
   646     SolutionStatus solve() { return _solve(); }
   647     
   648     ///@}
   649     
   650     ///\name Obtaining the solution LP
   651 
   652     ///@{
   653 
   654     ///\e
   655     Value primal(Col c) { return _getPrimal(cols.floatingId(c.id)); }
   656 
   657     ///@}
   658     
   659   };  
   660 
   661   ///\e
   662   
   663   ///\relates LpSolverBase::Expr
   664   ///
   665   inline LpSolverBase::Expr operator+(const LpSolverBase::Expr &a,
   666 				      const LpSolverBase::Expr &b) 
   667   {
   668     LpSolverBase::Expr tmp(a);
   669     tmp+=b; ///\todo Don't STL have some special 'merge' algorithm?
   670     return tmp;
   671   }
   672   ///\e
   673   
   674   ///\relates LpSolverBase::Expr
   675   ///
   676   inline LpSolverBase::Expr operator-(const LpSolverBase::Expr &a,
   677 				      const LpSolverBase::Expr &b) 
   678   {
   679     LpSolverBase::Expr tmp(a);
   680     tmp-=b; ///\todo Don't STL have some special 'merge' algorithm?
   681     return tmp;
   682   }
   683   ///\e
   684   
   685   ///\relates LpSolverBase::Expr
   686   ///
   687   inline LpSolverBase::Expr operator*(const LpSolverBase::Expr &a,
   688 				      const LpSolverBase::Value &b) 
   689   {
   690     LpSolverBase::Expr tmp(a);
   691     tmp*=b; ///\todo Don't STL have some special 'merge' algorithm?
   692     return tmp;
   693   }
   694   
   695   ///\e
   696   
   697   ///\relates LpSolverBase::Expr
   698   ///
   699   inline LpSolverBase::Expr operator*(const LpSolverBase::Value &a,
   700 				      const LpSolverBase::Expr &b) 
   701   {
   702     LpSolverBase::Expr tmp(b);
   703     tmp*=a; ///\todo Don't STL have some special 'merge' algorithm?
   704     return tmp;
   705   }
   706   ///\e
   707   
   708   ///\relates LpSolverBase::Expr
   709   ///
   710   inline LpSolverBase::Expr operator/(const LpSolverBase::Expr &a,
   711 				      const LpSolverBase::Value &b) 
   712   {
   713     LpSolverBase::Expr tmp(a);
   714     tmp/=b; ///\todo Don't STL have some special 'merge' algorithm?
   715     return tmp;
   716   }
   717   
   718   ///\e
   719   
   720   ///\relates LpSolverBase::Constr
   721   ///
   722   inline LpSolverBase::Constr operator<=(const LpSolverBase::Expr &e,
   723 					 const LpSolverBase::Expr &f) 
   724   {
   725     return LpSolverBase::Constr(-LpSolverBase::INF,e-f,0);
   726   }
   727 
   728   ///\e
   729   
   730   ///\relates LpSolverBase::Constr
   731   ///
   732   inline LpSolverBase::Constr operator<=(const LpSolverBase::Value &e,
   733 					 const LpSolverBase::Expr &f) 
   734   {
   735     return LpSolverBase::Constr(e,f);
   736   }
   737 
   738   ///\e
   739   
   740   ///\relates LpSolverBase::Constr
   741   ///
   742   inline LpSolverBase::Constr operator<=(const LpSolverBase::Expr &e,
   743 					 const LpSolverBase::Value &f) 
   744   {
   745     return LpSolverBase::Constr(e,f);
   746   }
   747 
   748   ///\e
   749   
   750   ///\relates LpSolverBase::Constr
   751   ///
   752   inline LpSolverBase::Constr operator>=(const LpSolverBase::Expr &e,
   753 					 const LpSolverBase::Expr &f) 
   754   {
   755     return LpSolverBase::Constr(-LpSolverBase::INF,f-e,0);
   756   }
   757 
   758 
   759   ///\e
   760   
   761   ///\relates LpSolverBase::Constr
   762   ///
   763   inline LpSolverBase::Constr operator>=(const LpSolverBase::Value &e,
   764 					 const LpSolverBase::Expr &f) 
   765   {
   766     return LpSolverBase::Constr(f,e);
   767   }
   768 
   769 
   770   ///\e
   771   
   772   ///\relates LpSolverBase::Constr
   773   ///
   774   inline LpSolverBase::Constr operator>=(const LpSolverBase::Expr &e,
   775 					 const LpSolverBase::Value &f) 
   776   {
   777     return LpSolverBase::Constr(f,e);
   778   }
   779 
   780   ///\e
   781   
   782   ///\relates LpSolverBase::Constr
   783   ///
   784   inline LpSolverBase::Constr operator==(const LpSolverBase::Expr &e,
   785 					 const LpSolverBase::Expr &f) 
   786   {
   787     return LpSolverBase::Constr(0,e-f,0);
   788   }
   789 
   790   ///\e
   791   
   792   ///\relates LpSolverBase::Constr
   793   ///
   794   inline LpSolverBase::Constr operator<=(const LpSolverBase::Value &n,
   795 					 const LpSolverBase::Constr&c) 
   796   {
   797     LpSolverBase::Constr tmp(c);
   798     ///\todo Create an own exception type.
   799     if(!isnan(tmp.lowerBound())) throw LogicError();
   800     else tmp.lowerBound()=n;
   801     return tmp;
   802   }
   803   ///\e
   804   
   805   ///\relates LpSolverBase::Constr
   806   ///
   807   inline LpSolverBase::Constr operator<=(const LpSolverBase::Constr& c,
   808 					 const LpSolverBase::Value &n)
   809   {
   810     LpSolverBase::Constr tmp(c);
   811     ///\todo Create an own exception type.
   812     if(!isnan(tmp.upperBound())) throw LogicError();
   813     else tmp.upperBound()=n;
   814     return tmp;
   815   }
   816 
   817   ///\e
   818   
   819   ///\relates LpSolverBase::Constr
   820   ///
   821   inline LpSolverBase::Constr operator>=(const LpSolverBase::Value &n,
   822 					 const LpSolverBase::Constr&c) 
   823   {
   824     LpSolverBase::Constr tmp(c);
   825     ///\todo Create an own exception type.
   826     if(!isnan(tmp.upperBound())) throw LogicError();
   827     else tmp.upperBound()=n;
   828     return tmp;
   829   }
   830   ///\e
   831   
   832   ///\relates LpSolverBase::Constr
   833   ///
   834   inline LpSolverBase::Constr operator>=(const LpSolverBase::Constr& c,
   835 					 const LpSolverBase::Value &n)
   836   {
   837     LpSolverBase::Constr tmp(c);
   838     ///\todo Create an own exception type.
   839     if(!isnan(tmp.lowerBound())) throw LogicError();
   840     else tmp.lowerBound()=n;
   841     return tmp;
   842   }
   843 
   844 
   845 } //namespace lemon
   846 
   847 #endif //LEMON_LP_BASE_H