Changeset 604:8c3112a66878 in lemon-main

Timestamp:

03/24/09 00:18:25 (16 years ago)

Author:

Peter Kovacs <kpeter@…>

Branch:

default

Phase:

public

Message:

Use XTI implementation instead of ATI in NetworkSimplex? (#234)

XTI (eXtended Threaded Index) is an imporved version of the widely
known ATI (Augmented Threaded Index) method for storing and updating
the spanning tree structure in Network Simplex algorithms.

In the ATI data structure three indices are stored for each node:
predecessor, thread and depth. In the XTI data structure depth is
replaced by the number of successors and the last successor
(according to the thread index).

File:

• lemon/network_simplex.h (modified) (7 diffs)

lemon/network_simplex.h

@@ -155 +155 @@
 
     // Data for storing the spanning tree structure
-    IntVector _depth;
     IntVector _parent;
     IntVector _pred;
     IntVector _thread;
+    IntVector _rev_thread;
+    IntVector _succ_num;
+    IntVector _last_succ;
+    IntVector _dirty_revs;
     BoolVector _forward;
     IntVector _state;
@@ -851 +854 @@
       _pi.resize(all_node_num, 0);
 
-      _depth.resize(all_node_num);
       _parent.resize(all_node_num);
       _pred.resize(all_node_num);
       _forward.resize(all_node_num);
       _thread.resize(all_node_num);
+      _rev_thread.resize(all_node_num);
+      _succ_num.resize(all_node_num);
+      _last_succ.resize(all_node_num);
       _state.resize(all_arc_num, STATE_LOWER);
 
@@ -882 +887 @@
       // Set data for the artificial root node
       _root = _node_num;
-      _depth[_root] = 0;
       _parent[_root] = -1;
       _pred[_root] = -1;
       _thread[_root] = 0;
+      _rev_thread[0] = _root;
+      _succ_num[_root] = all_node_num;
+      _last_succ[_root] = _root - 1;
       _supply[_root] = 0;
       _pi[_root] = 0;
@@ -923 +930 @@
       for (int u = 0, e = _arc_num; u != _node_num; ++u, ++e) {
         _thread[u] = u + 1;
-        _depth[u] = 1;
+        _rev_thread[u + 1] = u;
+        _succ_num[u] = 1;
+        _last_succ[u] = u;
         _parent[u] = _root;
         _pred[u] = e;
@@ -947 +956 @@
       int u = _source[in_arc];
       int v = _target[in_arc];
-      while (_depth[u] > _depth[v]) u = _parent[u];
-      while (_depth[v] > _depth[u]) v = _parent[v];
       while (u != v) {
-        u = _parent[u];
-        v = _parent[v];
+        if (_succ_num[u] < _succ_num[v]) {
+          u = _parent[u];
+        } else {
+          v = _parent[v];
+        }
       }
       join = u;
@@ -1027 +1037 @@
     }
 
-    // Update _thread and _parent vectors
-    void updateThreadParent() {
-      int u;
+    // Update the tree structure
+    void updateTreeStructure() {
+      int u, w;
+      int old_rev_thread = _rev_thread[u_out];
+      int old_succ_num = _succ_num[u_out];
+      int old_last_succ = _last_succ[u_out];
       v_out = _parent[u_out];
 
-      // Handle the case when join and v_out coincide
-      bool par_first = false;
-      if (join == v_out) {
-        for (u = join; u != u_in && u != v_in; u = _thread[u]) ;
-        if (u == v_in) {
-          par_first = true;
-          while (_thread[u] != u_out) u = _thread[u];
-          first = u;
-        }
-      }
-
-      // Find the last successor of u_in (u) and the node after it (right)
-      // according to the thread index
-      for (u = u_in; _depth[_thread[u]] > _depth[u_in]; u = _thread[u]) ;
-      right = _thread[u];
-      if (_thread[v_in] == u_out) {
-        for (last = u; _depth[last] > _depth[u_out]; last = _thread[last]) ;
-        if (last == u_out) last = _thread[last];
-      }
-      else last = _thread[v_in];
-
-      // Update stem nodes
+      u = _last_succ[u_in];  // the last successor of u_in
+      right = _thread[u];    // the node after it
+
+      // Handle the case when old_rev_thread equals to v_in
+      // (it also means that join and v_out coincide)
+      if (old_rev_thread == v_in) {
+        last = _thread[_last_succ[u_out]];
+      } else {
+        last = _thread[v_in];
+      }
+
+      // Update _thread and _parent along the stem nodes (i.e. the nodes
+      // between u_in and u_out, whose parent have to be changed)
       _thread[v_in] = stem = u_in;
+      _dirty_revs.clear();
+      _dirty_revs.push_back(v_in);
       par_stem = v_in;
       while (stem != u_out) {
-        _thread[u] = new_stem = _parent[stem];
-
-        // Find the node just before the stem node (u) according to
-        // the original thread index
-        for (u = new_stem; _thread[u] != stem; u = _thread[u]) ;
-        _thread[u] = right;
-
-        // Change the parent node of stem and shift stem and par_stem nodes
+        // Insert the next stem node into the thread list
+        new_stem = _parent[stem];
+        _thread[u] = new_stem;
+        _dirty_revs.push_back(u);
+
+        // Remove the subtree of stem from the thread list
+        w = _rev_thread[stem];
+        _thread[w] = right;
+        _rev_thread[right] = w;
+
+        // Change the parent node and shift stem nodes
         _parent[stem] = par_stem;
         par_stem = stem;
         stem = new_stem;
 
-        // Find the last successor of stem (u) and the node after it (right)
-        // according to the thread index
-        for (u = stem; _depth[_thread[u]] > _depth[stem]; u = _thread[u]) ;
+        // Update u and right
+        u = _last_succ[stem] == _last_succ[par_stem] ?
+          _rev_thread[par_stem] : _last_succ[stem];
         right = _thread[u];
       }
       _parent[u_out] = par_stem;
       _thread[u] = last;
-
-      if (join == v_out && par_first) {
-        if (first != v_in) _thread[first] = right;
-      } else {
-        for (u = v_out; _thread[u] != u_out; u = _thread[u]) ;
-        _thread[u] = right;
-      }
-    }
-
-    // Update _pred and _forward vectors
-    void updatePredArc() {
-      int u = u_out, v;
+      _rev_thread[last] = u;
+      _last_succ[u_out] = u;
+
+      // Remove the subtree of u_out from the thread list except for
+      // the case when old_rev_thread equals to v_in
+      // (it also means that join and v_out coincide)
+      if (old_rev_thread != v_in) {
+        _thread[old_rev_thread] = right;
+        _rev_thread[right] = old_rev_thread;
+      }
+
+      // Update _rev_thread using the new _thread values
+      for (int i = 0; i < int(_dirty_revs.size()); ++i) {
+        u = _dirty_revs[i];
+        _rev_thread[_thread[u]] = u;
+      }
+
+      // Update _pred, _forward, _last_succ and _succ_num for the
+      // stem nodes from u_out to u_in
+      int tmp_sc = 0, tmp_ls = _last_succ[u_out];
+      u = u_out;
       while (u != u_in) {
-        v = _parent[u];
-        _pred[u] = _pred[v];
-        _forward[u] = !_forward[v];
-        u = v;
+        w = _parent[u];
+        _pred[u] = _pred[w];
+        _forward[u] = !_forward[w];
+        tmp_sc += _succ_num[u] - _succ_num[w];
+        _succ_num[u] = tmp_sc;
+        _last_succ[w] = tmp_ls;
+        u = w;
       }
       _pred[u_in] = in_arc;
       _forward[u_in] = (u_in == _source[in_arc]);
-    }
-
-    // Update _depth and _potential vectors
-    void updateDepthPotential() {
-      _depth[u_in] = _depth[v_in] + 1;
+      _succ_num[u_in] = old_succ_num;
+
+      // Set limits for updating _last_succ form v_in and v_out
+      // towards the root
+      int up_limit_in = -1;
+      int up_limit_out = -1;
+      if (_last_succ[join] == v_in) {
+        up_limit_out = join;
+      } else {
+        up_limit_in = join;
+      }
+
+      // Update _last_succ from v_in towards the root
+      for (u = v_in; u != up_limit_in && _last_succ[u] == v_in;
+           u = _parent[u]) {
+        _last_succ[u] = _last_succ[u_out];
+      }
+      // Update _last_succ from v_out towards the root
+      if (join != old_rev_thread && v_in != old_rev_thread) {
+        for (u = v_out; u != up_limit_out && _last_succ[u] == old_last_succ;
+             u = _parent[u]) {
+          _last_succ[u] = old_rev_thread;
+        }
+      } else {
+        for (u = v_out; u != up_limit_out && _last_succ[u] == old_last_succ;
+             u = _parent[u]) {
+          _last_succ[u] = _last_succ[u_out];
+        }
+      }
+
+      // Update _succ_num from v_in to join
+      for (u = v_in; u != join; u = _parent[u]) {
+        _succ_num[u] += old_succ_num;
+      }
+      // Update _succ_num from v_out to join
+      for (u = v_out; u != join; u = _parent[u]) {
+        _succ_num[u] -= old_succ_num;
+      }
+    }
+
+    // Update potentials
+    void updatePotential() {
       Cost sigma = _forward[u_in] ?
         _pi[v_in] - _pi[u_in] - _cost[_pred[u_in]] :
         _pi[v_in] - _pi[u_in] + _cost[_pred[u_in]];
-      _pi[u_in] += sigma;
-      for(int u = _thread[u_in]; _parent[u] != -1; u = _thread[u]) {
-        _depth[u] = _depth[_parent[u]] + 1;
-        if (_depth[u] <= _depth[u_in]) break;
-        _pi[u] += sigma;
+      if (_succ_num[u_in] > _node_num / 2) {
+        // Update in the upper subtree (which contains the root)
+        int before = _rev_thread[u_in];
+        int after = _thread[_last_succ[u_in]];
+        _thread[before] = after;
+        _pi[_root] -= sigma;
+        for (int u = _thread[_root]; u != _root; u = _thread[u]) {
+          _pi[u] -= sigma;
+        }
+        _thread[before] = u_in;
+      } else {
+        // Update in the lower subtree (which has been moved)
+        int end = _thread[_last_succ[u_in]];
+        for (int u = u_in; u != end; u = _thread[u]) {
+          _pi[u] += sigma;
+        }
       }
     }
@@ -1140 +1209 @@
         changeFlow(change);
         if (change) {
-          updateThreadParent();
-          updatePredArc();
-          updateDepthPotential();
+          updateTreeStructure();
+          updatePotential();
         }
       }