No children
0
2
0
... | ... |
@@ -58,66 +58,66 @@ |
58 | 58 |
typedef BinHeap<Cost, RangeMap<int> > Heap; |
59 | 59 |
}; |
60 | 60 |
|
61 | 61 |
/// \addtogroup min_cost_flow_algs |
62 | 62 |
/// @{ |
63 | 63 |
|
64 | 64 |
/// \brief Implementation of the Capacity Scaling algorithm for |
65 | 65 |
/// finding a \ref min_cost_flow "minimum cost flow". |
66 | 66 |
/// |
67 | 67 |
/// \ref CapacityScaling implements the capacity scaling version |
68 | 68 |
/// of the successive shortest path algorithm for finding a |
69 | 69 |
/// \ref min_cost_flow "minimum cost flow" \ref amo93networkflows, |
70 | 70 |
/// \ref edmondskarp72theoretical. It is an efficient dual |
71 | 71 |
/// solution method. |
72 | 72 |
/// |
73 | 73 |
/// Most of the parameters of the problem (except for the digraph) |
74 | 74 |
/// can be given using separate functions, and the algorithm can be |
75 | 75 |
/// executed using the \ref run() function. If some parameters are not |
76 | 76 |
/// specified, then default values will be used. |
77 | 77 |
/// |
78 | 78 |
/// \tparam GR The digraph type the algorithm runs on. |
79 | 79 |
/// \tparam V The number type used for flow amounts, capacity bounds |
80 | 80 |
/// and supply values in the algorithm. By default, it is \c int. |
81 | 81 |
/// \tparam C The number type used for costs and potentials in the |
82 | 82 |
/// algorithm. By default, it is the same as \c V. |
83 | 83 |
/// \tparam TR The traits class that defines various types used by the |
84 | 84 |
/// algorithm. By default, it is \ref CapacityScalingDefaultTraits |
85 | 85 |
/// "CapacityScalingDefaultTraits<GR, V, C>". |
86 | 86 |
/// In most cases, this parameter should not be set directly, |
87 | 87 |
/// consider to use the named template parameters instead. |
88 | 88 |
/// |
89 | 89 |
/// \warning Both \c V and \c C must be signed number types. |
90 |
/// \warning All input data (capacities, supply values, and costs) must |
|
91 |
/// be integer. |
|
90 |
/// \warning Capacity bounds and supply values must be integer, but |
|
91 |
/// arc costs can be arbitrary real numbers. |
|
92 | 92 |
/// \warning This algorithm does not support negative costs for |
93 | 93 |
/// arcs having infinite upper bound. |
94 | 94 |
#ifdef DOXYGEN |
95 | 95 |
template <typename GR, typename V, typename C, typename TR> |
96 | 96 |
#else |
97 | 97 |
template < typename GR, typename V = int, typename C = V, |
98 | 98 |
typename TR = CapacityScalingDefaultTraits<GR, V, C> > |
99 | 99 |
#endif |
100 | 100 |
class CapacityScaling |
101 | 101 |
{ |
102 | 102 |
public: |
103 | 103 |
|
104 | 104 |
/// The type of the digraph |
105 | 105 |
typedef typename TR::Digraph Digraph; |
106 | 106 |
/// The type of the flow amounts, capacity bounds and supply values |
107 | 107 |
typedef typename TR::Value Value; |
108 | 108 |
/// The type of the arc costs |
109 | 109 |
typedef typename TR::Cost Cost; |
110 | 110 |
|
111 | 111 |
/// The type of the heap used for internal Dijkstra computations |
112 | 112 |
typedef typename TR::Heap Heap; |
113 | 113 |
|
114 | 114 |
/// The \ref CapacityScalingDefaultTraits "traits class" of the algorithm |
115 | 115 |
typedef TR Traits; |
116 | 116 |
|
117 | 117 |
public: |
118 | 118 |
|
119 | 119 |
/// \brief Problem type constants for the \c run() function. |
120 | 120 |
/// |
121 | 121 |
/// Enum type containing the problem type constants that can be |
122 | 122 |
/// returned by the \ref run() function of the algorithm. |
123 | 123 |
enum ProblemType { |
... | ... |
@@ -93,98 +93,101 @@ |
93 | 93 |
/// bounded), and the algorithm has found optimal flow and node |
94 | 94 |
/// potentials (primal and dual solutions). |
95 | 95 |
OPTIMAL, |
96 | 96 |
/// The objective function of the problem is unbounded, i.e. |
97 | 97 |
/// there is a directed cycle having negative total cost and |
98 | 98 |
/// infinite upper bound. |
99 | 99 |
UNBOUNDED |
100 | 100 |
}; |
101 | 101 |
|
102 | 102 |
/// \brief Constants for selecting the type of the supply constraints. |
103 | 103 |
/// |
104 | 104 |
/// Enum type containing constants for selecting the supply type, |
105 | 105 |
/// i.e. the direction of the inequalities in the supply/demand |
106 | 106 |
/// constraints of the \ref min_cost_flow "minimum cost flow problem". |
107 | 107 |
/// |
108 | 108 |
/// The default supply type is \c GEQ, the \c LEQ type can be |
109 | 109 |
/// selected using \ref supplyType(). |
110 | 110 |
/// The equality form is a special case of both supply types. |
111 | 111 |
enum SupplyType { |
112 | 112 |
/// This option means that there are <em>"greater or equal"</em> |
113 | 113 |
/// supply/demand constraints in the definition of the problem. |
114 | 114 |
GEQ, |
115 | 115 |
/// This option means that there are <em>"less or equal"</em> |
116 | 116 |
/// supply/demand constraints in the definition of the problem. |
117 | 117 |
LEQ |
118 | 118 |
}; |
119 | 119 |
|
120 | 120 |
/// \brief Constants for selecting the pivot rule. |
121 | 121 |
/// |
122 | 122 |
/// Enum type containing constants for selecting the pivot rule for |
123 | 123 |
/// the \ref run() function. |
124 | 124 |
/// |
125 |
/// \ref NetworkSimplex provides five different pivot rule |
|
126 |
/// implementations that significantly affect the running time |
|
125 |
/// \ref NetworkSimplex provides five different implementations for |
|
126 |
/// the pivot strategy that significantly affects the running time |
|
127 | 127 |
/// of the algorithm. |
128 |
/// By default, \ref BLOCK_SEARCH "Block Search" is used, which |
|
129 |
/// turend out to be the most efficient and the most robust on various |
|
130 |
/// test inputs. |
|
131 |
/// However, another pivot rule can be selected using the \ref run() |
|
132 |
/// |
|
128 |
/// According to experimental tests conducted on various problem |
|
129 |
/// instances, \ref BLOCK_SEARCH "Block Search" and |
|
130 |
/// \ref ALTERING_LIST "Altering Candidate List" rules turned out |
|
131 |
/// to be the most efficient. |
|
132 |
/// Since \ref BLOCK_SEARCH "Block Search" is a simpler strategy that |
|
133 |
/// seemed to be slightly more robust, it is used by default. |
|
134 |
/// However, another pivot rule can easily be selected using the |
|
135 |
/// \ref run() function with the proper parameter. |
|
133 | 136 |
enum PivotRule { |
134 | 137 |
|
135 | 138 |
/// The \e First \e Eligible pivot rule. |
136 | 139 |
/// The next eligible arc is selected in a wraparound fashion |
137 | 140 |
/// in every iteration. |
138 | 141 |
FIRST_ELIGIBLE, |
139 | 142 |
|
140 | 143 |
/// The \e Best \e Eligible pivot rule. |
141 | 144 |
/// The best eligible arc is selected in every iteration. |
142 | 145 |
BEST_ELIGIBLE, |
143 | 146 |
|
144 | 147 |
/// The \e Block \e Search pivot rule. |
145 | 148 |
/// A specified number of arcs are examined in every iteration |
146 | 149 |
/// in a wraparound fashion and the best eligible arc is selected |
147 | 150 |
/// from this block. |
148 | 151 |
BLOCK_SEARCH, |
149 | 152 |
|
150 | 153 |
/// The \e Candidate \e List pivot rule. |
151 | 154 |
/// In a major iteration a candidate list is built from eligible arcs |
152 | 155 |
/// in a wraparound fashion and in the following minor iterations |
153 | 156 |
/// the best eligible arc is selected from this list. |
154 | 157 |
CANDIDATE_LIST, |
155 | 158 |
|
156 | 159 |
/// The \e Altering \e Candidate \e List pivot rule. |
157 | 160 |
/// It is a modified version of the Candidate List method. |
158 |
/// It keeps only |
|
161 |
/// It keeps only a few of the best eligible arcs from the former |
|
159 | 162 |
/// candidate list and extends this list in every iteration. |
160 | 163 |
ALTERING_LIST |
161 | 164 |
}; |
162 | 165 |
|
163 | 166 |
private: |
164 | 167 |
|
165 | 168 |
TEMPLATE_DIGRAPH_TYPEDEFS(GR); |
166 | 169 |
|
167 | 170 |
typedef std::vector<int> IntVector; |
168 | 171 |
typedef std::vector<Value> ValueVector; |
169 | 172 |
typedef std::vector<Cost> CostVector; |
170 | 173 |
typedef std::vector<signed char> CharVector; |
171 | 174 |
// Note: vector<signed char> is used instead of vector<ArcState> and |
172 | 175 |
// vector<ArcDirection> for efficiency reasons |
173 | 176 |
|
174 | 177 |
// State constants for arcs |
175 | 178 |
enum ArcState { |
176 | 179 |
STATE_UPPER = -1, |
177 | 180 |
STATE_TREE = 0, |
178 | 181 |
STATE_LOWER = 1 |
179 | 182 |
}; |
180 | 183 |
|
181 | 184 |
// Direction constants for tree arcs |
182 | 185 |
enum ArcDirection { |
183 | 186 |
DIR_DOWN = -1, |
184 | 187 |
DIR_UP = 1 |
185 | 188 |
}; |
186 | 189 |
|
187 | 190 |
private: |
188 | 191 |
|
189 | 192 |
// Data related to the underlying digraph |
190 | 193 |
const GR &_graph; |
... | ... |
@@ -509,151 +512,151 @@ |
509 | 512 |
|
510 | 513 |
}; //class CandidateListPivotRule |
511 | 514 |
|
512 | 515 |
|
513 | 516 |
// Implementation of the Altering Candidate List pivot rule |
514 | 517 |
class AlteringListPivotRule |
515 | 518 |
{ |
516 | 519 |
private: |
517 | 520 |
|
518 | 521 |
// References to the NetworkSimplex class |
519 | 522 |
const IntVector &_source; |
520 | 523 |
const IntVector &_target; |
521 | 524 |
const CostVector &_cost; |
522 | 525 |
const CharVector &_state; |
523 | 526 |
const CostVector &_pi; |
524 | 527 |
int &_in_arc; |
525 | 528 |
int _search_arc_num; |
526 | 529 |
|
527 | 530 |
// Pivot rule data |
528 | 531 |
int _block_size, _head_length, _curr_length; |
529 | 532 |
int _next_arc; |
530 | 533 |
IntVector _candidates; |
531 | 534 |
CostVector _cand_cost; |
532 | 535 |
|
533 | 536 |
// Functor class to compare arcs during sort of the candidate list |
534 | 537 |
class SortFunc |
535 | 538 |
{ |
536 | 539 |
private: |
537 | 540 |
const CostVector &_map; |
538 | 541 |
public: |
539 | 542 |
SortFunc(const CostVector &map) : _map(map) {} |
540 | 543 |
bool operator()(int left, int right) { |
541 |
return _map[left] |
|
544 |
return _map[left] < _map[right]; |
|
542 | 545 |
} |
543 | 546 |
}; |
544 | 547 |
|
545 | 548 |
SortFunc _sort_func; |
546 | 549 |
|
547 | 550 |
public: |
548 | 551 |
|
549 | 552 |
// Constructor |
550 | 553 |
AlteringListPivotRule(NetworkSimplex &ns) : |
551 | 554 |
_source(ns._source), _target(ns._target), |
552 | 555 |
_cost(ns._cost), _state(ns._state), _pi(ns._pi), |
553 | 556 |
_in_arc(ns.in_arc), _search_arc_num(ns._search_arc_num), |
554 | 557 |
_next_arc(0), _cand_cost(ns._search_arc_num), _sort_func(_cand_cost) |
555 | 558 |
{ |
556 | 559 |
// The main parameters of the pivot rule |
557 | 560 |
const double BLOCK_SIZE_FACTOR = 1.0; |
558 | 561 |
const int MIN_BLOCK_SIZE = 10; |
559 |
const double HEAD_LENGTH_FACTOR = 0. |
|
562 |
const double HEAD_LENGTH_FACTOR = 0.01; |
|
560 | 563 |
const int MIN_HEAD_LENGTH = 3; |
561 | 564 |
|
562 | 565 |
_block_size = std::max( int(BLOCK_SIZE_FACTOR * |
563 | 566 |
std::sqrt(double(_search_arc_num))), |
564 | 567 |
MIN_BLOCK_SIZE ); |
565 | 568 |
_head_length = std::max( int(HEAD_LENGTH_FACTOR * _block_size), |
566 | 569 |
MIN_HEAD_LENGTH ); |
567 | 570 |
_candidates.resize(_head_length + _block_size); |
568 | 571 |
_curr_length = 0; |
569 | 572 |
} |
570 | 573 |
|
571 | 574 |
// Find next entering arc |
572 | 575 |
bool findEnteringArc() { |
573 | 576 |
// Check the current candidate list |
574 | 577 |
int e; |
575 | 578 |
Cost c; |
576 | 579 |
for (int i = 0; i != _curr_length; ++i) { |
577 | 580 |
e = _candidates[i]; |
578 | 581 |
c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); |
579 | 582 |
if (c < 0) { |
580 | 583 |
_cand_cost[e] = c; |
581 | 584 |
} else { |
582 | 585 |
_candidates[i--] = _candidates[--_curr_length]; |
583 | 586 |
} |
584 | 587 |
} |
585 | 588 |
|
586 | 589 |
// Extend the list |
587 | 590 |
int cnt = _block_size; |
588 | 591 |
int limit = _head_length; |
589 | 592 |
|
590 | 593 |
for (e = _next_arc; e != _search_arc_num; ++e) { |
591 | 594 |
c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); |
592 | 595 |
if (c < 0) { |
593 | 596 |
_cand_cost[e] = c; |
594 | 597 |
_candidates[_curr_length++] = e; |
595 | 598 |
} |
596 | 599 |
if (--cnt == 0) { |
597 | 600 |
if (_curr_length > limit) goto search_end; |
598 | 601 |
limit = 0; |
599 | 602 |
cnt = _block_size; |
600 | 603 |
} |
601 | 604 |
} |
602 | 605 |
for (e = 0; e != _next_arc; ++e) { |
603 |
_cand_cost[e] = _state[e] * |
|
604 |
(_cost[e] + _pi[_source[e]] - _pi[_target[e]]); |
|
605 |
|
|
606 |
c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); |
|
607 |
if (c < 0) { |
|
608 |
_cand_cost[e] = c; |
|
606 | 609 |
_candidates[_curr_length++] = e; |
607 | 610 |
} |
608 | 611 |
if (--cnt == 0) { |
609 | 612 |
if (_curr_length > limit) goto search_end; |
610 | 613 |
limit = 0; |
611 | 614 |
cnt = _block_size; |
612 | 615 |
} |
613 | 616 |
} |
614 | 617 |
if (_curr_length == 0) return false; |
615 | 618 |
|
616 | 619 |
search_end: |
617 | 620 |
|
618 |
// Make heap of the candidate list (approximating a partial sort) |
|
619 |
make_heap( _candidates.begin(), _candidates.begin() + _curr_length, |
|
620 |
|
|
621 |
// Perform partial sort operation on the candidate list |
|
622 |
int new_length = std::min(_head_length + 1, _curr_length); |
|
623 |
std::partial_sort(_candidates.begin(), _candidates.begin() + new_length, |
|
624 |
_candidates.begin() + _curr_length, _sort_func); |
|
621 | 625 |
|
622 |
// |
|
626 |
// Select the entering arc and remove it from the list |
|
623 | 627 |
_in_arc = _candidates[0]; |
624 | 628 |
_next_arc = e; |
625 |
pop_heap( _candidates.begin(), _candidates.begin() + _curr_length, |
|
626 |
_sort_func ); |
|
627 |
|
|
629 |
_candidates[0] = _candidates[new_length - 1]; |
|
630 |
_curr_length = new_length - 1; |
|
628 | 631 |
return true; |
629 | 632 |
} |
630 | 633 |
|
631 | 634 |
}; //class AlteringListPivotRule |
632 | 635 |
|
633 | 636 |
public: |
634 | 637 |
|
635 | 638 |
/// \brief Constructor. |
636 | 639 |
/// |
637 | 640 |
/// The constructor of the class. |
638 | 641 |
/// |
639 | 642 |
/// \param graph The digraph the algorithm runs on. |
640 | 643 |
/// \param arc_mixing Indicate if the arcs will be stored in a |
641 | 644 |
/// mixed order in the internal data structure. |
642 | 645 |
/// In general, it leads to similar performance as using the original |
643 | 646 |
/// arc order, but it makes the algorithm more robust and in special |
644 | 647 |
/// cases, even significantly faster. Therefore, it is enabled by default. |
645 | 648 |
NetworkSimplex(const GR& graph, bool arc_mixing = true) : |
646 | 649 |
_graph(graph), _node_id(graph), _arc_id(graph), |
647 | 650 |
_arc_mixing(arc_mixing), |
648 | 651 |
MAX(std::numeric_limits<Value>::max()), |
649 | 652 |
INF(std::numeric_limits<Value>::has_infinity ? |
650 | 653 |
std::numeric_limits<Value>::infinity() : MAX) |
651 | 654 |
{ |
652 | 655 |
// Check the number types |
653 | 656 |
LEMON_ASSERT(std::numeric_limits<Value>::is_signed, |
654 | 657 |
"The flow type of NetworkSimplex must be signed"); |
655 | 658 |
LEMON_ASSERT(std::numeric_limits<Cost>::is_signed, |
656 | 659 |
"The cost type of NetworkSimplex must be signed"); |
657 | 660 |
|
658 | 661 |
// Reset data structures |
659 | 662 |
reset(); |
0 comments (0 inline)