408 using namespace std; |
408 using namespace std; |
409 return finite(_ub); |
409 return finite(_ub); |
410 } |
410 } |
411 }; |
411 }; |
412 |
412 |
|
413 ///Linear expression of rows |
|
414 |
|
415 ///This data structure represents a column of the matrix, |
|
416 ///thas is it strores a linear expression of the dual variables |
|
417 ///(\ref Row "Row"s). |
|
418 /// |
|
419 ///There are several ways to access and modify the contents of this |
|
420 ///container. |
|
421 ///- Its it fully compatible with \c std::map<Row,double>, so for expamle |
|
422 ///if \c e is an DualExpr and \c v |
|
423 ///and \c w are of type \ref Row, then you can |
|
424 ///read and modify the coefficients like |
|
425 ///these. |
|
426 ///\code |
|
427 ///e[v]=5; |
|
428 ///e[v]+=12; |
|
429 ///e.erase(v); |
|
430 ///\endcode |
|
431 ///or you can also iterate through its elements. |
|
432 ///\code |
|
433 ///double s=0; |
|
434 ///for(LpSolverBase::DualExpr::iterator i=e.begin();i!=e.end();++i) |
|
435 /// s+=i->second; |
|
436 ///\endcode |
|
437 ///(This code computes the sum of all coefficients). |
|
438 ///- Numbers (<tt>double</tt>'s) |
|
439 ///and variables (\ref Row "Row"s) directly convert to an |
|
440 ///\ref DualExpr and the usual linear operations are defined so |
|
441 ///\code |
|
442 ///v+w |
|
443 ///2*v-3.12*(v-w/2) |
|
444 ///v*2.1+(3*v+(v*12+w)*3)/2 |
|
445 ///\endcode |
|
446 ///are valid \ref DualExpr "DualExpr"essions. |
|
447 ///The usual assignment operations are also defined. |
|
448 ///\code |
|
449 ///e=v+w; |
|
450 ///e+=2*v-3.12*(v-w/2); |
|
451 ///e*=3.4; |
|
452 ///e/=5; |
|
453 ///\endcode |
|
454 /// |
|
455 ///\sa Expr |
|
456 /// |
|
457 class DualExpr : public std::map<Row,Value> |
|
458 { |
|
459 public: |
|
460 typedef LpSolverBase::Row Key; |
|
461 typedef LpSolverBase::Value Value; |
|
462 |
|
463 protected: |
|
464 typedef std::map<Row,Value> Base; |
|
465 |
|
466 public: |
|
467 typedef True IsLinExpression; |
|
468 ///\e |
|
469 DualExpr() : Base() { } |
|
470 ///\e |
|
471 DualExpr(const Key &v) { |
|
472 Base::insert(std::make_pair(v, 1)); |
|
473 } |
|
474 ///\e |
|
475 DualExpr(const Value &v) {} |
|
476 ///\e |
|
477 void set(const Key &v,const Value &c) { |
|
478 Base::insert(std::make_pair(v, c)); |
|
479 } |
|
480 |
|
481 ///Removes the components with zero coefficient. |
|
482 void simplify() { |
|
483 for (Base::iterator i=Base::begin(); i!=Base::end();) { |
|
484 Base::iterator j=i; |
|
485 ++j; |
|
486 if ((*i).second==0) Base::erase(i); |
|
487 j=i; |
|
488 } |
|
489 } |
|
490 |
|
491 ///Sets all coefficients to 0. |
|
492 void clear() { |
|
493 Base::clear(); |
|
494 } |
|
495 |
|
496 ///\e |
|
497 DualExpr &operator+=(const DualExpr &e) { |
|
498 for (Base::const_iterator j=e.begin(); j!=e.end(); ++j) |
|
499 (*this)[j->first]+=j->second; |
|
500 ///\todo it might be speeded up using "hints" |
|
501 return *this; |
|
502 } |
|
503 ///\e |
|
504 DualExpr &operator-=(const DualExpr &e) { |
|
505 for (Base::const_iterator j=e.begin(); j!=e.end(); ++j) |
|
506 (*this)[j->first]-=j->second; |
|
507 return *this; |
|
508 } |
|
509 ///\e |
|
510 DualExpr &operator*=(const Value &c) { |
|
511 for (Base::iterator j=Base::begin(); j!=Base::end(); ++j) |
|
512 j->second*=c; |
|
513 return *this; |
|
514 } |
|
515 ///\e |
|
516 DualExpr &operator/=(const Value &c) { |
|
517 for (Base::iterator j=Base::begin(); j!=Base::end(); ++j) |
|
518 j->second/=c; |
|
519 return *this; |
|
520 } |
|
521 }; |
|
522 |
413 |
523 |
414 protected: |
524 protected: |
415 _FixId rows; |
525 _FixId rows; |
416 _FixId cols; |
526 _FixId cols; |
417 |
527 |
545 } |
655 } |
546 return s; |
656 return s; |
547 } |
657 } |
548 #endif |
658 #endif |
549 |
659 |
550 ///Add a new empty row (i.e a new constaint) to the LP |
660 ///Set a column (i.e a dual constraint) of the LP |
551 |
661 |
552 ///This function adds a new empty row (i.e a new constaint) to the LP. |
662 ///\param c is the column to be modified |
|
663 ///\param e is a dual linear expression (see \ref DualExpr) |
|
664 ///\bug This is a temportary function. The interface will change to |
|
665 ///a better one. |
|
666 void setCol(Col c,const DualExpr &e) { |
|
667 std::vector<int> indices; |
|
668 std::vector<Value> values; |
|
669 indices.push_back(0); |
|
670 values.push_back(0); |
|
671 for(DualExpr::const_iterator i=e.begin(); i!=e.end(); ++i) |
|
672 if((*i).second!=0) { ///\bug EPSILON would be necessary here!!! |
|
673 indices.push_back(cols.floatingId((*i).first.id)); |
|
674 values.push_back((*i).second); |
|
675 } |
|
676 _setColCoeffs(cols.floatingId(c.id),indices.size()-1, |
|
677 &indices[0],&values[0]); |
|
678 } |
|
679 |
|
680 ///Add a new column to the LP |
|
681 |
|
682 ///\param e is a dual linear expression (see \ref DualExpr) |
|
683 ///\param obj is the corresponding component of the objective |
|
684 ///function. It is 0 by default. |
|
685 ///\return The created column. |
|
686 ///\bug This is a temportary function. The interface will change to |
|
687 ///a better one. |
|
688 Col addCol(Value l,const DualExpr &e, Value obj=0) { |
|
689 Col c=addCol(); |
|
690 setCol(c,e); |
|
691 objCoeff(c,0); |
|
692 return c; |
|
693 } |
|
694 |
|
695 ///Add a new empty row (i.e a new constraint) to the LP |
|
696 |
|
697 ///This function adds a new empty row (i.e a new constraint) to the LP. |
553 ///\return The created row |
698 ///\return The created row |
554 Row addRow() { Row r; r.id=rows.insert(_addRow()); return r;} |
699 Row addRow() { Row r; r.id=rows.insert(_addRow()); return r;} |
555 |
700 |
556 ///Set a row (i.e a constaint) of the LP |
701 ///\brief Adds several new row |
|
702 ///(i.e a variables) at once |
|
703 /// |
|
704 ///This magic function takes a container as its argument |
|
705 ///and fills its elements |
|
706 ///with new row (i.e. variables) |
|
707 ///\param t can be |
|
708 ///- a standard STL compatible iterable container with |
|
709 ///\ref Row as its \c values_type |
|
710 ///like |
|
711 ///\code |
|
712 ///std::vector<LpSolverBase::Row> |
|
713 ///std::list<LpSolverBase::Row> |
|
714 ///\endcode |
|
715 ///- a standard STL compatible iterable container with |
|
716 ///\ref Row as its \c mapped_type |
|
717 ///like |
|
718 ///\code |
|
719 ///std::map<AnyType,LpSolverBase::Row> |
|
720 ///\endcode |
|
721 ///- an iterable lemon \ref concept::WriteMap "write map" like |
|
722 ///\code |
|
723 ///ListGraph::NodeMap<LpSolverBase::Row> |
|
724 ///ListGraph::EdgeMap<LpSolverBase::Row> |
|
725 ///\endcode |
|
726 ///\return The number of rows created. |
|
727 #ifdef DOXYGEN |
|
728 template<class T> |
|
729 int addRowSet(T &t) { return 0;} |
|
730 #else |
|
731 template<class T> |
|
732 typename enable_if<typename T::value_type::LpSolverRow,int>::type |
|
733 addRowSet(T &t,dummy<0> = 0) { |
|
734 int s=0; |
|
735 for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addRow();s++;} |
|
736 return s; |
|
737 } |
|
738 template<class T> |
|
739 typename enable_if<typename T::value_type::second_type::LpSolverRow, |
|
740 int>::type |
|
741 addRowSet(T &t,dummy<1> = 1) { |
|
742 int s=0; |
|
743 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
744 i->second=addRow(); |
|
745 s++; |
|
746 } |
|
747 return s; |
|
748 } |
|
749 template<class T> |
|
750 typename enable_if<typename T::ValueSet::value_type::LpSolverRow, |
|
751 int>::type |
|
752 addRowSet(T &t,dummy<2> = 2) { |
|
753 ///\bug <tt>return addRowSet(t.valueSet());</tt> should also work. |
|
754 int s=0; |
|
755 for(typename T::ValueSet::iterator i=t.valueSet().begin(); |
|
756 i!=t.valueSet().end(); |
|
757 ++i) |
|
758 { |
|
759 *i=addRow(); |
|
760 s++; |
|
761 } |
|
762 return s; |
|
763 } |
|
764 #endif |
|
765 |
|
766 ///Set a row (i.e a constraint) of the LP |
557 |
767 |
558 ///\param r is the row to be modified |
768 ///\param r is the row to be modified |
559 ///\param l is lower bound (-\ref INF means no bound) |
769 ///\param l is lower bound (-\ref INF means no bound) |
560 ///\param e is a linear expression (see \ref Expr) |
770 ///\param e is a linear expression (see \ref Expr) |
561 ///\param u is the upper bound (\ref INF means no bound) |
771 ///\param u is the upper bound (\ref INF means no bound) |
578 // _setRowLowerBound(rows.floatingId(r.id),l-e.constComp()); |
788 // _setRowLowerBound(rows.floatingId(r.id),l-e.constComp()); |
579 // _setRowUpperBound(rows.floatingId(r.id),u-e.constComp()); |
789 // _setRowUpperBound(rows.floatingId(r.id),u-e.constComp()); |
580 _setRowBounds(rows.floatingId(r.id),l-e.constComp(),u-e.constComp()); |
790 _setRowBounds(rows.floatingId(r.id),l-e.constComp(),u-e.constComp()); |
581 } |
791 } |
582 |
792 |
583 ///Set a row (i.e a constaint) of the LP |
793 ///Set a row (i.e a constraint) of the LP |
584 |
794 |
585 ///\param r is the row to be modified |
795 ///\param r is the row to be modified |
586 ///\param c is a linear expression (see \ref Constr) |
796 ///\param c is a linear expression (see \ref Constr) |
587 void setRow(Row r, const Constr &c) { |
797 void setRow(Row r, const Constr &c) { |
588 setRow(r, |
798 setRow(r, |
589 c.lowerBounded()?c.lowerBound():-INF, |
799 c.lowerBounded()?c.lowerBound():-INF, |
590 c.expr(), |
800 c.expr(), |
591 c.upperBounded()?c.upperBound():INF); |
801 c.upperBounded()?c.upperBound():INF); |
592 } |
802 } |
593 |
803 |
594 ///Add a new row (i.e a new constaint) to the LP |
804 ///Add a new row (i.e a new constraint) to the LP |
595 |
805 |
596 ///\param l is the lower bound (-\ref INF means no bound) |
806 ///\param l is the lower bound (-\ref INF means no bound) |
597 ///\param e is a linear expression (see \ref Expr) |
807 ///\param e is a linear expression (see \ref Expr) |
598 ///\param u is the upper bound (\ref INF means no bound) |
808 ///\param u is the upper bound (\ref INF means no bound) |
599 ///\return The created row. |
809 ///\return The created row. |
915 if(!isnan(tmp.lowerBound())) throw LogicError(); |
1125 if(!isnan(tmp.lowerBound())) throw LogicError(); |
916 else tmp.lowerBound()=n; |
1126 else tmp.lowerBound()=n; |
917 return tmp; |
1127 return tmp; |
918 } |
1128 } |
919 |
1129 |
|
1130 ///\e |
|
1131 |
|
1132 ///\relates LpSolverBase::DualExpr |
|
1133 /// |
|
1134 inline LpSolverBase::DualExpr operator+(const LpSolverBase::DualExpr &a, |
|
1135 const LpSolverBase::DualExpr &b) |
|
1136 { |
|
1137 LpSolverBase::DualExpr tmp(a); |
|
1138 tmp+=b; ///\todo Doesn't STL have some special 'merge' algorithm? |
|
1139 return tmp; |
|
1140 } |
|
1141 ///\e |
|
1142 |
|
1143 ///\relates LpSolverBase::DualExpr |
|
1144 /// |
|
1145 inline LpSolverBase::DualExpr operator-(const LpSolverBase::DualExpr &a, |
|
1146 const LpSolverBase::DualExpr &b) |
|
1147 { |
|
1148 LpSolverBase::DualExpr tmp(a); |
|
1149 tmp-=b; ///\todo Doesn't STL have some special 'merge' algorithm? |
|
1150 return tmp; |
|
1151 } |
|
1152 ///\e |
|
1153 |
|
1154 ///\relates LpSolverBase::DualExpr |
|
1155 /// |
|
1156 inline LpSolverBase::DualExpr operator*(const LpSolverBase::DualExpr &a, |
|
1157 const LpSolverBase::Value &b) |
|
1158 { |
|
1159 LpSolverBase::DualExpr tmp(a); |
|
1160 tmp*=b; ///\todo Doesn't STL have some special 'merge' algorithm? |
|
1161 return tmp; |
|
1162 } |
|
1163 |
|
1164 ///\e |
|
1165 |
|
1166 ///\relates LpSolverBase::DualExpr |
|
1167 /// |
|
1168 inline LpSolverBase::DualExpr operator*(const LpSolverBase::Value &a, |
|
1169 const LpSolverBase::DualExpr &b) |
|
1170 { |
|
1171 LpSolverBase::DualExpr tmp(b); |
|
1172 tmp*=a; ///\todo Doesn't STL have some special 'merge' algorithm? |
|
1173 return tmp; |
|
1174 } |
|
1175 ///\e |
|
1176 |
|
1177 ///\relates LpSolverBase::DualExpr |
|
1178 /// |
|
1179 inline LpSolverBase::DualExpr operator/(const LpSolverBase::DualExpr &a, |
|
1180 const LpSolverBase::Value &b) |
|
1181 { |
|
1182 LpSolverBase::DualExpr tmp(a); |
|
1183 tmp/=b; ///\todo Doesn't STL have some special 'merge' algorithm? |
|
1184 return tmp; |
|
1185 } |
|
1186 |
920 |
1187 |
921 } //namespace lemon |
1188 } //namespace lemon |
922 |
1189 |
923 #endif //LEMON_LP_BASE_H |
1190 #endif //LEMON_LP_BASE_H |