LEMON/LEMON-main Changeset - r811:fe80a8145653

Repositories » LEMON » LEMON-main

Summary
Changelog
Files
Switch To
- loading...
Options
- Lightweight changelog
- Search
Pull Requests

Changeset - r811:fe80a8145653

« 810:3b5349

812:4b1b37 »

r811:fe80a8145653

Thu, 12 Nov 2009 23:45:15

[Not Reviewed]

0 comments (0 inline)

kpeter (Peter Kovacs)
kpeter@inf.elte.hu

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.

0 2 0

default

2 files changed with 61 insertions and 50 deletions:

lemon/capacity_scaling.h

lemon/network_simplex.h

↑ Collapse diff ↑

lemon/capacity_scaling.h

Show inline comments

@@ -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;

lemon/network_simplex.h

Show inline comments

@@ -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 @@
 		      while (pivot.findEnteringArc()) {
 		        findJoinNode();
 		        bool change = findLeavingArc();
-		        if (delta >= INF) return UNBOUNDED;
+		        if (delta >= MAX) return UNBOUNDED;
 		        changeFlow(change);
 		        if (change) {
 		          updateTreeStructure();

0 comments (0 inline)

RhodeCode

Login to your account