diff -r bafe730864da -r 61bf01404ede lemon/network_simplex.h --- a/lemon/network_simplex.h Mon Jun 01 16:53:59 2009 +0000 +++ b/lemon/network_simplex.h Thu Jun 04 01:19:06 2009 +0000 @@ -241,10 +241,10 @@ BlockSearchPivotRule(NetworkSimplex &ns) : _edge(ns._edge), _source(ns._source), _target(ns._target), _cost(ns._cost), _state(ns._state), _pi(ns._pi), - _in_edge(ns._in_edge), _edge_num(ns._edge_num + ns._node_num), _next_edge(0) + _in_edge(ns._in_edge), _edge_num(ns._edge_num), _next_edge(0) { // The main parameters of the pivot rule - const double BLOCK_SIZE_FACTOR = 2.0; + const double BLOCK_SIZE_FACTOR = 0.5; const int MIN_BLOCK_SIZE = 10; _block_size = std::max( int(BLOCK_SIZE_FACTOR * sqrt(_edge_num)), @@ -255,33 +255,32 @@ bool findEnteringEdge() { Cost c, min = 0; int cnt = _block_size; - int e, min_edge = _next_edge; + int e; for (e = _next_edge; e < _edge_num; ++e) { c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); if (c < min) { min = c; - min_edge = e; + _in_edge = e; } if (--cnt == 0) { - if (min < 0) break; + if (min < 0) goto search_end; cnt = _block_size; } } - if (min == 0 || cnt > 0) { - for (e = 0; e < _next_edge; ++e) { - c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); - if (c < min) { - min = c; - min_edge = e; - } - if (--cnt == 0) { - if (min < 0) break; - cnt = _block_size; - } + for (e = 0; e < _next_edge; ++e) { + c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); + if (c < min) { + min = c; + _in_edge = e; + } + if (--cnt == 0) { + if (min < 0) goto search_end; + cnt = _block_size; } } if (min >= 0) return false; - _in_edge = min_edge; + + search_end: _next_edge = e; return true; } @@ -325,7 +324,7 @@ _in_edge(ns._in_edge), _edge_num(ns._edge_num), _next_edge(0) { // The main parameters of the pivot rule - const double LIST_LENGTH_FACTOR = 1.0; + const double LIST_LENGTH_FACTOR = 0.25; const int MIN_LIST_LENGTH = 10; const double MINOR_LIMIT_FACTOR = 0.1; const int MIN_MINOR_LIMIT = 3; @@ -341,7 +340,7 @@ /// Find next entering edge bool findEnteringEdge() { Cost min, c; - int e, min_edge = _next_edge; + int e; if (_curr_length > 0 && _minor_count < _minor_limit) { // Minor iteration: select the best eligible edge from the // current candidate list @@ -352,16 +351,13 @@ c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); if (c < min) { min = c; - min_edge = e; + _in_edge = e; } - if (c >= 0) { + else if (c >= 0) { _candidates[i--] = _candidates[--_curr_length]; } } - if (min < 0) { - _in_edge = min_edge; - return true; - } + if (min < 0) return true; } // Major iteration: build a new candidate list @@ -373,27 +369,26 @@ _candidates[_curr_length++] = e; if (c < min) { min = c; - min_edge = e; + _in_edge = e; } - if (_curr_length == _list_length) break; + if (_curr_length == _list_length) goto search_end; } } - if (_curr_length < _list_length) { - for (e = 0; e < _next_edge; ++e) { - c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); - if (c < 0) { - _candidates[_curr_length++] = e; - if (c < min) { - min = c; - min_edge = e; - } - if (_curr_length == _list_length) break; + for (e = 0; e < _next_edge; ++e) { + c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); + if (c < 0) { + _candidates[_curr_length++] = e; + if (c < min) { + min = c; + _in_edge = e; } + if (_curr_length == _list_length) goto search_end; } } if (_curr_length == 0) return false; + + search_end: _minor_count = 1; - _in_edge = min_edge; _next_edge = e; return true; } @@ -451,7 +446,7 @@ _next_edge(0), _cand_cost(ns._edge_num),_sort_func(_cand_cost) { // The main parameters of the pivot rule - const double BLOCK_SIZE_FACTOR = 1.5; + const double BLOCK_SIZE_FACTOR = 1.0; const int MIN_BLOCK_SIZE = 10; const double HEAD_LENGTH_FACTOR = 0.1; const int MIN_HEAD_LENGTH = 3; @@ -479,39 +474,35 @@ // Extend the list int cnt = _block_size; - int last_edge = 0; int limit = _head_length; - for (int e = _next_edge; e < _edge_num; ++e) { + for (e = _next_edge; e < _edge_num; ++e) { _cand_cost[e] = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); if (_cand_cost[e] < 0) { _candidates[_curr_length++] = e; - last_edge = e; } if (--cnt == 0) { - if (_curr_length > limit) break; + if (_curr_length > limit) goto search_end; limit = 0; cnt = _block_size; } } - if (_curr_length <= limit) { - for (int e = 0; e < _next_edge; ++e) { - _cand_cost[e] = _state[e] * - (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); - if (_cand_cost[e] < 0) { - _candidates[_curr_length++] = e; - last_edge = e; - } - if (--cnt == 0) { - if (_curr_length > limit) break; - limit = 0; - cnt = _block_size; - } + for (e = 0; e < _next_edge; ++e) { + _cand_cost[e] = _state[e] * + (_cost[e] + _pi[_source[e]] - _pi[_target[e]]); + if (_cand_cost[e] < 0) { + _candidates[_curr_length++] = e; + } + if (--cnt == 0) { + if (_curr_length > limit) goto search_end; + limit = 0; + cnt = _block_size; } } if (_curr_length == 0) return false; - _next_edge = last_edge + 1; + + search_end: // Make heap of the candidate list (approximating a partial sort) make_heap( _candidates.begin(), _candidates.begin() + _curr_length, @@ -519,6 +510,7 @@ // Pop the first element of the heap _in_edge = _candidates[0]; + _next_edge = e; pop_heap( _candidates.begin(), _candidates.begin() + _curr_length, _sort_func ); _curr_length = std::min(_head_length, _curr_length - 1); @@ -931,12 +923,10 @@ _supply[_root] = 0; _pi[_root] = 0; - // Store the edges in a mixed order - int k = std::max(int(sqrt(_edge_num)), 10); + // Store the edges int i = 0; for (EdgeIt e(_orig_graph); e != INVALID; ++e) { - _edge[i] = e; - if ((i += k) >= _edge_num) i = (i % k) + 1; + _edge[i++] = e; } // Initialize edge maps @@ -961,7 +951,7 @@ } // Add artificial edges and initialize the spanning tree data structure - Cost max_cost = std::numeric_limits::max() / 4; + Cost max_cost = std::numeric_limits::max() / 2 + 1; Capacity max_cap = std::numeric_limits::max(); for (int u = 0, e = _edge_num; u != _node_num; ++u, ++e) { _parent[u] = _root;