... |
... |
@@ -504,197 +504,197 @@
|
504 |
504 |
operator Col() const {
|
505 |
505 |
return colFromId(_it->first);
|
506 |
506 |
}
|
507 |
507 |
|
508 |
508 |
/// Returns the coefficient of the term
|
509 |
509 |
const Value& operator*() const { return _it->second; }
|
510 |
510 |
|
511 |
511 |
/// Next term
|
512 |
512 |
|
513 |
513 |
/// Assign the iterator to the next term.
|
514 |
514 |
///
|
515 |
515 |
ConstCoeffIt& operator++() { ++_it; return *this; }
|
516 |
516 |
|
517 |
517 |
/// Equality operator
|
518 |
518 |
bool operator==(Invalid) const { return _it == _end; }
|
519 |
519 |
/// Inequality operator
|
520 |
520 |
bool operator!=(Invalid) const { return _it != _end; }
|
521 |
521 |
};
|
522 |
522 |
|
523 |
523 |
};
|
524 |
524 |
|
525 |
525 |
///Linear constraint
|
526 |
526 |
|
527 |
527 |
///This data stucture represents a linear constraint in the LP.
|
528 |
528 |
///Basically it is a linear expression with a lower or an upper bound
|
529 |
529 |
///(or both). These parts of the constraint can be obtained by the member
|
530 |
530 |
///functions \ref expr(), \ref lowerBound() and \ref upperBound(),
|
531 |
531 |
///respectively.
|
532 |
532 |
///There are two ways to construct a constraint.
|
533 |
533 |
///- You can set the linear expression and the bounds directly
|
534 |
534 |
/// by the functions above.
|
535 |
535 |
///- The operators <tt>\<=</tt>, <tt>==</tt> and <tt>\>=</tt>
|
536 |
536 |
/// are defined between expressions, or even between constraints whenever
|
537 |
537 |
/// it makes sense. Therefore if \c e and \c f are linear expressions and
|
538 |
538 |
/// \c s and \c t are numbers, then the followings are valid expressions
|
539 |
539 |
/// and thus they can be used directly e.g. in \ref addRow() whenever
|
540 |
540 |
/// it makes sense.
|
541 |
541 |
///\code
|
542 |
542 |
/// e<=s
|
543 |
543 |
/// e<=f
|
544 |
544 |
/// e==f
|
545 |
545 |
/// s<=e<=t
|
546 |
546 |
/// e>=t
|
547 |
547 |
///\endcode
|
548 |
548 |
///\warning The validity of a constraint is checked only at run
|
549 |
549 |
///time, so e.g. \ref addRow(<tt>x[1]\<=x[2]<=5</tt>) will
|
550 |
550 |
///compile, but will fail an assertion.
|
551 |
551 |
class Constr
|
552 |
552 |
{
|
553 |
553 |
public:
|
554 |
554 |
typedef LpBase::Expr Expr;
|
555 |
555 |
typedef Expr::Key Key;
|
556 |
556 |
typedef Expr::Value Value;
|
557 |
557 |
|
558 |
558 |
protected:
|
559 |
559 |
Expr _expr;
|
560 |
560 |
Value _lb,_ub;
|
561 |
561 |
public:
|
562 |
562 |
///\e
|
563 |
563 |
Constr() : _expr(), _lb(NaN), _ub(NaN) {}
|
564 |
564 |
///\e
|
565 |
565 |
Constr(Value lb, const Expr &e, Value ub) :
|
566 |
566 |
_expr(e), _lb(lb), _ub(ub) {}
|
567 |
567 |
Constr(const Expr &e) :
|
568 |
568 |
_expr(e), _lb(NaN), _ub(NaN) {}
|
569 |
569 |
///\e
|
570 |
570 |
void clear()
|
571 |
571 |
{
|
572 |
572 |
_expr.clear();
|
573 |
573 |
_lb=_ub=NaN;
|
574 |
574 |
}
|
575 |
575 |
|
576 |
576 |
///Reference to the linear expression
|
577 |
577 |
Expr &expr() { return _expr; }
|
578 |
578 |
///Cont reference to the linear expression
|
579 |
579 |
const Expr &expr() const { return _expr; }
|
580 |
580 |
///Reference to the lower bound.
|
581 |
581 |
|
582 |
582 |
///\return
|
583 |
583 |
///- \ref INF "INF": the constraint is lower unbounded.
|
584 |
584 |
///- \ref NaN "NaN": lower bound has not been set.
|
585 |
585 |
///- finite number: the lower bound
|
586 |
586 |
Value &lowerBound() { return _lb; }
|
587 |
587 |
///The const version of \ref lowerBound()
|
588 |
588 |
const Value &lowerBound() const { return _lb; }
|
589 |
589 |
///Reference to the upper bound.
|
590 |
590 |
|
591 |
591 |
///\return
|
592 |
592 |
///- \ref INF "INF": the constraint is upper unbounded.
|
593 |
593 |
///- \ref NaN "NaN": upper bound has not been set.
|
594 |
594 |
///- finite number: the upper bound
|
595 |
595 |
Value &upperBound() { return _ub; }
|
596 |
596 |
///The const version of \ref upperBound()
|
597 |
597 |
const Value &upperBound() const { return _ub; }
|
598 |
598 |
///Is the constraint lower bounded?
|
599 |
599 |
bool lowerBounded() const {
|
600 |
|
return _lb != -INF && !isnan(_lb);
|
|
600 |
return _lb != -INF && !isNaN(_lb);
|
601 |
601 |
}
|
602 |
602 |
///Is the constraint upper bounded?
|
603 |
603 |
bool upperBounded() const {
|
604 |
|
return _ub != INF && !isnan(_ub);
|
|
604 |
return _ub != INF && !isNaN(_ub);
|
605 |
605 |
}
|
606 |
606 |
|
607 |
607 |
};
|
608 |
608 |
|
609 |
609 |
///Linear expression of rows
|
610 |
610 |
|
611 |
611 |
///This data structure represents a column of the matrix,
|
612 |
612 |
///thas is it strores a linear expression of the dual variables
|
613 |
613 |
///(\ref Row "Row"s).
|
614 |
614 |
///
|
615 |
615 |
///There are several ways to access and modify the contents of this
|
616 |
616 |
///container.
|
617 |
617 |
///\code
|
618 |
618 |
///e[v]=5;
|
619 |
619 |
///e[v]+=12;
|
620 |
620 |
///e.erase(v);
|
621 |
621 |
///\endcode
|
622 |
622 |
///or you can also iterate through its elements.
|
623 |
623 |
///\code
|
624 |
624 |
///double s=0;
|
625 |
625 |
///for(LpBase::DualExpr::ConstCoeffIt i(e);i!=INVALID;++i)
|
626 |
626 |
/// s+=*i;
|
627 |
627 |
///\endcode
|
628 |
628 |
///(This code computes the sum of all coefficients).
|
629 |
629 |
///- Numbers (<tt>double</tt>'s)
|
630 |
630 |
///and variables (\ref Row "Row"s) directly convert to an
|
631 |
631 |
///\ref DualExpr and the usual linear operations are defined, so
|
632 |
632 |
///\code
|
633 |
633 |
///v+w
|
634 |
634 |
///2*v-3.12*(v-w/2)
|
635 |
635 |
///v*2.1+(3*v+(v*12+w)*3)/2
|
636 |
636 |
///\endcode
|
637 |
637 |
///are valid \ref DualExpr dual expressions.
|
638 |
638 |
///The usual assignment operations are also defined.
|
639 |
639 |
///\code
|
640 |
640 |
///e=v+w;
|
641 |
641 |
///e+=2*v-3.12*(v-w/2);
|
642 |
642 |
///e*=3.4;
|
643 |
643 |
///e/=5;
|
644 |
644 |
///\endcode
|
645 |
645 |
///
|
646 |
646 |
///\sa Expr
|
647 |
647 |
class DualExpr {
|
648 |
648 |
friend class LpBase;
|
649 |
649 |
public:
|
650 |
650 |
/// The key type of the expression
|
651 |
651 |
typedef LpBase::Row Key;
|
652 |
652 |
/// The value type of the expression
|
653 |
653 |
typedef LpBase::Value Value;
|
654 |
654 |
|
655 |
655 |
protected:
|
656 |
656 |
std::map<int, Value> comps;
|
657 |
657 |
|
658 |
658 |
public:
|
659 |
659 |
typedef True SolverExpr;
|
660 |
660 |
/// Default constructor
|
661 |
661 |
|
662 |
662 |
/// Construct an empty expression, the coefficients are
|
663 |
663 |
/// initialized to zero.
|
664 |
664 |
DualExpr() {}
|
665 |
665 |
/// Construct an expression from a row
|
666 |
666 |
|
667 |
667 |
/// Construct an expression, which has a term with \c r dual
|
668 |
668 |
/// variable and 1.0 coefficient.
|
669 |
669 |
DualExpr(const Row &r) {
|
670 |
670 |
typedef std::map<int, Value>::value_type pair_type;
|
671 |
671 |
comps.insert(pair_type(id(r), 1));
|
672 |
672 |
}
|
673 |
673 |
/// Returns the coefficient of the row
|
674 |
674 |
Value operator[](const Row& r) const {
|
675 |
675 |
std::map<int, Value>::const_iterator it = comps.find(id(r));
|
676 |
676 |
if (it != comps.end()) {
|
677 |
677 |
return it->second;
|
678 |
678 |
} else {
|
679 |
679 |
return 0;
|
680 |
680 |
}
|
681 |
681 |
}
|
682 |
682 |
/// Returns the coefficient of the row
|
683 |
683 |
Value& operator[](const Row& r) {
|
684 |
684 |
return comps[id(r)];
|
685 |
685 |
}
|
686 |
686 |
/// Sets the coefficient of the row
|
687 |
687 |
void set(const Row &r, const Value &v) {
|
688 |
688 |
if (v != 0.0) {
|
689 |
689 |
typedef std::map<int, Value>::value_type pair_type;
|
690 |
690 |
comps.insert(pair_type(id(r), v));
|
691 |
691 |
} else {
|
692 |
692 |
comps.erase(id(r));
|
693 |
693 |
}
|
694 |
694 |
}
|
695 |
695 |
/// \brief Removes the coefficients which's absolute value does
|
696 |
696 |
/// not exceed \c epsilon.
|
697 |
697 |
void simplify(Value epsilon = 0.0) {
|
698 |
698 |
std::map<int, Value>::iterator it=comps.begin();
|
699 |
699 |
while (it != comps.end()) {
|
700 |
700 |
std::map<int, Value>::iterator jt=it;
|
... |
... |
@@ -1573,229 +1573,229 @@
|
1573 |
1573 |
inline LpBase::Expr operator*(const LpBase::Value &a, const LpBase::Expr &b) {
|
1574 |
1574 |
LpBase::Expr tmp(b);
|
1575 |
1575 |
tmp*=a;
|
1576 |
1576 |
return tmp;
|
1577 |
1577 |
}
|
1578 |
1578 |
///Divide with constant
|
1579 |
1579 |
|
1580 |
1580 |
///\relates LpBase::Expr
|
1581 |
1581 |
///
|
1582 |
1582 |
inline LpBase::Expr operator/(const LpBase::Expr &a, const LpBase::Value &b) {
|
1583 |
1583 |
LpBase::Expr tmp(a);
|
1584 |
1584 |
tmp/=b;
|
1585 |
1585 |
return tmp;
|
1586 |
1586 |
}
|
1587 |
1587 |
|
1588 |
1588 |
///Create constraint
|
1589 |
1589 |
|
1590 |
1590 |
///\relates LpBase::Constr
|
1591 |
1591 |
///
|
1592 |
1592 |
inline LpBase::Constr operator<=(const LpBase::Expr &e,
|
1593 |
1593 |
const LpBase::Expr &f) {
|
1594 |
1594 |
return LpBase::Constr(0, f - e, LpBase::INF);
|
1595 |
1595 |
}
|
1596 |
1596 |
|
1597 |
1597 |
///Create constraint
|
1598 |
1598 |
|
1599 |
1599 |
///\relates LpBase::Constr
|
1600 |
1600 |
///
|
1601 |
1601 |
inline LpBase::Constr operator<=(const LpBase::Value &e,
|
1602 |
1602 |
const LpBase::Expr &f) {
|
1603 |
1603 |
return LpBase::Constr(e, f, LpBase::NaN);
|
1604 |
1604 |
}
|
1605 |
1605 |
|
1606 |
1606 |
///Create constraint
|
1607 |
1607 |
|
1608 |
1608 |
///\relates LpBase::Constr
|
1609 |
1609 |
///
|
1610 |
1610 |
inline LpBase::Constr operator<=(const LpBase::Expr &e,
|
1611 |
1611 |
const LpBase::Value &f) {
|
1612 |
1612 |
return LpBase::Constr(- LpBase::INF, e, f);
|
1613 |
1613 |
}
|
1614 |
1614 |
|
1615 |
1615 |
///Create constraint
|
1616 |
1616 |
|
1617 |
1617 |
///\relates LpBase::Constr
|
1618 |
1618 |
///
|
1619 |
1619 |
inline LpBase::Constr operator>=(const LpBase::Expr &e,
|
1620 |
1620 |
const LpBase::Expr &f) {
|
1621 |
1621 |
return LpBase::Constr(0, e - f, LpBase::INF);
|
1622 |
1622 |
}
|
1623 |
1623 |
|
1624 |
1624 |
|
1625 |
1625 |
///Create constraint
|
1626 |
1626 |
|
1627 |
1627 |
///\relates LpBase::Constr
|
1628 |
1628 |
///
|
1629 |
1629 |
inline LpBase::Constr operator>=(const LpBase::Value &e,
|
1630 |
1630 |
const LpBase::Expr &f) {
|
1631 |
1631 |
return LpBase::Constr(LpBase::NaN, f, e);
|
1632 |
1632 |
}
|
1633 |
1633 |
|
1634 |
1634 |
|
1635 |
1635 |
///Create constraint
|
1636 |
1636 |
|
1637 |
1637 |
///\relates LpBase::Constr
|
1638 |
1638 |
///
|
1639 |
1639 |
inline LpBase::Constr operator>=(const LpBase::Expr &e,
|
1640 |
1640 |
const LpBase::Value &f) {
|
1641 |
1641 |
return LpBase::Constr(f, e, LpBase::INF);
|
1642 |
1642 |
}
|
1643 |
1643 |
|
1644 |
1644 |
///Create constraint
|
1645 |
1645 |
|
1646 |
1646 |
///\relates LpBase::Constr
|
1647 |
1647 |
///
|
1648 |
1648 |
inline LpBase::Constr operator==(const LpBase::Expr &e,
|
1649 |
1649 |
const LpBase::Value &f) {
|
1650 |
1650 |
return LpBase::Constr(f, e, f);
|
1651 |
1651 |
}
|
1652 |
1652 |
|
1653 |
1653 |
///Create constraint
|
1654 |
1654 |
|
1655 |
1655 |
///\relates LpBase::Constr
|
1656 |
1656 |
///
|
1657 |
1657 |
inline LpBase::Constr operator==(const LpBase::Expr &e,
|
1658 |
1658 |
const LpBase::Expr &f) {
|
1659 |
1659 |
return LpBase::Constr(0, f - e, 0);
|
1660 |
1660 |
}
|
1661 |
1661 |
|
1662 |
1662 |
///Create constraint
|
1663 |
1663 |
|
1664 |
1664 |
///\relates LpBase::Constr
|
1665 |
1665 |
///
|
1666 |
1666 |
inline LpBase::Constr operator<=(const LpBase::Value &n,
|
1667 |
1667 |
const LpBase::Constr &c) {
|
1668 |
1668 |
LpBase::Constr tmp(c);
|
1669 |
|
LEMON_ASSERT(isnan(tmp.lowerBound()), "Wrong LP constraint");
|
|
1669 |
LEMON_ASSERT(isNaN(tmp.lowerBound()), "Wrong LP constraint");
|
1670 |
1670 |
tmp.lowerBound()=n;
|
1671 |
1671 |
return tmp;
|
1672 |
1672 |
}
|
1673 |
1673 |
///Create constraint
|
1674 |
1674 |
|
1675 |
1675 |
///\relates LpBase::Constr
|
1676 |
1676 |
///
|
1677 |
1677 |
inline LpBase::Constr operator<=(const LpBase::Constr &c,
|
1678 |
1678 |
const LpBase::Value &n)
|
1679 |
1679 |
{
|
1680 |
1680 |
LpBase::Constr tmp(c);
|
1681 |
|
LEMON_ASSERT(isnan(tmp.upperBound()), "Wrong LP constraint");
|
|
1681 |
LEMON_ASSERT(isNaN(tmp.upperBound()), "Wrong LP constraint");
|
1682 |
1682 |
tmp.upperBound()=n;
|
1683 |
1683 |
return tmp;
|
1684 |
1684 |
}
|
1685 |
1685 |
|
1686 |
1686 |
///Create constraint
|
1687 |
1687 |
|
1688 |
1688 |
///\relates LpBase::Constr
|
1689 |
1689 |
///
|
1690 |
1690 |
inline LpBase::Constr operator>=(const LpBase::Value &n,
|
1691 |
1691 |
const LpBase::Constr &c) {
|
1692 |
1692 |
LpBase::Constr tmp(c);
|
1693 |
|
LEMON_ASSERT(isnan(tmp.upperBound()), "Wrong LP constraint");
|
|
1693 |
LEMON_ASSERT(isNaN(tmp.upperBound()), "Wrong LP constraint");
|
1694 |
1694 |
tmp.upperBound()=n;
|
1695 |
1695 |
return tmp;
|
1696 |
1696 |
}
|
1697 |
1697 |
///Create constraint
|
1698 |
1698 |
|
1699 |
1699 |
///\relates LpBase::Constr
|
1700 |
1700 |
///
|
1701 |
1701 |
inline LpBase::Constr operator>=(const LpBase::Constr &c,
|
1702 |
1702 |
const LpBase::Value &n)
|
1703 |
1703 |
{
|
1704 |
1704 |
LpBase::Constr tmp(c);
|
1705 |
|
LEMON_ASSERT(isnan(tmp.lowerBound()), "Wrong LP constraint");
|
|
1705 |
LEMON_ASSERT(isNaN(tmp.lowerBound()), "Wrong LP constraint");
|
1706 |
1706 |
tmp.lowerBound()=n;
|
1707 |
1707 |
return tmp;
|
1708 |
1708 |
}
|
1709 |
1709 |
|
1710 |
1710 |
///Addition
|
1711 |
1711 |
|
1712 |
1712 |
///\relates LpBase::DualExpr
|
1713 |
1713 |
///
|
1714 |
1714 |
inline LpBase::DualExpr operator+(const LpBase::DualExpr &a,
|
1715 |
1715 |
const LpBase::DualExpr &b) {
|
1716 |
1716 |
LpBase::DualExpr tmp(a);
|
1717 |
1717 |
tmp+=b;
|
1718 |
1718 |
return tmp;
|
1719 |
1719 |
}
|
1720 |
1720 |
///Substraction
|
1721 |
1721 |
|
1722 |
1722 |
///\relates LpBase::DualExpr
|
1723 |
1723 |
///
|
1724 |
1724 |
inline LpBase::DualExpr operator-(const LpBase::DualExpr &a,
|
1725 |
1725 |
const LpBase::DualExpr &b) {
|
1726 |
1726 |
LpBase::DualExpr tmp(a);
|
1727 |
1727 |
tmp-=b;
|
1728 |
1728 |
return tmp;
|
1729 |
1729 |
}
|
1730 |
1730 |
///Multiply with constant
|
1731 |
1731 |
|
1732 |
1732 |
///\relates LpBase::DualExpr
|
1733 |
1733 |
///
|
1734 |
1734 |
inline LpBase::DualExpr operator*(const LpBase::DualExpr &a,
|
1735 |
1735 |
const LpBase::Value &b) {
|
1736 |
1736 |
LpBase::DualExpr tmp(a);
|
1737 |
1737 |
tmp*=b;
|
1738 |
1738 |
return tmp;
|
1739 |
1739 |
}
|
1740 |
1740 |
|
1741 |
1741 |
///Multiply with constant
|
1742 |
1742 |
|
1743 |
1743 |
///\relates LpBase::DualExpr
|
1744 |
1744 |
///
|
1745 |
1745 |
inline LpBase::DualExpr operator*(const LpBase::Value &a,
|
1746 |
1746 |
const LpBase::DualExpr &b) {
|
1747 |
1747 |
LpBase::DualExpr tmp(b);
|
1748 |
1748 |
tmp*=a;
|
1749 |
1749 |
return tmp;
|
1750 |
1750 |
}
|
1751 |
1751 |
///Divide with constant
|
1752 |
1752 |
|
1753 |
1753 |
///\relates LpBase::DualExpr
|
1754 |
1754 |
///
|
1755 |
1755 |
inline LpBase::DualExpr operator/(const LpBase::DualExpr &a,
|
1756 |
1756 |
const LpBase::Value &b) {
|
1757 |
1757 |
LpBase::DualExpr tmp(a);
|
1758 |
1758 |
tmp/=b;
|
1759 |
1759 |
return tmp;
|
1760 |
1760 |
}
|
1761 |
1761 |
|
1762 |
1762 |
/// \ingroup lp_group
|
1763 |
1763 |
///
|
1764 |
1764 |
/// \brief Common base class for LP solvers
|
1765 |
1765 |
///
|
1766 |
1766 |
/// This class is an abstract base class for LP solvers. This class
|
1767 |
1767 |
/// provides a full interface for set and modify an LP problem,
|
1768 |
1768 |
/// solve it and retrieve the solution. You can use one of the
|
1769 |
1769 |
/// descendants as a concrete implementation, or the \c Lp
|
1770 |
1770 |
/// default LP solver. However, if you would like to handle LP
|
1771 |
1771 |
/// solvers as reference or pointer in a generic way, you can use
|
1772 |
1772 |
/// this class directly.
|
1773 |
1773 |
class LpSolver : virtual public LpBase {
|
1774 |
1774 |
public:
|
1775 |
1775 |
|
1776 |
1776 |
/// The problem types for primal and dual problems
|
1777 |
1777 |
enum ProblemType {
|
1778 |
1778 |
///Feasible solution hasn't been found (but may exist).
|
1779 |
1779 |
UNDEFINED = 0,
|
1780 |
1780 |
///The problem has no feasible solution
|
1781 |
1781 |
INFEASIBLE = 1,
|
1782 |
1782 |
///Feasible solution found
|
1783 |
1783 |
FEASIBLE = 2,
|
1784 |
1784 |
///Optimal solution exists and found
|
1785 |
1785 |
OPTIMAL = 3,
|
1786 |
1786 |
///The cost function is unbounded
|
1787 |
1787 |
UNBOUNDED = 4
|
1788 |
1788 |
};
|
1789 |
1789 |
|
1790 |
1790 |
///The basis status of variables
|
1791 |
1791 |
enum VarStatus {
|
1792 |
1792 |
/// The variable is in the basis
|
1793 |
1793 |
BASIC,
|
1794 |
1794 |
/// The variable is free, but not basic
|
1795 |
1795 |
FREE,
|
1796 |
1796 |
/// The variable has active lower bound
|
1797 |
1797 |
LOWER,
|
1798 |
1798 |
/// The variable has active upper bound
|
1799 |
1799 |
UPPER,
|
1800 |
1800 |
/// The variable is non-basic and fixed
|
1801 |
1801 |
FIXED
|