... | ... |
@@ -51,13 +51,14 @@ |
51 | 51 |
/// \tparam GR The digraph type the algorithm runs on. |
52 | 52 |
/// \tparam F The value type used for flow amounts, capacity bounds |
53 | 53 |
/// and supply values in the algorithm. By default it is \c int. |
54 | 54 |
/// \tparam C The value type used for costs and potentials in the |
55 | 55 |
/// algorithm. By default it is the same as \c F. |
56 | 56 |
/// |
57 |
/// \warning Both value types must be signed |
|
57 |
/// \warning Both value types must be signed and all input data must |
|
58 |
/// be integer. |
|
58 | 59 |
/// |
59 | 60 |
/// \note %NetworkSimplex provides five different pivot rule |
60 | 61 |
/// implementations. For more information see \ref PivotRule. |
61 | 62 |
template <typename GR, typename F = int, typename C = F> |
62 | 63 |
class NetworkSimplex |
63 | 64 |
{ |
... | ... |
@@ -1041,14 +1042,16 @@ |
1041 | 1042 |
for (ArcIt e(_graph); e != INVALID; ++e) { |
1042 | 1043 |
_arc_ref[i] = e; |
1043 | 1044 |
if ((i += k) >= _arc_num) i = (i % k) + 1; |
1044 | 1045 |
} |
1045 | 1046 |
|
1046 | 1047 |
// Initialize arc maps |
1047 |
Flow max_cap = std::numeric_limits<Flow>::max(); |
|
1048 |
Cost max_cost = std::numeric_limits<Cost>::max() / 4; |
|
1048 |
Flow inf_cap = |
|
1049 |
std::numeric_limits<Flow>::has_infinity ? |
|
1050 |
std::numeric_limits<Flow>::infinity() : |
|
1051 |
std::numeric_limits<Flow>::max(); |
|
1049 | 1052 |
if (_pupper && _pcost) { |
1050 | 1053 |
for (int i = 0; i != _arc_num; ++i) { |
1051 | 1054 |
Arc e = _arc_ref[i]; |
1052 | 1055 |
_source[i] = _node_id[_graph.source(e)]; |
1053 | 1056 |
_target[i] = _node_id[_graph.target(e)]; |
1054 | 1057 |
_cap[i] = (*_pupper)[e]; |
... | ... |
@@ -1066,22 +1069,34 @@ |
1066 | 1069 |
} |
1067 | 1070 |
if (_pupper) { |
1068 | 1071 |
for (int i = 0; i != _arc_num; ++i) |
1069 | 1072 |
_cap[i] = (*_pupper)[_arc_ref[i]]; |
1070 | 1073 |
} else { |
1071 | 1074 |
for (int i = 0; i != _arc_num; ++i) |
1072 |
_cap[i] = |
|
1075 |
_cap[i] = inf_cap; |
|
1073 | 1076 |
} |
1074 | 1077 |
if (_pcost) { |
1075 | 1078 |
for (int i = 0; i != _arc_num; ++i) |
1076 | 1079 |
_cost[i] = (*_pcost)[_arc_ref[i]]; |
1077 | 1080 |
} else { |
1078 | 1081 |
for (int i = 0; i != _arc_num; ++i) |
1079 | 1082 |
_cost[i] = 1; |
1080 | 1083 |
} |
1081 | 1084 |
} |
1085 |
|
|
1086 |
// Initialize artifical cost |
|
1087 |
Cost art_cost; |
|
1088 |
if (std::numeric_limits<Cost>::is_exact) { |
|
1089 |
art_cost = std::numeric_limits<Cost>::max() / 4 + 1; |
|
1090 |
} else { |
|
1091 |
art_cost = std::numeric_limits<Cost>::min(); |
|
1092 |
for (int i = 0; i != _arc_num; ++i) { |
|
1093 |
if (_cost[i] > art_cost) art_cost = _cost[i]; |
|
1094 |
} |
|
1095 |
art_cost = (art_cost + 1) * _node_num; |
|
1096 |
} |
|
1082 | 1097 |
|
1083 | 1098 |
// Remove non-zero lower bounds |
1084 | 1099 |
if (_plower) { |
1085 | 1100 |
for (int i = 0; i != _arc_num; ++i) { |
1086 | 1101 |
Flow c = (*_plower)[_arc_ref[i]]; |
1087 | 1102 |
if (c != 0) { |
... | ... |
@@ -1097,23 +1112,23 @@ |
1097 | 1112 |
_thread[u] = u + 1; |
1098 | 1113 |
_rev_thread[u + 1] = u; |
1099 | 1114 |
_succ_num[u] = 1; |
1100 | 1115 |
_last_succ[u] = u; |
1101 | 1116 |
_parent[u] = _root; |
1102 | 1117 |
_pred[u] = e; |
1103 |
_cost[e] = max_cost; |
|
1104 |
_cap[e] = max_cap; |
|
1118 |
_cost[e] = art_cost; |
|
1119 |
_cap[e] = inf_cap; |
|
1105 | 1120 |
_state[e] = STATE_TREE; |
1106 | 1121 |
if (_supply[u] >= 0) { |
1107 | 1122 |
_flow[e] = _supply[u]; |
1108 | 1123 |
_forward[u] = true; |
1109 |
_pi[u] = - |
|
1124 |
_pi[u] = -art_cost; |
|
1110 | 1125 |
} else { |
1111 | 1126 |
_flow[e] = -_supply[u]; |
1112 | 1127 |
_forward[u] = false; |
1113 |
_pi[u] = |
|
1128 |
_pi[u] = art_cost; |
|
1114 | 1129 |
} |
1115 | 1130 |
} |
1116 | 1131 |
|
1117 | 1132 |
return true; |
1118 | 1133 |
} |
1119 | 1134 |
|
... | ... |
@@ -1324,28 +1339,16 @@ |
1324 | 1339 |
|
1325 | 1340 |
// Update potentials |
1326 | 1341 |
void updatePotential() { |
1327 | 1342 |
Cost sigma = _forward[u_in] ? |
1328 | 1343 |
_pi[v_in] - _pi[u_in] - _cost[_pred[u_in]] : |
1329 | 1344 |
_pi[v_in] - _pi[u_in] + _cost[_pred[u_in]]; |
1330 |
if (_succ_num[u_in] > _node_num / 2) { |
|
1331 |
// Update in the upper subtree (which contains the root) |
|
1332 |
int before = _rev_thread[u_in]; |
|
1333 |
int after = _thread[_last_succ[u_in]]; |
|
1334 |
_thread[before] = after; |
|
1335 |
_pi[_root] -= sigma; |
|
1336 |
for (int u = _thread[_root]; u != _root; u = _thread[u]) { |
|
1337 |
_pi[u] -= sigma; |
|
1338 |
} |
|
1339 |
_thread[before] = u_in; |
|
1340 |
} else { |
|
1341 |
// Update in the lower subtree (which has been moved) |
|
1342 |
int end = _thread[_last_succ[u_in]]; |
|
1343 |
for (int u = u_in; u != end; u = _thread[u]) { |
|
1344 |
_pi[u] += sigma; |
|
1345 |
} |
|
1345 |
// Update potentials in the subtree, which has been moved |
|
1346 |
int end = _thread[_last_succ[u_in]]; |
|
1347 |
for (int u = u_in; u != end; u = _thread[u]) { |
|
1348 |
_pi[u] += sigma; |
|
1346 | 1349 |
} |
1347 | 1350 |
} |
1348 | 1351 |
|
1349 | 1352 |
// Execute the algorithm |
1350 | 1353 |
bool start(PivotRule pivot_rule) { |
1351 | 1354 |
// Select the pivot rule implementation |
0 comments (0 inline)