COIN-OR::LEMON - Graph Library

source: lemon-0.x/lemon/lp_base.h @ 1458:7a483c1d38b5

Last change on this file since 1458:7a483c1d38b5 was 1458:7a483c1d38b5, checked in by athos, 15 years ago

Not ready, but I commit it for simplicity.

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