0
2
0
| ... | ... |
@@ -132,9 +132,9 @@ |
| 132 | 132 |
|
| 133 | 133 |
TEMPLATE_DIGRAPH_TYPEDEFS(GR); |
| 134 | 134 |
|
| 135 | 135 |
typedef std::vector<int> IntVector; |
| 136 |
typedef std::vector< |
|
| 136 |
typedef std::vector<char> BoolVector; |
|
| 137 | 137 |
typedef std::vector<Value> ValueVector; |
| 138 | 138 |
typedef std::vector<Cost> CostVector; |
| 139 | 139 |
|
| 140 | 140 |
private: |
| ... | ... |
@@ -195,8 +195,9 @@ |
| 195 | 195 |
{
|
| 196 | 196 |
private: |
| 197 | 197 |
|
| 198 | 198 |
int _node_num; |
| 199 |
bool _geq; |
|
| 199 | 200 |
const IntVector &_first_out; |
| 200 | 201 |
const IntVector &_target; |
| 201 | 202 |
const CostVector &_cost; |
| 202 | 203 |
const ValueVector &_res_cap; |
| ... | ... |
@@ -209,12 +210,12 @@ |
| 209 | 210 |
|
| 210 | 211 |
public: |
| 211 | 212 |
|
| 212 | 213 |
ResidualDijkstra(CapacityScaling& cs) : |
| 213 |
_node_num(cs._node_num), _first_out(cs._first_out), |
|
| 214 |
_target(cs._target), _cost(cs._cost), _res_cap(cs._res_cap), |
|
| 215 |
_excess(cs._excess), _pi(cs._pi), _pred(cs._pred), |
|
| 216 |
_dist(cs._node_num) |
|
| 214 |
_node_num(cs._node_num), _geq(cs._sum_supply < 0), |
|
| 215 |
_first_out(cs._first_out), _target(cs._target), _cost(cs._cost), |
|
| 216 |
_res_cap(cs._res_cap), _excess(cs._excess), _pi(cs._pi), |
|
| 217 |
_pred(cs._pred), _dist(cs._node_num) |
|
| 217 | 218 |
{}
|
| 218 | 219 |
|
| 219 | 220 |
int run(int s, Value delta = 1) {
|
| 220 | 221 |
RangeMap<int> heap_cross_ref(_node_num, Heap::PRE_HEAP); |
| ... | ... |
@@ -231,9 +232,10 @@ |
| 231 | 232 |
_proc_nodes.push_back(u); |
| 232 | 233 |
heap.pop(); |
| 233 | 234 |
|
| 234 | 235 |
// Traverse outgoing residual arcs |
| 235 |
|
|
| 236 |
int last_out = _geq ? _first_out[u+1] : _first_out[u+1] - 1; |
|
| 237 |
for (int a = _first_out[u]; a != last_out; ++a) {
|
|
| 236 | 238 |
if (_res_cap[a] < delta) continue; |
| 237 | 239 |
v = _target[a]; |
| 238 | 240 |
switch (heap.state(v)) {
|
| 239 | 241 |
case Heap::PRE_HEAP: |
| ... | ... |
@@ -686,24 +688,27 @@ |
| 686 | 688 |
_sum_supply += _supply[i]; |
| 687 | 689 |
} |
| 688 | 690 |
if (_sum_supply > 0) return INFEASIBLE; |
| 689 | 691 |
|
| 690 |
// Initialize |
|
| 692 |
// Initialize vectors |
|
| 691 | 693 |
for (int i = 0; i != _root; ++i) {
|
| 692 | 694 |
_pi[i] = 0; |
| 693 | 695 |
_excess[i] = _supply[i]; |
| 694 | 696 |
} |
| 695 | 697 |
|
| 696 | 698 |
// Remove non-zero lower bounds |
| 699 |
const Value MAX = std::numeric_limits<Value>::max(); |
|
| 700 |
int last_out; |
|
| 697 | 701 |
if (_have_lower) {
|
| 698 | 702 |
for (int i = 0; i != _root; ++i) {
|
| 699 |
|
|
| 703 |
last_out = _first_out[i+1]; |
|
| 704 |
for (int j = _first_out[i]; j != last_out; ++j) {
|
|
| 700 | 705 |
if (_forward[j]) {
|
| 701 | 706 |
Value c = _lower[j]; |
| 702 | 707 |
if (c >= 0) {
|
| 703 |
_res_cap[j] = _upper[j] < |
|
| 708 |
_res_cap[j] = _upper[j] < MAX ? _upper[j] - c : INF; |
|
| 704 | 709 |
} else {
|
| 705 |
_res_cap[j] = _upper[j] < |
|
| 710 |
_res_cap[j] = _upper[j] < MAX + c ? _upper[j] - c : INF; |
|
| 706 | 711 |
} |
| 707 | 712 |
_excess[i] -= c; |
| 708 | 713 |
_excess[_target[j]] += c; |
| 709 | 714 |
} else {
|
| ... | ... |
@@ -717,17 +722,18 @@ |
| 717 | 722 |
} |
| 718 | 723 |
} |
| 719 | 724 |
|
| 720 | 725 |
// Handle negative costs |
| 721 |
for (int u = 0; u != _root; ++u) {
|
|
| 722 |
for (int a = _first_out[u]; a != _first_out[u+1]; ++a) {
|
|
| 723 |
Value rc = _res_cap[a]; |
|
| 724 |
if (_cost[a] < 0 && rc > 0) {
|
|
| 725 |
if (rc == INF) return UNBOUNDED; |
|
| 726 |
_excess[u] -= rc; |
|
| 727 |
_excess[_target[a]] += rc; |
|
| 728 |
_res_cap[a] = 0; |
|
| 729 |
|
|
| 726 |
for (int i = 0; i != _root; ++i) {
|
|
| 727 |
last_out = _first_out[i+1] - 1; |
|
| 728 |
for (int j = _first_out[i]; j != last_out; ++j) {
|
|
| 729 |
Value rc = _res_cap[j]; |
|
| 730 |
if (_cost[j] < 0 && rc > 0) {
|
|
| 731 |
if (rc >= MAX) return UNBOUNDED; |
|
| 732 |
_excess[i] -= rc; |
|
| 733 |
_excess[_target[j]] += rc; |
|
| 734 |
_res_cap[j] = 0; |
|
| 735 |
_res_cap[_reverse[j]] += rc; |
|
| 730 | 736 |
} |
| 731 | 737 |
} |
| 732 | 738 |
} |
| 733 | 739 |
|
| ... | ... |
@@ -735,36 +741,34 @@ |
| 735 | 741 |
if (_sum_supply < 0) {
|
| 736 | 742 |
_pi[_root] = 0; |
| 737 | 743 |
_excess[_root] = -_sum_supply; |
| 738 | 744 |
for (int a = _first_out[_root]; a != _res_arc_num; ++a) {
|
| 739 |
int u = _target[a]; |
|
| 740 |
if (_excess[u] < 0) {
|
|
| 741 |
_res_cap[a] = -_excess[u] + 1; |
|
| 742 |
} else {
|
|
| 743 |
_res_cap[a] = 1; |
|
| 744 |
} |
|
| 745 |
|
|
| 745 |
int ra = _reverse[a]; |
|
| 746 |
_res_cap[a] = -_sum_supply + 1; |
|
| 747 |
_res_cap[ra] = 0; |
|
| 746 | 748 |
_cost[a] = 0; |
| 747 |
_cost[ |
|
| 749 |
_cost[ra] = 0; |
|
| 748 | 750 |
} |
| 749 | 751 |
} else {
|
| 750 | 752 |
_pi[_root] = 0; |
| 751 | 753 |
_excess[_root] = 0; |
| 752 | 754 |
for (int a = _first_out[_root]; a != _res_arc_num; ++a) {
|
| 755 |
int ra = _reverse[a]; |
|
| 753 | 756 |
_res_cap[a] = 1; |
| 754 |
_res_cap[ |
|
| 757 |
_res_cap[ra] = 0; |
|
| 755 | 758 |
_cost[a] = 0; |
| 756 |
_cost[ |
|
| 759 |
_cost[ra] = 0; |
|
| 757 | 760 |
} |
| 758 | 761 |
} |
| 759 | 762 |
|
| 760 | 763 |
// Initialize delta value |
| 761 | 764 |
if (_factor > 1) {
|
| 762 | 765 |
// With scaling |
| 763 | 766 |
Value max_sup = 0, max_dem = 0; |
| 764 | 767 |
for (int i = 0; i != _node_num; ++i) {
|
| 765 |
if ( _excess[i] > max_sup) max_sup = _excess[i]; |
|
| 766 |
if (-_excess[i] > max_dem) max_dem = -_excess[i]; |
|
| 768 |
Value ex = _excess[i]; |
|
| 769 |
if ( ex > max_sup) max_sup = ex; |
|
| 770 |
if (-ex > max_dem) max_dem = -ex; |
|
| 767 | 771 |
} |
| 768 | 772 |
Value max_cap = 0; |
| 769 | 773 |
for (int j = 0; j != _res_arc_num; ++j) {
|
| 770 | 774 |
if (_res_cap[j] > max_cap) max_cap = _res_cap[j]; |
| ... | ... |
@@ -788,9 +792,10 @@ |
| 788 | 792 |
pt = startWithoutScaling(); |
| 789 | 793 |
|
| 790 | 794 |
// Handle non-zero lower bounds |
| 791 | 795 |
if (_have_lower) {
|
| 792 |
|
|
| 796 |
int limit = _first_out[_root]; |
|
| 797 |
for (int j = 0; j != limit; ++j) {
|
|
| 793 | 798 |
if (!_forward[j]) _res_cap[j] += _lower[j]; |
| 794 | 799 |
} |
| 795 | 800 |
} |
| 796 | 801 |
|
| ... | ... |
@@ -811,10 +816,13 @@ |
| 811 | 816 |
int s, t; |
| 812 | 817 |
ResidualDijkstra _dijkstra(*this); |
| 813 | 818 |
while (true) {
|
| 814 | 819 |
// Saturate all arcs not satisfying the optimality condition |
| 820 |
int last_out; |
|
| 815 | 821 |
for (int u = 0; u != _node_num; ++u) {
|
| 816 |
|
|
| 822 |
last_out = _sum_supply < 0 ? |
|
| 823 |
_first_out[u+1] : _first_out[u+1] - 1; |
|
| 824 |
for (int a = _first_out[u]; a != last_out; ++a) {
|
|
| 817 | 825 |
int v = _target[a]; |
| 818 | 826 |
Cost c = _cost[a] + _pi[u] - _pi[v]; |
| 819 | 827 |
Value rc = _res_cap[a]; |
| 820 | 828 |
if (c < 0 && rc >= _delta) {
|
| ... | ... |
@@ -829,10 +837,11 @@ |
| 829 | 837 |
// Find excess nodes and deficit nodes |
| 830 | 838 |
_excess_nodes.clear(); |
| 831 | 839 |
_deficit_nodes.clear(); |
| 832 | 840 |
for (int u = 0; u != _node_num; ++u) {
|
| 833 |
if (_excess[u] >= _delta) _excess_nodes.push_back(u); |
|
| 834 |
if (_excess[u] <= -_delta) _deficit_nodes.push_back(u); |
|
| 841 |
Value ex = _excess[u]; |
|
| 842 |
if (ex >= _delta) _excess_nodes.push_back(u); |
|
| 843 |
if (ex <= -_delta) _deficit_nodes.push_back(u); |
|
| 835 | 844 |
} |
| 836 | 845 |
int next_node = 0, next_def_node = 0; |
| 837 | 846 |
|
| 838 | 847 |
// Find augmenting shortest paths |
| ... | ... |
@@ -163,9 +163,9 @@ |
| 163 | 163 |
|
| 164 | 164 |
TEMPLATE_DIGRAPH_TYPEDEFS(GR); |
| 165 | 165 |
|
| 166 | 166 |
typedef std::vector<int> IntVector; |
| 167 |
typedef std::vector< |
|
| 167 |
typedef std::vector<char> CharVector; |
|
| 168 | 168 |
typedef std::vector<Value> ValueVector; |
| 169 | 169 |
typedef std::vector<Cost> CostVector; |
| 170 | 170 |
|
| 171 | 171 |
// State constants for arcs |
| ... | ... |
@@ -211,17 +211,19 @@ |
| 211 | 211 |
IntVector _rev_thread; |
| 212 | 212 |
IntVector _succ_num; |
| 213 | 213 |
IntVector _last_succ; |
| 214 | 214 |
IntVector _dirty_revs; |
| 215 |
BoolVector _forward; |
|
| 216 |
IntVector _state; |
|
| 215 |
CharVector _forward; |
|
| 216 |
CharVector _state; |
|
| 217 | 217 |
int _root; |
| 218 | 218 |
|
| 219 | 219 |
// Temporary data used in the current pivot iteration |
| 220 | 220 |
int in_arc, join, u_in, v_in, u_out, v_out; |
| 221 | 221 |
int first, second, right, last; |
| 222 | 222 |
int stem, par_stem, new_stem; |
| 223 | 223 |
Value delta; |
| 224 |
|
|
| 225 |
const Value MAX; |
|
| 224 | 226 |
|
| 225 | 227 |
public: |
| 226 | 228 |
|
| 227 | 229 |
/// \brief Constant for infinite upper bounds (capacities). |
| ... | ... |
@@ -241,9 +243,9 @@ |
| 241 | 243 |
// References to the NetworkSimplex class |
| 242 | 244 |
const IntVector &_source; |
| 243 | 245 |
const IntVector &_target; |
| 244 | 246 |
const CostVector &_cost; |
| 245 |
const |
|
| 247 |
const CharVector &_state; |
|
| 246 | 248 |
const CostVector &_pi; |
| 247 | 249 |
int &_in_arc; |
| 248 | 250 |
int _search_arc_num; |
| 249 | 251 |
|
| ... | ... |
@@ -293,9 +295,9 @@ |
| 293 | 295 |
// References to the NetworkSimplex class |
| 294 | 296 |
const IntVector &_source; |
| 295 | 297 |
const IntVector &_target; |
| 296 | 298 |
const CostVector &_cost; |
| 297 |
const |
|
| 299 |
const CharVector &_state; |
|
| 298 | 300 |
const CostVector &_pi; |
| 299 | 301 |
int &_in_arc; |
| 300 | 302 |
int _search_arc_num; |
| 301 | 303 |
|
| ... | ... |
@@ -332,9 +334,9 @@ |
| 332 | 334 |
// References to the NetworkSimplex class |
| 333 | 335 |
const IntVector &_source; |
| 334 | 336 |
const IntVector &_target; |
| 335 | 337 |
const CostVector &_cost; |
| 336 |
const |
|
| 338 |
const CharVector &_state; |
|
| 337 | 339 |
const CostVector &_pi; |
| 338 | 340 |
int &_in_arc; |
| 339 | 341 |
int _search_arc_num; |
| 340 | 342 |
|
| ... | ... |
@@ -405,9 +407,9 @@ |
| 405 | 407 |
// References to the NetworkSimplex class |
| 406 | 408 |
const IntVector &_source; |
| 407 | 409 |
const IntVector &_target; |
| 408 | 410 |
const CostVector &_cost; |
| 409 |
const |
|
| 411 |
const CharVector &_state; |
|
| 410 | 412 |
const CostVector &_pi; |
| 411 | 413 |
int &_in_arc; |
| 412 | 414 |
int _search_arc_num; |
| 413 | 415 |
|
| ... | ... |
@@ -508,9 +510,9 @@ |
| 508 | 510 |
// References to the NetworkSimplex class |
| 509 | 511 |
const IntVector &_source; |
| 510 | 512 |
const IntVector &_target; |
| 511 | 513 |
const CostVector &_cost; |
| 512 |
const |
|
| 514 |
const CharVector &_state; |
|
| 513 | 515 |
const CostVector &_pi; |
| 514 | 516 |
int &_in_arc; |
| 515 | 517 |
int _search_arc_num; |
| 516 | 518 |
|
| ... | ... |
@@ -630,11 +632,11 @@ |
| 630 | 632 |
/// In special cases, it could lead to better overall performance, |
| 631 | 633 |
/// but it is usually slower. Therefore it is disabled by default. |
| 632 | 634 |
NetworkSimplex(const GR& graph, bool arc_mixing = false) : |
| 633 | 635 |
_graph(graph), _node_id(graph), _arc_id(graph), |
| 636 |
MAX(std::numeric_limits<Value>::max()), |
|
| 634 | 637 |
INF(std::numeric_limits<Value>::has_infinity ? |
| 635 |
std::numeric_limits<Value>::infinity() : |
|
| 636 |
std::numeric_limits<Value>::max()) |
|
| 638 |
std::numeric_limits<Value>::infinity() : MAX) |
|
| 637 | 639 |
{
|
| 638 | 640 |
// Check the value types |
| 639 | 641 |
LEMON_ASSERT(std::numeric_limits<Value>::is_signed, |
| 640 | 642 |
"The flow type of NetworkSimplex must be signed"); |
| ... | ... |
@@ -1019,11 +1021,11 @@ |
| 1019 | 1021 |
if (_have_lower) {
|
| 1020 | 1022 |
for (int i = 0; i != _arc_num; ++i) {
|
| 1021 | 1023 |
Value c = _lower[i]; |
| 1022 | 1024 |
if (c >= 0) {
|
| 1023 |
_cap[i] = _upper[i] < |
|
| 1025 |
_cap[i] = _upper[i] < MAX ? _upper[i] - c : INF; |
|
| 1024 | 1026 |
} else {
|
| 1025 |
_cap[i] = _upper[i] < |
|
| 1027 |
_cap[i] = _upper[i] < MAX + c ? _upper[i] - c : INF; |
|
| 1026 | 1028 |
} |
| 1027 | 1029 |
_supply[_source[i]] -= c; |
| 1028 | 1030 |
_supply[_target[i]] += c; |
| 1029 | 1031 |
} |
| ... | ... |
@@ -1213,9 +1215,9 @@ |
| 1213 | 1215 |
// Search the cycle along the path form the first node to the root |
| 1214 | 1216 |
for (int u = first; u != join; u = _parent[u]) {
|
| 1215 | 1217 |
e = _pred[u]; |
| 1216 | 1218 |
d = _forward[u] ? |
| 1217 |
_flow[e] : (_cap[e] |
|
| 1219 |
_flow[e] : (_cap[e] >= MAX ? INF : _cap[e] - _flow[e]); |
|
| 1218 | 1220 |
if (d < delta) {
|
| 1219 | 1221 |
delta = d; |
| 1220 | 1222 |
u_out = u; |
| 1221 | 1223 |
result = 1; |
| ... | ... |
@@ -1224,9 +1226,9 @@ |
| 1224 | 1226 |
// Search the cycle along the path form the second node to the root |
| 1225 | 1227 |
for (int u = second; u != join; u = _parent[u]) {
|
| 1226 | 1228 |
e = _pred[u]; |
| 1227 | 1229 |
d = _forward[u] ? |
| 1228 |
(_cap[e] |
|
| 1230 |
(_cap[e] >= MAX ? INF : _cap[e] - _flow[e]) : _flow[e]; |
|
| 1229 | 1231 |
if (d <= delta) {
|
| 1230 | 1232 |
delta = d; |
| 1231 | 1233 |
u_out = u; |
| 1232 | 1234 |
result = 2; |
| ... | ... |
@@ -1423,9 +1425,9 @@ |
| 1423 | 1425 |
// Execute the Network Simplex algorithm |
| 1424 | 1426 |
while (pivot.findEnteringArc()) {
|
| 1425 | 1427 |
findJoinNode(); |
| 1426 | 1428 |
bool change = findLeavingArc(); |
| 1427 |
if (delta >= |
|
| 1429 |
if (delta >= MAX) return UNBOUNDED; |
|
| 1428 | 1430 |
changeFlow(change); |
| 1429 | 1431 |
if (change) {
|
| 1430 | 1432 |
updateTreeStructure(); |
| 1431 | 1433 |
updatePotential(); |
0 comments (0 inline)