lemon/lp_base.h
changeset 2312 07e46cbb7d85
parent 2309 468a525d5b45
child 2324 18fc834761d9
     1.1 --- a/lemon/lp_base.h	Tue Nov 28 17:25:22 2006 +0000
     1.2 +++ b/lemon/lp_base.h	Wed Nov 29 15:01:13 2006 +0000
     1.3 @@ -32,7 +32,8 @@
     1.4  ///\brief The interface of the LP solver interface.
     1.5  ///\ingroup gen_opt_group
     1.6  namespace lemon {
     1.7 -  
     1.8 +
     1.9 +
    1.10    ///Internal data structure to convert floating id's to fix one's
    1.11      
    1.12    ///\todo This might be implemented to be also usable in other places.
    1.13 @@ -109,7 +110,7 @@
    1.14      ///or -1 if no index has been inserted before.
    1.15      int firstIndex() {return _first_index;}
    1.16    };
    1.17 -    
    1.18 +
    1.19    ///Common base class for LP solvers
    1.20    
    1.21    ///\todo Much more docs
    1.22 @@ -128,7 +129,8 @@
    1.23        ///an optimal solution has been found or infeasibility/unboundedness
    1.24        ///has been proved.
    1.25        SOLVED = 0,
    1.26 -      ///Any other case (including the case when some user specified limit has been exceeded)
    1.27 +      ///Any other case (including the case when some user specified
    1.28 +      ///limit has been exceeded)
    1.29        UNSOLVED = 1
    1.30      };
    1.31        
    1.32 @@ -221,6 +223,9 @@
    1.33  	return *this;
    1.34        }
    1.35      };
    1.36 +
    1.37 +    static int id(const Col& col) { return col.id; }
    1.38 + 
    1.39        
    1.40      ///Refer to a row of the LP.
    1.41  
    1.42 @@ -245,7 +250,22 @@
    1.43        bool operator> (Row c) const  {return id> c.id;}
    1.44        bool operator==(Row c) const  {return id==c.id;}
    1.45        bool operator!=(Row c) const  {return id!=c.id;} 
    1.46 -   };
    1.47 +    };
    1.48 +
    1.49 +    static int id(const Row& row) { return row.id; }
    1.50 +
    1.51 +  protected:
    1.52 +
    1.53 +    int _lpId(const Col& col) const {
    1.54 +      return cols.floatingId(id(col));
    1.55 +    }
    1.56 +
    1.57 +    int _lpId(const Row& row) const {
    1.58 +      return rows.floatingId(id(row));
    1.59 +    }
    1.60 +
    1.61 +
    1.62 +  public:
    1.63      
    1.64      ///Linear expression of variables and a constant component
    1.65      
    1.66 @@ -336,6 +356,10 @@
    1.67  	}
    1.68        }
    1.69  
    1.70 +      void simplify() const {
    1.71 +        const_cast<Expr*>(this)->simplify();
    1.72 +      }
    1.73 +
    1.74        ///Removes the coefficients closer to zero than \c tolerance.
    1.75        void simplify(double &tolerance) {
    1.76  	for (Base::iterator i=Base::begin(); i!=Base::end();) {
    1.77 @@ -415,9 +439,6 @@
    1.78        typedef Expr::Key Key;
    1.79        typedef Expr::Value Value;
    1.80        
    1.81 -//       static const Value INF;
    1.82 -//       static const Value NaN;
    1.83 -
    1.84      protected:
    1.85        Expr _expr;
    1.86        Value _lb,_ub;
    1.87 @@ -553,6 +574,10 @@
    1.88  	}
    1.89        }
    1.90  
    1.91 +      void simplify() const {
    1.92 +        const_cast<DualExpr*>(this)->simplify();
    1.93 +      }
    1.94 +
    1.95        ///Removes the coefficients closer to zero than \c tolerance.
    1.96        void simplify(double &tolerance) {
    1.97  	for (Base::iterator i=Base::begin(); i!=Base::end();) {
    1.98 @@ -563,7 +588,6 @@
    1.99  	}
   1.100        }
   1.101  
   1.102 -
   1.103        ///Sets all coefficients to 0.
   1.104        void clear() {
   1.105  	Base::clear();
   1.106 @@ -596,12 +620,73 @@
   1.107      };
   1.108      
   1.109  
   1.110 +  private:
   1.111 +
   1.112 +    template <typename _Base>
   1.113 +    class MappedIterator {
   1.114 +    public:
   1.115 +
   1.116 +      typedef _Base Base;
   1.117 +
   1.118 +      typedef typename Base::iterator_category iterator_category;
   1.119 +      typedef typename Base::difference_type difference_type;
   1.120 +      typedef const std::pair<int, Value> value_type;
   1.121 +      typedef value_type reference;
   1.122 +      class pointer {
   1.123 +      public:
   1.124 +        pointer(value_type& _value) : value(_value) {}
   1.125 +        value_type* operator->() { return &value; }
   1.126 +      private:
   1.127 +        value_type value;
   1.128 +      };
   1.129 +
   1.130 +      MappedIterator(const Base& _base, const LpSolverBase& _lp) 
   1.131 +        : base(_base), lp(_lp) {}
   1.132 +
   1.133 +      reference operator*() {
   1.134 +        return std::make_pair(lp._lpId(base->first), base->second);
   1.135 +      }
   1.136 +
   1.137 +      pointer operator->() {
   1.138 +        return pointer(operator*());
   1.139 +      }
   1.140 +
   1.141 +      MappedIterator& operator++() {
   1.142 +        ++base;
   1.143 +        return *this;
   1.144 +      }
   1.145 +
   1.146 +      MappedIterator& operator++(int) {
   1.147 +        MappedIterator tmp(*this);
   1.148 +        ++base;
   1.149 +        return tmp;
   1.150 +      }
   1.151 +
   1.152 +      bool operator==(const MappedIterator& it) const {
   1.153 +        return base == it.base;
   1.154 +      }
   1.155 +
   1.156 +      bool operator!=(const MappedIterator& it) const {
   1.157 +        return base != it.base;
   1.158 +      }
   1.159 +
   1.160 +    private:
   1.161 +      Base base;
   1.162 +      const LpSolverBase& lp;
   1.163 +    };
   1.164 +
   1.165    protected:
   1.166  
   1.167 +    /// STL compatible iterator for lp col
   1.168 +    typedef MappedIterator<Expr::const_iterator> LpRowIterator;
   1.169 +    /// STL compatible iterator for lp row
   1.170 +    typedef MappedIterator<DualExpr::const_iterator> LpColIterator;
   1.171 +
   1.172      //Abstract virtual functions
   1.173      virtual LpSolverBase &_newLp() = 0;
   1.174      virtual LpSolverBase &_copyLp(){
   1.175 -      ///\todo This should be implemented here, too,  when we have problem retrieving routines. It can be overriden.
   1.176 +      ///\todo This should be implemented here, too, when we have
   1.177 +      ///problem retrieving routines. It can be overriden.
   1.178  
   1.179        //Starting:
   1.180        LpSolverBase & newlp(_newLp());
   1.181 @@ -613,16 +698,10 @@
   1.182      virtual int _addRow() = 0; 
   1.183      virtual void _eraseCol(int col) = 0;
   1.184      virtual void _eraseRow(int row) = 0;
   1.185 -    virtual void _getColName(int col,       std::string & name) = 0;
   1.186 +    virtual void _getColName(int col, std::string & name) = 0;
   1.187      virtual void _setColName(int col, const std::string & name) = 0;
   1.188 -    virtual void _setRowCoeffs(int i, 
   1.189 -			       int length,
   1.190 -                               int  const * indices, 
   1.191 -                               Value  const * values ) = 0;
   1.192 -    virtual void _setColCoeffs(int i, 
   1.193 -			       int length,
   1.194 -                               int  const * indices, 
   1.195 -                               Value  const * values ) = 0;
   1.196 +    virtual void _setRowCoeffs(int i, LpRowIterator b, LpRowIterator e) = 0;
   1.197 +    virtual void _setColCoeffs(int i, LpColIterator b, LpColIterator e) = 0;
   1.198      virtual void _setCoeff(int row, int col, Value value) = 0;
   1.199      virtual void _setColLowerBound(int i, Value value) = 0;
   1.200      virtual void _setColUpperBound(int i, Value value) = 0;
   1.201 @@ -631,9 +710,7 @@
   1.202      virtual void _setRowBounds(int i, Value lower, Value upper) = 0;
   1.203      virtual void _setObjCoeff(int i, Value obj_coef) = 0;
   1.204      virtual void _clearObj()=0;
   1.205 -//     virtual void _setObj(int length,
   1.206 -//                          int  const * indices, 
   1.207 -//                          Value  const * values ) = 0;
   1.208 +
   1.209      virtual SolveExitStatus _solve() = 0;
   1.210      virtual Value _getPrimal(int i) = 0;
   1.211      virtual Value _getDual(int i) = 0;
   1.212 @@ -652,10 +729,7 @@
   1.213      
   1.214      //Constant component of the objective function
   1.215      Value obj_const_comp;
   1.216 -    
   1.217 -
   1.218 -
   1.219 -    
   1.220 +        
   1.221    public:
   1.222  
   1.223      ///\e
   1.224 @@ -744,17 +818,9 @@
   1.225      ///\param e is a dual linear expression (see \ref DualExpr)
   1.226      ///a better one.
   1.227      void col(Col c,const DualExpr &e) {
   1.228 -      std::vector<int> indices;
   1.229 -      std::vector<Value> values;
   1.230 -      indices.push_back(0);
   1.231 -      values.push_back(0);
   1.232 -      for(DualExpr::const_iterator i=e.begin(); i!=e.end(); ++i)
   1.233 -	if((*i).second!=0) {
   1.234 -	  indices.push_back(rows.floatingId((*i).first.id));
   1.235 -	  values.push_back((*i).second);
   1.236 -	}
   1.237 -      _setColCoeffs(cols.floatingId(c.id),indices.size()-1,
   1.238 -		    &indices[0],&values[0]);
   1.239 +      e.simplify();
   1.240 +      _setColCoeffs(_lpId(c), LpColIterator(e.begin(), *this), 
   1.241 +                    LpColIterator(e.end(), *this));
   1.242      }
   1.243  
   1.244      ///Add a new column to the LP
   1.245 @@ -849,20 +915,12 @@
   1.246      ///\todo Option to control whether a constraint with a single variable is
   1.247      ///added or not.
   1.248      void row(Row r, Value l,const Expr &e, Value u) {
   1.249 -      std::vector<int> indices;
   1.250 -      std::vector<Value> values;
   1.251 -      indices.push_back(0);
   1.252 -      values.push_back(0);
   1.253 -      for(Expr::const_iterator i=e.begin(); i!=e.end(); ++i)
   1.254 -	if((*i).second!=0) { ///\bug EPSILON would be necessary here!!!
   1.255 -	  indices.push_back(cols.floatingId((*i).first.id));
   1.256 -	  values.push_back((*i).second);
   1.257 -	}
   1.258 -      _setRowCoeffs(rows.floatingId(r.id),indices.size()-1,
   1.259 -		    &indices[0],&values[0]);
   1.260 -//       _setRowLowerBound(rows.floatingId(r.id),l-e.constComp());
   1.261 -//       _setRowUpperBound(rows.floatingId(r.id),u-e.constComp());
   1.262 -       _setRowBounds(rows.floatingId(r.id),l-e.constComp(),u-e.constComp());
   1.263 +      e.simplify();
   1.264 +      _setRowCoeffs(_lpId(r), LpRowIterator(e.begin(), *this),
   1.265 +                    LpRowIterator(e.end(), *this));
   1.266 +//       _setRowLowerBound(_lpId(r),l-e.constComp());
   1.267 +//       _setRowUpperBound(_lpId(r),u-e.constComp());
   1.268 +       _setRowBounds(_lpId(r),l-e.constComp(),u-e.constComp());
   1.269      }
   1.270  
   1.271      ///Set a row (i.e a constraint) of the LP
   1.272 @@ -870,10 +928,8 @@
   1.273      ///\param r is the row to be modified
   1.274      ///\param c is a linear expression (see \ref Constr)
   1.275      void row(Row r, const Constr &c) {
   1.276 -      row(r,
   1.277 -	     c.lowerBounded()?c.lowerBound():-INF,
   1.278 -	     c.expr(),
   1.279 -	     c.upperBounded()?c.upperBound():INF);
   1.280 +      row(r, c.lowerBounded()?c.lowerBound():-INF,
   1.281 +          c.expr(), c.upperBounded()?c.upperBound():INF);
   1.282      }
   1.283  
   1.284      ///Add a new row (i.e a new constraint) to the LP
   1.285 @@ -904,7 +960,7 @@
   1.286      ///\param c is the coloumn to be deleted
   1.287      ///\todo Please check this
   1.288      void eraseCol(Col c) {
   1.289 -      _eraseCol(cols.floatingId(c.id));
   1.290 +      _eraseCol(_lpId(c));
   1.291        cols.erase(c.id);
   1.292      }
   1.293      ///Erase a  row (i.e a constraint) from the LP
   1.294 @@ -912,7 +968,7 @@
   1.295      ///\param r is the row to be deleted
   1.296      ///\todo Please check this
   1.297      void eraseRow(Row r) {
   1.298 -      _eraseRow(rows.floatingId(r.id));
   1.299 +      _eraseRow(_lpId(r));
   1.300        rows.erase(r.id);
   1.301      }
   1.302  
   1.303 @@ -922,7 +978,7 @@
   1.304      ///\return The name of the colunm
   1.305      std::string colName(Col c){
   1.306        std::string name;
   1.307 -      _getColName(cols.floatingId(c.id), name);
   1.308 +      _getColName(_lpId(c), name);
   1.309        return name;
   1.310      }
   1.311      
   1.312 @@ -930,8 +986,8 @@
   1.313      
   1.314      ///\param c is the coresponding coloumn 
   1.315      ///\param name The name to be given
   1.316 -    void colName(Col c, const std::string & name){
   1.317 -      _setColName(cols.floatingId(c.id), name);
   1.318 +    void colName(Col c, const std::string& name){
   1.319 +      _setColName(_lpId(c), name);
   1.320      }
   1.321      
   1.322      /// Set an element of the coefficient matrix of the LP
   1.323 @@ -941,7 +997,7 @@
   1.324      ///\param val is the new value of the coefficient
   1.325  
   1.326      void coeff(Row r, Col c, Value val){
   1.327 -      _setCoeff(rows.floatingId(r.id),cols.floatingId(c.id), val);
   1.328 +      _setCoeff(_lpId(r),_lpId(c), val);
   1.329      }
   1.330  
   1.331      /// Set the lower bound of a column (i.e a variable)
   1.332 @@ -950,7 +1006,7 @@
   1.333      /// extended number of type Value, i.e. a finite number of type 
   1.334      /// Value or -\ref INF.
   1.335      void colLowerBound(Col c, Value value) {
   1.336 -      _setColLowerBound(cols.floatingId(c.id),value);
   1.337 +      _setColLowerBound(_lpId(c),value);
   1.338      }
   1.339      
   1.340      ///\brief Set the lower bound of  several columns
   1.341 @@ -996,7 +1052,7 @@
   1.342      /// extended number of type Value, i.e. a finite number of type 
   1.343      /// Value or \ref INF.
   1.344      void colUpperBound(Col c, Value value) {
   1.345 -      _setColUpperBound(cols.floatingId(c.id),value);
   1.346 +      _setColUpperBound(_lpId(c),value);
   1.347      };
   1.348  
   1.349      ///\brief Set the lower bound of  several columns
   1.350 @@ -1043,8 +1099,8 @@
   1.351      /// extended number of type Value, i.e. a finite number of type 
   1.352      /// Value, -\ref INF or \ref INF.
   1.353      void colBounds(Col c, Value lower, Value upper) {
   1.354 -      _setColLowerBound(cols.floatingId(c.id),lower);
   1.355 -      _setColUpperBound(cols.floatingId(c.id),upper);
   1.356 +      _setColLowerBound(_lpId(c),lower);
   1.357 +      _setColUpperBound(_lpId(c),upper);
   1.358      }
   1.359      
   1.360      ///\brief Set the lower and the upper bound of several columns
   1.361 @@ -1091,7 +1147,7 @@
   1.362  //     /// extended number of type Value, i.e. a finite number of type 
   1.363  //     /// Value or -\ref INF.
   1.364  //     void rowLowerBound(Row r, Value value) {
   1.365 -//       _setRowLowerBound(rows.floatingId(r.id),value);
   1.366 +//       _setRowLowerBound(_lpId(r),value);
   1.367  //     };
   1.368  //     /// Set the upper bound of a row (i.e a constraint)
   1.369  
   1.370 @@ -1099,7 +1155,7 @@
   1.371  //     /// extended number of type Value, i.e. a finite number of type 
   1.372  //     /// Value or \ref INF.
   1.373  //     void rowUpperBound(Row r, Value value) {
   1.374 -//       _setRowUpperBound(rows.floatingId(r.id),value);
   1.375 +//       _setRowUpperBound(_lpId(r),value);
   1.376  //     };
   1.377  
   1.378      /// Set the lower and the upper bounds of a row (i.e a constraint)
   1.379 @@ -1109,12 +1165,12 @@
   1.380      /// extended number of type Value, i.e. a finite number of type 
   1.381      /// Value, -\ref INF or \ref INF.
   1.382      void rowBounds(Row c, Value lower, Value upper) {
   1.383 -      _setRowBounds(rows.floatingId(c.id),lower, upper);
   1.384 -      // _setRowUpperBound(rows.floatingId(c.id),upper);
   1.385 +      _setRowBounds(_lpId(c),lower, upper);
   1.386 +      // _setRowUpperBound(_lpId(c),upper);
   1.387      }
   1.388      
   1.389      ///Set an element of the objective function
   1.390 -    void objCoeff(Col c, Value v) {_setObjCoeff(cols.floatingId(c.id),v); };
   1.391 +    void objCoeff(Col c, Value v) {_setObjCoeff(_lpId(c),v); };
   1.392      ///Set the objective function
   1.393      
   1.394      ///\param e is a linear expression of type \ref Expr.
   1.395 @@ -1170,13 +1226,13 @@
   1.396      }
   1.397  
   1.398      ///\e
   1.399 -    Value primal(Col c) { return _getPrimal(cols.floatingId(c.id)); }
   1.400 +    Value primal(Col c) { return _getPrimal(_lpId(c)); }
   1.401  
   1.402      ///\e
   1.403 -    Value dual(Row r) { return _getDual(rows.floatingId(r.id)); }
   1.404 +    Value dual(Row r) { return _getDual(_lpId(r)); }
   1.405  
   1.406      ///\e
   1.407 -    bool isBasicCol(Col c) { return _isBasicCol(cols.floatingId(c.id)); }
   1.408 +    bool isBasicCol(Col c) { return _isBasicCol(_lpId(c)); }
   1.409  
   1.410      ///\e
   1.411  
   1.412 @@ -1213,14 +1269,14 @@
   1.413      ///
   1.414      ///Sets the type of the given coloumn to the given type.
   1.415      void colType(Col c, ColTypes col_type) {
   1.416 -      _colType(cols.floatingId(c.id),col_type);
   1.417 +      _colType(_lpId(c),col_type);
   1.418      }
   1.419  
   1.420      ///Gives back the type of the column.
   1.421      ///
   1.422      ///Gives back the type of the column.
   1.423      ColTypes colType(Col c){
   1.424 -      return _colType(cols.floatingId(c.id));
   1.425 +      return _colType(_lpId(c));
   1.426      }
   1.427  
   1.428      ///Sets the type of the given Col to integer or remove that property.
   1.429 @@ -1440,7 +1496,7 @@
   1.430    ///\relates LpSolverBase::DualExpr
   1.431    ///
   1.432    inline LpSolverBase::DualExpr operator+(const LpSolverBase::DualExpr &a,
   1.433 -				      const LpSolverBase::DualExpr &b) 
   1.434 +                                          const LpSolverBase::DualExpr &b) 
   1.435    {
   1.436      LpSolverBase::DualExpr tmp(a);
   1.437      tmp+=b;
   1.438 @@ -1451,7 +1507,7 @@
   1.439    ///\relates LpSolverBase::DualExpr
   1.440    ///
   1.441    inline LpSolverBase::DualExpr operator-(const LpSolverBase::DualExpr &a,
   1.442 -				      const LpSolverBase::DualExpr &b) 
   1.443 +                                          const LpSolverBase::DualExpr &b) 
   1.444    {
   1.445      LpSolverBase::DualExpr tmp(a);
   1.446      tmp-=b;
   1.447 @@ -1462,7 +1518,7 @@
   1.448    ///\relates LpSolverBase::DualExpr
   1.449    ///
   1.450    inline LpSolverBase::DualExpr operator*(const LpSolverBase::DualExpr &a,
   1.451 -				      const LpSolverBase::Value &b) 
   1.452 +                                          const LpSolverBase::Value &b) 
   1.453    {
   1.454      LpSolverBase::DualExpr tmp(a);
   1.455      tmp*=b;
   1.456 @@ -1474,7 +1530,7 @@
   1.457    ///\relates LpSolverBase::DualExpr
   1.458    ///
   1.459    inline LpSolverBase::DualExpr operator*(const LpSolverBase::Value &a,
   1.460 -				      const LpSolverBase::DualExpr &b) 
   1.461 +                                          const LpSolverBase::DualExpr &b) 
   1.462    {
   1.463      LpSolverBase::DualExpr tmp(b);
   1.464      tmp*=a;
   1.465 @@ -1485,7 +1541,7 @@
   1.466    ///\relates LpSolverBase::DualExpr
   1.467    ///
   1.468    inline LpSolverBase::DualExpr operator/(const LpSolverBase::DualExpr &a,
   1.469 -				      const LpSolverBase::Value &b) 
   1.470 +                                          const LpSolverBase::Value &b) 
   1.471    {
   1.472      LpSolverBase::DualExpr tmp(a);
   1.473      tmp/=b;