[Lemon-commits] Peter Kovacs: Small implementation improvements ...
Lemon HG
hg at lemon.cs.elte.hu
Mon Dec 14 06:17:45 CET 2009
details: http://lemon.cs.elte.hu/hg/lemon/rev/fe80a8145653
changeset: 877:fe80a8145653
user: Peter Kovacs <kpeter [at] inf.elte.hu>
date: Thu Nov 12 23:45:15 2009 +0100
description:
Small implementation improvements in MCF algorithms (#180)
- Handle max() as infinite value (not only infinity()).
- Better GEQ handling in CapacityScaling.
- Skip the unnecessary saturating operations in the first phase in
CapacityScaling.
- Use vector<char> instead of vector<bool> and vector<int> if it is
possible and it proved to be usually faster.
diffstat:
lemon/capacity_scaling.h | 79 ++++++++++++++++++++++-----------------
lemon/network_simplex.h | 32 ++++++++-------
2 files changed, 61 insertions(+), 50 deletions(-)
diffs (truncated from 308 to 300 lines):
diff --git a/lemon/capacity_scaling.h b/lemon/capacity_scaling.h
--- a/lemon/capacity_scaling.h
+++ b/lemon/capacity_scaling.h
@@ -133,7 +133,7 @@
TEMPLATE_DIGRAPH_TYPEDEFS(GR);
typedef std::vector<int> IntVector;
- typedef std::vector<bool> BoolVector;
+ typedef std::vector<char> BoolVector;
typedef std::vector<Value> ValueVector;
typedef std::vector<Cost> CostVector;
@@ -196,6 +196,7 @@
private:
int _node_num;
+ bool _geq;
const IntVector &_first_out;
const IntVector &_target;
const CostVector &_cost;
@@ -210,10 +211,10 @@
public:
ResidualDijkstra(CapacityScaling& cs) :
- _node_num(cs._node_num), _first_out(cs._first_out),
- _target(cs._target), _cost(cs._cost), _res_cap(cs._res_cap),
- _excess(cs._excess), _pi(cs._pi), _pred(cs._pred),
- _dist(cs._node_num)
+ _node_num(cs._node_num), _geq(cs._sum_supply < 0),
+ _first_out(cs._first_out), _target(cs._target), _cost(cs._cost),
+ _res_cap(cs._res_cap), _excess(cs._excess), _pi(cs._pi),
+ _pred(cs._pred), _dist(cs._node_num)
{}
int run(int s, Value delta = 1) {
@@ -232,7 +233,8 @@
heap.pop();
// Traverse outgoing residual arcs
- for (int a = _first_out[u]; a != _first_out[u+1]; ++a) {
+ int last_out = _geq ? _first_out[u+1] : _first_out[u+1] - 1;
+ for (int a = _first_out[u]; a != last_out; ++a) {
if (_res_cap[a] < delta) continue;
v = _target[a];
switch (heap.state(v)) {
@@ -687,22 +689,25 @@
}
if (_sum_supply > 0) return INFEASIBLE;
- // Initialize maps
+ // Initialize vectors
for (int i = 0; i != _root; ++i) {
_pi[i] = 0;
_excess[i] = _supply[i];
}
// Remove non-zero lower bounds
+ const Value MAX = std::numeric_limits<Value>::max();
+ int last_out;
if (_have_lower) {
for (int i = 0; i != _root; ++i) {
- for (int j = _first_out[i]; j != _first_out[i+1]; ++j) {
+ last_out = _first_out[i+1];
+ for (int j = _first_out[i]; j != last_out; ++j) {
if (_forward[j]) {
Value c = _lower[j];
if (c >= 0) {
- _res_cap[j] = _upper[j] < INF ? _upper[j] - c : INF;
+ _res_cap[j] = _upper[j] < MAX ? _upper[j] - c : INF;
} else {
- _res_cap[j] = _upper[j] < INF + c ? _upper[j] - c : INF;
+ _res_cap[j] = _upper[j] < MAX + c ? _upper[j] - c : INF;
}
_excess[i] -= c;
_excess[_target[j]] += c;
@@ -718,15 +723,16 @@
}
// Handle negative costs
- for (int u = 0; u != _root; ++u) {
- for (int a = _first_out[u]; a != _first_out[u+1]; ++a) {
- Value rc = _res_cap[a];
- if (_cost[a] < 0 && rc > 0) {
- if (rc == INF) return UNBOUNDED;
- _excess[u] -= rc;
- _excess[_target[a]] += rc;
- _res_cap[a] = 0;
- _res_cap[_reverse[a]] += rc;
+ for (int i = 0; i != _root; ++i) {
+ last_out = _first_out[i+1] - 1;
+ for (int j = _first_out[i]; j != last_out; ++j) {
+ Value rc = _res_cap[j];
+ if (_cost[j] < 0 && rc > 0) {
+ if (rc >= MAX) return UNBOUNDED;
+ _excess[i] -= rc;
+ _excess[_target[j]] += rc;
+ _res_cap[j] = 0;
+ _res_cap[_reverse[j]] += rc;
}
}
}
@@ -736,24 +742,21 @@
_pi[_root] = 0;
_excess[_root] = -_sum_supply;
for (int a = _first_out[_root]; a != _res_arc_num; ++a) {
- int u = _target[a];
- if (_excess[u] < 0) {
- _res_cap[a] = -_excess[u] + 1;
- } else {
- _res_cap[a] = 1;
- }
- _res_cap[_reverse[a]] = 0;
+ int ra = _reverse[a];
+ _res_cap[a] = -_sum_supply + 1;
+ _res_cap[ra] = 0;
_cost[a] = 0;
- _cost[_reverse[a]] = 0;
+ _cost[ra] = 0;
}
} else {
_pi[_root] = 0;
_excess[_root] = 0;
for (int a = _first_out[_root]; a != _res_arc_num; ++a) {
+ int ra = _reverse[a];
_res_cap[a] = 1;
- _res_cap[_reverse[a]] = 0;
+ _res_cap[ra] = 0;
_cost[a] = 0;
- _cost[_reverse[a]] = 0;
+ _cost[ra] = 0;
}
}
@@ -762,8 +765,9 @@
// With scaling
Value max_sup = 0, max_dem = 0;
for (int i = 0; i != _node_num; ++i) {
- if ( _excess[i] > max_sup) max_sup = _excess[i];
- if (-_excess[i] > max_dem) max_dem = -_excess[i];
+ Value ex = _excess[i];
+ if ( ex > max_sup) max_sup = ex;
+ if (-ex > max_dem) max_dem = -ex;
}
Value max_cap = 0;
for (int j = 0; j != _res_arc_num; ++j) {
@@ -789,7 +793,8 @@
// Handle non-zero lower bounds
if (_have_lower) {
- for (int j = 0; j != _res_arc_num - _node_num + 1; ++j) {
+ int limit = _first_out[_root];
+ for (int j = 0; j != limit; ++j) {
if (!_forward[j]) _res_cap[j] += _lower[j];
}
}
@@ -812,8 +817,11 @@
ResidualDijkstra _dijkstra(*this);
while (true) {
// Saturate all arcs not satisfying the optimality condition
+ int last_out;
for (int u = 0; u != _node_num; ++u) {
- for (int a = _first_out[u]; a != _first_out[u+1]; ++a) {
+ last_out = _sum_supply < 0 ?
+ _first_out[u+1] : _first_out[u+1] - 1;
+ for (int a = _first_out[u]; a != last_out; ++a) {
int v = _target[a];
Cost c = _cost[a] + _pi[u] - _pi[v];
Value rc = _res_cap[a];
@@ -830,8 +838,9 @@
_excess_nodes.clear();
_deficit_nodes.clear();
for (int u = 0; u != _node_num; ++u) {
- if (_excess[u] >= _delta) _excess_nodes.push_back(u);
- if (_excess[u] <= -_delta) _deficit_nodes.push_back(u);
+ Value ex = _excess[u];
+ if (ex >= _delta) _excess_nodes.push_back(u);
+ if (ex <= -_delta) _deficit_nodes.push_back(u);
}
int next_node = 0, next_def_node = 0;
diff --git a/lemon/network_simplex.h b/lemon/network_simplex.h
--- a/lemon/network_simplex.h
+++ b/lemon/network_simplex.h
@@ -164,7 +164,7 @@
TEMPLATE_DIGRAPH_TYPEDEFS(GR);
typedef std::vector<int> IntVector;
- typedef std::vector<bool> BoolVector;
+ typedef std::vector<char> CharVector;
typedef std::vector<Value> ValueVector;
typedef std::vector<Cost> CostVector;
@@ -212,8 +212,8 @@
IntVector _succ_num;
IntVector _last_succ;
IntVector _dirty_revs;
- BoolVector _forward;
- IntVector _state;
+ CharVector _forward;
+ CharVector _state;
int _root;
// Temporary data used in the current pivot iteration
@@ -221,6 +221,8 @@
int first, second, right, last;
int stem, par_stem, new_stem;
Value delta;
+
+ const Value MAX;
public:
@@ -242,7 +244,7 @@
const IntVector &_source;
const IntVector &_target;
const CostVector &_cost;
- const IntVector &_state;
+ const CharVector &_state;
const CostVector &_pi;
int &_in_arc;
int _search_arc_num;
@@ -294,7 +296,7 @@
const IntVector &_source;
const IntVector &_target;
const CostVector &_cost;
- const IntVector &_state;
+ const CharVector &_state;
const CostVector &_pi;
int &_in_arc;
int _search_arc_num;
@@ -333,7 +335,7 @@
const IntVector &_source;
const IntVector &_target;
const CostVector &_cost;
- const IntVector &_state;
+ const CharVector &_state;
const CostVector &_pi;
int &_in_arc;
int _search_arc_num;
@@ -406,7 +408,7 @@
const IntVector &_source;
const IntVector &_target;
const CostVector &_cost;
- const IntVector &_state;
+ const CharVector &_state;
const CostVector &_pi;
int &_in_arc;
int _search_arc_num;
@@ -509,7 +511,7 @@
const IntVector &_source;
const IntVector &_target;
const CostVector &_cost;
- const IntVector &_state;
+ const CharVector &_state;
const CostVector &_pi;
int &_in_arc;
int _search_arc_num;
@@ -631,9 +633,9 @@
/// but it is usually slower. Therefore it is disabled by default.
NetworkSimplex(const GR& graph, bool arc_mixing = false) :
_graph(graph), _node_id(graph), _arc_id(graph),
+ MAX(std::numeric_limits<Value>::max()),
INF(std::numeric_limits<Value>::has_infinity ?
- std::numeric_limits<Value>::infinity() :
- std::numeric_limits<Value>::max())
+ std::numeric_limits<Value>::infinity() : MAX)
{
// Check the value types
LEMON_ASSERT(std::numeric_limits<Value>::is_signed,
@@ -1020,9 +1022,9 @@
for (int i = 0; i != _arc_num; ++i) {
Value c = _lower[i];
if (c >= 0) {
- _cap[i] = _upper[i] < INF ? _upper[i] - c : INF;
+ _cap[i] = _upper[i] < MAX ? _upper[i] - c : INF;
} else {
- _cap[i] = _upper[i] < INF + c ? _upper[i] - c : INF;
+ _cap[i] = _upper[i] < MAX + c ? _upper[i] - c : INF;
}
_supply[_source[i]] -= c;
_supply[_target[i]] += c;
@@ -1214,7 +1216,7 @@
for (int u = first; u != join; u = _parent[u]) {
e = _pred[u];
d = _forward[u] ?
- _flow[e] : (_cap[e] == INF ? INF : _cap[e] - _flow[e]);
+ _flow[e] : (_cap[e] >= MAX ? INF : _cap[e] - _flow[e]);
if (d < delta) {
delta = d;
u_out = u;
@@ -1225,7 +1227,7 @@
for (int u = second; u != join; u = _parent[u]) {
e = _pred[u];
d = _forward[u] ?
- (_cap[e] == INF ? INF : _cap[e] - _flow[e]) : _flow[e];
+ (_cap[e] >= MAX ? INF : _cap[e] - _flow[e]) : _flow[e];
if (d <= delta) {
delta = d;
u_out = u;
@@ -1424,7 +1426,7 @@
More information about the Lemon-commits
mailing list