Merge branch 'opt' of https://git01.codeplex.com/z3 into opt

src/smt/network_flow.h

@@ -51,7 +51,7 @@ namespace smt {
         typedef typename Ext::fin_numeral fin_numeral;
 
         graph m_graph;
-        thread_spanning_tree<Ext> tree;
+        spanning_tree_base m_tree;
 
         // Denote supply/demand b_i on node i
         vector<fin_numeral> m_balances;
@@ -66,17 +66,13 @@ namespace smt {
 
         unsigned m_step;
 
-        edge_id m_entering_edge;
+        edge_id m_enter_id, m_leave_id;
-        edge_id m_leaving_edge;
-        node m_join_node;
         optional<numeral> m_delta;
-        bool m_in_edge_dir;
+        bool m_is_swap_enter, m_is_swap_leave;
 
         // Initialize the network with a feasible spanning tree
         void initialize();
 
-        edge_id get_edge_id(dl_var source, dl_var target) const;
-
         void update_potentials();
 
         void update_flows();
@@ -94,7 +90,6 @@ namespace smt {
         std::string display_spanning_tree();
 
         bool edge_in_tree(edge_id id) const;
-        bool edge_in_tree(node src, node dst) const;
 
         bool check_well_formed();
         bool check_optimal();

src/smt/network_flow_def.h

@@ -40,30 +40,26 @@ namespace smt {
                 m_graph.add_edge(e.get_target(), e.get_source(), e.get_weight(), explanation());
             }
         }
-        
+        m_tree = thread_spanning_tree<Ext>(m_graph);
-        unsigned num_nodes = m_graph.get_num_nodes() + 1;
-        unsigned num_edges = m_graph.get_num_edges();
-
-        m_balances.resize(num_nodes);
-        m_potentials.resize(num_nodes);      
-
-        tree = thread_spanning_tree<Ext>();
         m_step = 0;
     }
 
     template<typename Ext>
     void network_flow<Ext>::initialize() {
         TRACE("network_flow", tout << "initialize...\n";);
-        // Create an artificial root node to construct initial spanning tree
+        // Create an artificial root node to construct initial spanning m_tree
         unsigned num_nodes = m_graph.get_num_nodes();
         unsigned num_edges = m_graph.get_num_edges();
         
         node root = num_nodes;
         m_graph.init_var(root);
+
+        m_potentials.resize(num_nodes + 1);     
         m_potentials[root] = numeral::zero();
 
+        m_balances.resize(num_nodes + 1);
         fin_numeral sum_supply = fin_numeral::zero();
-        for (unsigned i = 0; i < m_balances.size(); ++i) {
+        for (unsigned i = 0; i < num_nodes; ++i) {
             sum_supply += m_balances[i];
         }        
         m_balances[root] = -sum_supply;
@@ -73,54 +69,37 @@ namespace smt {
         m_states.resize(num_nodes + num_edges);
         m_states.fill(LOWER);
 
-        // Create artificial edges from/to root node to/from other nodes and initialize the spanning tree
+        // Create artificial edges from/to root node to/from other nodes and initialize the spanning m_tree
         svector<bool> upwards(num_nodes, false);
         for (unsigned i = 0; i < num_nodes; ++i) {
             upwards[i] = !m_balances[i].is_neg();
             m_states[num_edges + i] = BASIS;
             node src = upwards[i] ? i : root;
             node tgt = upwards[i] ? root : i;            
-            m_flows[num_edges + i] = upwards[i] ? m_balances[i] : -m_balances[i];            
+            m_flows[num_edges + i] = upwards[i] ? m_balances[i] : -m_balances[i];
+            m_potentials[i] = upwards[i] ? numeral::one() : -numeral::one();
             m_graph.add_edge(src, tgt, numeral::one(), explanation());
         }
 
-        tree.initialize(upwards, num_nodes);
+        m_tree.initialize(upwards);
-
-        // Compute initial potentials
-        svector<node> descendants;
-        tree.get_descendants(root, descendants);
-        // Skip root node
-        for (unsigned i = 1; i < descendants.size(); ++i) {
-            node u = descendants[i];
-            node v = tree.get_parent(u);   
-            edge_id e_id = get_edge_id(u, v);            
-            m_potentials[u] = m_potentials[v] + (tree.get_arc_direction(u) ? - m_graph.get_weight(e_id) : m_graph.get_weight(e_id));
-        }
 
         TRACE("network_flow", {
                 tout << pp_vector("Potentials", m_potentials, true) << pp_vector("Flows", m_flows);
             });
-        TRACE("network_flow", tout << "Spanning tree:\n" << display_spanning_tree(););
+        TRACE("network_flow", tout << "Spanning m_tree:\n" << display_spanning_tree(););
         SASSERT(check_well_formed());
     }
 
     template<typename Ext>
-    edge_id network_flow<Ext>::get_edge_id(dl_var source, dl_var target) const {
+    void network_flow<Ext>::update_potentials() {        
-        // tree.get_arc_direction(source) decides which node is the real source
+        node src = m_graph.get_source(m_enter_id);
-        edge_id id;
+        node tgt = m_graph.get_target(m_enter_id); 
-        VERIFY(tree.get_arc_direction(source) ? m_graph.get_edge_id(source, target, id) : m_graph.get_edge_id(target, source, id));
+        numeral cost = m_potentials[src] - m_potentials[tgt] - m_graph.get_weight(m_enter_id);
-        return id;
+        numeral change = m_is_swap_leave ? -cost : cost;
-    }
-
-    template<typename Ext>
-    void network_flow<Ext>::update_potentials() {
-        TRACE("network_flow", tout << "update_potentials...\n";);
-        node src = m_graph.get_source(m_entering_edge);
-        node tgt = m_graph.get_target(m_entering_edge); 
-        numeral cost = m_graph.get_weight(m_entering_edge);
-        numeral change = m_potentials[tgt] - m_potentials[src] + (tree.get_arc_direction(src) ? -cost : cost);
         svector<node> descendants;
-        tree.get_descendants(src, descendants);
+        node start = m_is_swap_enter ? src : tgt;
+        TRACE("network_flow", tout << "update_potentials of T_" << start << " with delta = " << change << "...\n";);
+        m_tree.get_descendants(start, descendants);
         for (unsigned i = 0; i < descendants.size(); ++i) {
             node u = descendants[i];
             m_potentials[u] += change;
@@ -131,23 +110,16 @@ namespace smt {
     template<typename Ext>
     void network_flow<Ext>::update_flows() {
         TRACE("network_flow", tout << "update_flows...\n";);
-        numeral val = fin_numeral(m_states[m_entering_edge]) * (*m_delta);
+        numeral val = *m_delta;
-        m_flows[m_entering_edge] += val;
+        m_flows[m_enter_id] += val;
-        node source = m_graph.get_source(m_entering_edge);
+        node src = m_graph.get_source(m_enter_id);
-        svector<node> ancestors;
+        node tgt = m_graph.get_target(m_enter_id); 
-        tree.get_ancestors(source, ancestors);
+        svector<edge_id> path;
-        for (unsigned i = 0; i < ancestors.size() && ancestors[i] != m_join_node; ++i) {
+        svector<bool> against;
-            node u = ancestors[i];
+        m_tree.get_path(src, tgt, path, against);
-            edge_id e_id = get_edge_id(u, tree.get_parent(u));
+        for (unsigned i = 0; i < path.size(); ++i) {
-            m_flows[e_id] += tree.get_arc_direction(u) ? -val : val;
+            edge_id e_id = path[i];
-        }
+            m_flows[e_id] += against[i] ? -val : val;
-
-        node target = m_graph.get_target(m_entering_edge);
-        tree.get_ancestors(target, ancestors);
-        for (unsigned i = 0; i < ancestors.size() && ancestors[i] != m_join_node; ++i) {
-            node u = ancestors[i];
-            edge_id e_id = get_edge_id(u, tree.get_parent(u));
-            m_flows[e_id] += tree.get_arc_direction(u) ? val : -val;
         }
         TRACE("network_flow", tout << pp_vector("Flows", m_flows, true););
     }
@@ -160,11 +132,10 @@ namespace smt {
             node src = m_graph.get_source(i);
             node tgt = m_graph.get_target(i);
             if (m_states[i] != BASIS) {
-                numeral change = fin_numeral(m_states[i]) * (m_graph.get_weight(i) + m_potentials[src] - m_potentials[tgt]);
+                numeral cost = m_potentials[src] - m_potentials[tgt] - m_graph.get_weight(i);
-                // Choose the first negative-cost edge to be the violating edge
                 // TODO: add multiple pivoting strategies
-                if (change.is_neg()) {
+                if (cost.is_pos()) {
-                    m_entering_edge = i;
+                    m_enter_id = i;
                     TRACE("network_flow", {
                         tout << "Found entering edge " << i << " between node ";
                         tout << src << " and node " << tgt << "...\n";
@@ -180,52 +151,28 @@ namespace smt {
     template<typename Ext>
     bool network_flow<Ext>::choose_leaving_edge() {
         TRACE("network_flow", tout << "choose_leaving_edge...\n";);
-        node source = m_graph.get_source(m_entering_edge);
+        node src = m_graph.get_source(m_enter_id);
-        node target = m_graph.get_target(m_entering_edge);
+        node tgt = m_graph.get_target(m_enter_id); 
-        if (m_states[m_entering_edge] == UPPER) {
-            std::swap(source, target);
-        }
-        
-        m_join_node = tree.get_common_ancestor(source, target);
-
-        TRACE("network_flow", tout << "Found join node " << m_join_node << std::endl;);    
         m_delta.set_invalid();
-        node src, tgt;
+        edge_id leave_id;
-
+        svector<edge_id> path;
-        // Send flows along the path from source to the ancestor
+        svector<bool> against;
-        svector<node> ancestors;
+        m_tree.get_path(src, tgt, path, against);
-        tree.get_ancestors(source, ancestors);
+        for (unsigned i = 0; i < path.size(); ++i) {
-        for (unsigned i = 0; i < ancestors.size() && ancestors[i] != m_join_node; ++i) {
+            edge_id e_id = path[i];
-            node u = ancestors[i];
+            if (against[i] && (!m_delta || m_flows[e_id] <= *m_delta)) {
-            edge_id e_id = get_edge_id(u, tree.get_parent(u));            
-            if (tree.get_arc_direction(u) && (!m_delta || m_flows[e_id] < *m_delta)) {
                 m_delta = m_flows[e_id];
-                src = u;
+                leave_id = e_id;
-                tgt = tree.get_parent(u);
-                SASSERT(edge_in_tree(src,tgt));
-                m_in_edge_dir = true;
-            }
-        }
-
-        // Send flows along the path from target to the ancestor
-        tree.get_ancestors(target, ancestors);
-        for (unsigned i = 0; i < ancestors.size() && ancestors[i] != m_join_node; ++i) {
-            node u = ancestors[i];
-            edge_id e_id = get_edge_id(u, tree.get_parent(u)); 
-            if (!tree.get_arc_direction(u) && (!m_delta || m_flows[e_id] <= *m_delta)) {
-                m_delta = m_flows[e_id];
-                src = u;
-                tgt = tree.get_parent(u);
-                SASSERT(edge_in_tree(src,tgt));
-                m_in_edge_dir = false;
             }
         }
 
         if (m_delta) {
-            m_leaving_edge = get_edge_id(src, tgt);
+            m_leave_id = leave_id;
+
             TRACE("network_flow", { 
-                tout << "Found leaving edge " << m_leaving_edge;
+                tout << "Found leaving edge " << m_leave_id;
-                tout << " between node " << src << " and node " << tgt << "...\n";
+                tout << " between node " << m_graph.get_source(m_leave_id);
+                tout << " and node " << m_graph.get_target(m_leave_id) << " with delta = " << *m_delta << "...\n";
                 });
             return true;
         }
@@ -234,13 +181,8 @@ namespace smt {
     }
 
     template<typename Ext>
-    void network_flow<Ext>::update_spanning_tree() {        
+    void network_flow<Ext>::update_spanning_tree() {  
-        node p = m_graph.get_source(m_entering_edge);
+        m_tree.update(m_enter_id, m_leave_id, m_is_swap_enter, m_is_swap_leave);
-        node q = m_graph.get_target(m_entering_edge);        
-        node u = m_graph.get_source(m_leaving_edge);
-        node v = m_graph.get_target(m_leaving_edge);
-       
-        tree.update(p, q, u, v);
     }
 
     // Minimize cost flows
@@ -252,18 +194,18 @@ namespace smt {
             bool bounded = choose_leaving_edge();
             if (!bounded) return false;
             update_flows();
-            if (m_entering_edge != m_leaving_edge) {
+            if (m_enter_id != m_leave_id) {
-				SASSERT(edge_in_tree(m_leaving_edge));
+				SASSERT(edge_in_tree(m_leave_id));
-				SASSERT(!edge_in_tree(m_entering_edge));
+				SASSERT(!edge_in_tree(m_enter_id));
-                m_states[m_entering_edge] = BASIS;
+                m_states[m_enter_id] = BASIS;
-                m_states[m_leaving_edge] = (m_flows[m_leaving_edge].is_zero()) ? LOWER : UPPER;
+                m_states[m_leave_id] = (m_flows[m_leave_id].is_zero()) ? LOWER : UPPER;
                 update_spanning_tree();
                 update_potentials();                
-                TRACE("network_flow", tout << "Spanning tree:\n" << display_spanning_tree(););
+                TRACE("network_flow", tout << "Spanning m_tree:\n" << display_spanning_tree(););
                 SASSERT(check_well_formed());
             } 
             else {
-                m_states[m_leaving_edge] = m_states[m_leaving_edge] == LOWER ? UPPER : LOWER;
+                m_states[m_leave_id] = m_states[m_leave_id] == LOWER ? UPPER : LOWER;
             }            
         }
         TRACE("network_flow", tout << "Found optimal solution.\n";);
@@ -297,27 +239,24 @@ namespace smt {
     bool network_flow<Ext>::edge_in_tree(edge_id id) const {
         return m_states[id] == BASIS;
     }
-
-    template<typename Ext>
-    bool network_flow<Ext>::edge_in_tree(node src, node dst) const {
-        return edge_in_tree(get_edge_id(src, dst));
-    }
     
     template<typename Ext>
     bool network_flow<Ext>::check_well_formed() {
-        SASSERT(tree.check_well_formed());
+        SASSERT(m_tree.check_well_formed());
 
         // m_flows are zero on non-basic edges
         for (unsigned i = 0; i < m_flows.size(); ++i) {
+            SASSERT(!m_flows[i].is_neg());
             SASSERT(m_states[i] == BASIS || m_flows[i].is_zero());
         }
 
-        // m_upwards show correct direction
+        vector<edge> const & es = m_graph.get_all_edges();
-        for (unsigned i = 0; i < m_potentials.size(); ++i) {
+        for (unsigned i = 0; i < es.size(); ++i) {
-            node p = tree.get_parent(i);
+            edge const & e = es[i];
-            edge_id id;
+            if (m_states[i] == BASIS) {
-            SASSERT(!tree.get_arc_direction(i) || m_graph.get_edge_id(i, p, id));            
+                SASSERT(m_potentials[e.get_source()] - m_potentials[e.get_target()] == e.get_weight());
-        }    
+            }
+        }
 
         return true;
     }
@@ -328,8 +267,7 @@ namespace smt {
         vector<edge> const & es = m_graph.get_all_edges();
         for (unsigned i = 0; i < es.size(); ++i) {
             edge const & e = es[i];
-            if (m_states[i] == BASIS) 
+            if (m_states[i] == BASIS) {
-            {
                 total_cost += e.get_weight().get_rational() * m_flows[i];
             }
         }

src/smt/spanning_tree.h

@@ -44,25 +44,27 @@ namespace smt {
         // (i, m_pred[i]) points upwards (pointing toward the root node)
         svector<bool> m_upwards;
 
+        graph & m_graph;
+
         void swap_order(node q, node v);
         node find_rev_thread(node n) const;
         void fix_depth(node start, node end);
         node get_final(int start);
-        bool is_preorder_traversal(node start, node end);
+        bool is_preorder_traversal(node start, node end);        
-        bool is_ancestor_of(node ancestor, node child);
+        edge_id get_edge_to_parent(node start) const;
+        node get_common_ancestor(node u, node v);
 
     public:
+        thread_spanning_tree(graph & g);
 
-        void initialize(svector<bool> const & upwards, int num_nodes);
+        void initialize(svector<bool> const & upwards);
         void get_descendants(node start, svector<node> & descendants);
-        void get_ancestors(node start, svector<node> & ancestors);
+        
-        node get_common_ancestor(node u, node v);
+        void update(edge_id enter_id, edge_id leave_id, bool & is_swap_enter, bool & is_swap_leave);
-        void update(node p, node q, node u, node v);
         bool check_well_formed();
-
+        void get_path(node start, node end, svector<edge_id> & path, svector<bool> & against);
-        // TODO: remove these two unnatural functions
+        bool is_forward_edge(edge_id e_id) const;
-        bool get_arc_direction(node start) const;
+        bool is_ancestor_of(node ancestor, node child);
-        node get_parent(node start);
     };
 
 }

src/smt/spanning_tree_base.h

@@ -47,24 +47,14 @@ namespace smt {
         typedef int node;
 
     public:
-        virtual void initialize(svector<bool> const & upwards, int num_nodes) {};
+        virtual void initialize(svector<bool> const & upwards) {};
-
-        /**
-            \brief Get all descendants of a node including itself
-        */
         virtual void get_descendants(node start, svector<node> & descendants) {};
-        /**
+        
-            \brief Get all ancestors of a node including itself
+        virtual void update(edge_id enter_id, edge_id leave_id, bool & is_swap_enter, bool & is_swap_leave) {};
-        */
-        virtual void get_ancestors(node start, svector<node> & ancestors) {};
-
-        virtual node get_common_ancestor(node u, node v) {UNREACHABLE(); return -1;};
-        virtual void update(node p, node q, node u, node v) {};
         virtual bool check_well_formed() {UNREACHABLE(); return false;};
-
+        virtual void get_path(node start, node end, svector<edge_id> & path, svector<bool> & against) {};
-        // TODO: remove these two unnatural functions
+        virtual bool is_forward_edge(edge_id e_id) const {UNREACHABLE(); return false;};
-        virtual bool get_arc_direction(node start) const {UNREACHABLE(); return false;};
+        virtual bool is_ancestor_of(node ancestor, node child) {UNREACHABLE(); return false;};
-        virtual node get_parent(node start) {UNREACHABLE(); return -1;};
     };
 }
 

src/smt/spanning_tree_def.h

@@ -22,7 +22,280 @@ Notes:
 #include "spanning_tree.h"
 
 namespace smt {
+
+    template<typename Ext>
+    thread_spanning_tree<Ext>::thread_spanning_tree(graph & g) :
+        m_graph(g) {
+    }
+
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::initialize(svector<bool> const & upwards) {
+        unsigned num_nodes = m_graph.get_num_nodes();
+        m_pred.resize(num_nodes);
+        m_depth.resize(num_nodes);
+        m_thread.resize(num_nodes);
+        m_upwards.resize(num_nodes);
+
+        node root = m_graph.get_num_nodes() - 1;
+        m_pred[root] = -1;
+        m_depth[root] = 0;
+        m_thread[root] = 0;
+            
+        // Create artificial edges from/to root node to/from other nodes and initialize the spanning tree
+        for (int i = 0; i < root; ++i) {
+            m_pred[i] = root;
+            m_depth[i] = 1;
+            m_thread[i] = i + 1;
+            m_upwards[i] = upwards[i];
+        }
+
+        TRACE("network_flow", {
+            tout << pp_vector("Predecessors", m_pred, true) << pp_vector("Threads", m_thread); 
+            tout << pp_vector("Depths", m_depth) << pp_vector("Upwards", m_upwards);
+        });
+    }
+
+    template<typename Ext>
+    typename thread_spanning_tree<Ext>::node thread_spanning_tree<Ext>::get_common_ancestor(node u, node v) {
+        while (u != v) {
+            if (m_depth[u] > m_depth[v])
+                u = m_pred[u];
+            else 
+                v = m_pred[v];                    
+        }
+        return u;
+    }
+
+    template<typename Ext>
+    edge_id thread_spanning_tree<Ext>::get_edge_to_parent(node start) const {
+        SASSERT(m_pred[start] != -1);
+        edge_id id;
+        node end = m_pred[start];
+        VERIFY(m_upwards[start] ? m_graph.get_edge_id(start, end, id) : m_graph.get_edge_id(end, start, id));
+        return id;
+    }
+
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::get_path(node start, node end, svector<edge_id> & path, svector<bool> & against) {
+        node join = get_common_ancestor(start, end);
+        path.reset();
+        while (start != join) {
+            edge_id e_id = get_edge_to_parent(start);
+            path.push_back(e_id);
+            against.push_back(is_forward_edge(e_id));
+            start = m_pred[start];
+        }
+        while (end != join) {
+            edge_id e_id = get_edge_to_parent(end);
+            path.push_back(e_id);
+            against.push_back(!is_forward_edge(e_id));
+            end = m_pred[end];
+        }
+    }
+
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::is_forward_edge(edge_id e_id) const {
+        node start = m_graph.get_source(e_id);
+        node end = m_graph.get_target(e_id);
+        SASSERT(m_pred[start] == end || m_pred[end] == start);
+        return m_pred[start] == end;
+    }
+
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::get_descendants(node start, svector<node> & descendants) {
+        descendants.reset();
+        node u = start;
+        while (m_depth[m_thread[u]] > m_depth[start]) {
+            descendants.push_back(u);
+            u = m_thread[u];
+        }
+    }
+
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::is_ancestor_of(node ancestor, node child) {
+        for (node n = child; n != -1; n = m_pred[n]) {
+            if (n == ancestor) {
+                return true;
+            }
+        }
+        return false;
+    }
+        
     /**
+        \brief add entering_edge, remove leaving_edge from spanning tree.
+
+        Old tree:                    New tree:
+                       root                      root
+                     /      \                  /      \
+                    x       y                  x       y
+                  /   \    / \               /   \    / \
+                       u     s                   u      s
+                       |    /                          /
+                       v    w                    v    w
+                      / \    \                  / \    \
+                         z   p                     z   p
+                          \                         \ /
+                           q                         q
+        */
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::update(edge_id enter_id, edge_id leave_id, bool & is_swap_enter, bool & is_swap_leave) {    
+        node p = m_graph.get_source(enter_id);
+        node q = m_graph.get_target(enter_id);
+        node u = m_graph.get_source(leave_id);
+        node v = m_graph.get_target(leave_id);
+        
+        if (m_pred[u] == v) {
+            std::swap(u, v);
+            is_swap_leave = true;
+        }
+        else {
+            is_swap_leave = false;
+        }
+        SASSERT(m_pred[v] == u);         
+
+        bool prev_upwards = false;   
+        if (is_ancestor_of(v, p)) {
+            std::swap(p, q);
+            prev_upwards = true;            
+        }
+
+        is_swap_enter = prev_upwards;
+        SASSERT(is_ancestor_of(v, q));
+
+        TRACE("network_flow", { 
+            tout << "update_spanning_tree: (" << p << ", " << q << ") enters, (";
+            tout << u << ", " << v << ") leaves\n";
+        });
+
+        // Update m_pred (for nodes in the stem from q to v)
+        // Note: m_pred[v] == u
+        // Initialize m_upwards[q] = q_upwards
+
+        node old_pred = m_pred[q];        
+        if (q != v) {
+            for (node n = q; n != u; ) {
+                SASSERT(old_pred != u || n == v); // the last processed node is v
+                SASSERT(-1 != m_pred[old_pred]);
+                int next_old_pred = m_pred[old_pred];  
+                swap_order(n, old_pred);
+                std::swap(m_upwards[n], prev_upwards);
+                prev_upwards = !prev_upwards; // flip previous version of upwards.
+                n = old_pred;
+                old_pred = next_old_pred;
+            }     
+        }
+        m_pred[q] = p;
+
+        // m_thread were updated.
+        // update the depth.
+
+        fix_depth(q, get_final(q));
+
+        TRACE("network_flow", {
+            tout << pp_vector("Predecessors", m_pred, true) << pp_vector("Threads", m_thread); 
+            tout << pp_vector("Depths", m_depth) << pp_vector("Upwards", m_upwards);
+            });
+    }
+
+    /**
+        \brief Check invariants of main data-structures.
+
+        Spanning tree of m_graph + root is represented using:
+        
+        svector<edge_state> m_states;      edge_id |-> edge_state
+        svector<bool> m_upwards;           node |-> bool
+        svector<node> m_pred;              node |-> node
+        svector<int>  m_depth;             node |-> int
+        svector<node> m_thread;            node |-> node
+
+        Tree is determined by m_pred:
+        - m_pred[root] == -1
+        - m_pred[n] = m != n     for each node n, acyclic until reaching root.
+        - m_depth[m_pred[n]] + 1 == m_depth[n] for each n != root        
+
+        m_thread is a linked list traversing all nodes.
+        Furthermore, the nodes linked in m_thread follows a 
+        depth-first traversal order.
+
+        m_upwards direction of edge from i to m_pred[i] m_graph
+                
+    */
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::check_well_formed() {
+        node root = m_pred.size()-1;
+
+        // Check that m_thread traverses each node.
+        // This gets checked using union-find as well.
+        svector<bool> found(m_thread.size(), false);
+        found[root] = true;
+        for (node x = m_thread[root]; x != root; x = m_thread[x]) {
+            SASSERT(x != m_thread[x]);
+            found[x] = true;
+        }
+        for (unsigned i = 0; i < found.size(); ++i) {
+            SASSERT(found[i]);
+        }
+
+        // m_pred is acyclic, and points to root.
+        SASSERT(m_pred[root] == -1);
+        SASSERT(m_depth[root] == 0);
+        for (node i = 0; i < root; ++i) {
+            SASSERT(m_depth[m_pred[i]] < m_depth[i]);
+        }
+
+        // m_depth[x] denotes distance from x to the root node
+        for (node x = m_thread[root]; x != root; x = m_thread[x]) {
+            SASSERT(m_depth[x] > 0);
+            SASSERT(m_depth[x] == m_depth[m_pred[x]] + 1);
+        }
+
+        // m_thread forms a spanning tree over [0..root]
+        // Union-find structure
+        svector<int> roots(m_pred.size(), -1);
+                      
+        for (node x = m_thread[root]; x != root; x = m_thread[x]) {
+            node y = m_pred[x];
+            // We are now going to check the edge between x and y
+            SASSERT(find(roots, x) != find(roots, y));
+            merge(roots, x, y);
+        }
+
+        // All nodes belong to the same spanning tree
+        for (unsigned i = 0; i < roots.size(); ++i) {            
+            SASSERT(roots[i] + roots.size() == 0 || roots[i] >= 0);            
+        }        
+
+        return true;
+    }
+
+    static unsigned find(svector<int>& roots, unsigned x) {
+        unsigned old_x = x;
+        while (roots[x] >= 0) {
+            x = roots[x];
+        }
+        SASSERT(roots[x] < 0);
+        if (old_x != x) {
+            roots[old_x] = x;
+        }
+        return x;
+    }
+
+    static void merge(svector<int>& roots, unsigned x, unsigned y) {
+        x = find(roots, x);
+        y = find(roots, y);
+        SASSERT(roots[x] < 0 && roots[y] < 0);
+        if (x == y) {
+            return;
+        }
+        if (roots[x] > roots[y]) {
+            std::swap(x, y);
+        }
+        SASSERT(roots[x] <= roots[y]);
+        roots[x] += roots[y];
+        roots[y] = x;
+    }
+
+        /**
         swap v and q in tree.
         - fixup m_thread
         - fixup m_pred
@@ -119,251 +392,6 @@ namespace smt {
         SASSERT(children.empty());
         return true;
     }
-
-    template<typename Ext>
-    bool thread_spanning_tree<Ext>::is_ancestor_of(node ancestor, node child) {
-        for (node n = child; n != -1; n = m_pred[n]) {
-            if (n == ancestor) {
-                return true;
-            }
-        }
-        return false;
-    }
-
-    static unsigned find(svector<int>& roots, unsigned x) {
-        unsigned old_x = x;
-        while (roots[x] >= 0) {
-            x = roots[x];
-        }
-        SASSERT(roots[x] < 0);
-        if (old_x != x) {
-            roots[old_x] = x;
-        }
-        return x;
-    }
-
-    static void merge(svector<int>& roots, unsigned x, unsigned y) {
-        x = find(roots, x);
-        y = find(roots, y);
-        SASSERT(roots[x] < 0 && roots[y] < 0);
-        if (x == y) {
-            return;
-        }
-        if (roots[x] > roots[y]) {
-            std::swap(x, y);
-        }
-        SASSERT(roots[x] <= roots[y]);
-        roots[x] += roots[y];
-        roots[y] = x;
-    }
-
-    template<typename Ext>
-    void thread_spanning_tree<Ext>::initialize(svector<bool> const & upwards, int num_nodes) {
-        m_pred.resize(num_nodes + 1);
-        m_depth.resize(num_nodes + 1);
-        m_thread.resize(num_nodes + 1);
-        m_upwards.resize(num_nodes + 1);
-
-        node root = num_nodes;
-        m_pred[root] = -1;
-        m_depth[root] = 0;
-        m_thread[root] = 0;
-            
-        // Create artificial edges from/to root node to/from other nodes and initialize the spanning tree
-        for (int i = 0; i < num_nodes; ++i) {
-            m_pred[i] = root;
-            m_depth[i] = 1;
-            m_thread[i] = i + 1;
-            m_upwards[i] = upwards[i];
-        }
-
-        TRACE("network_flow", {
-            tout << pp_vector("Predecessors", m_pred, true) << pp_vector("Threads", m_thread); 
-            tout << pp_vector("Depths", m_depth) << pp_vector("Upwards", m_upwards);
-        });
-    }
-
-    template<typename Ext>
-    typename thread_spanning_tree<Ext>::node thread_spanning_tree<Ext>::get_common_ancestor(node u, node v) {
-        while (u != v) {
-            if (m_depth[u] > m_depth[v])
-                u = m_pred[u];
-            else 
-                v = m_pred[v];                    
-        }
-        return u;
-    }
-
-    template<typename Ext>
-    void thread_spanning_tree<Ext>::get_descendants(node start, svector<node>& descendants) {
-        descendants.reset();
-        node u = start;
-        while (m_depth[m_thread[u]] > m_depth[start]) {
-            descendants.push_back(u);
-            u = m_thread[u];
-        }
-    }
-
-    template<typename Ext>
-    void thread_spanning_tree<Ext>::get_ancestors(node start, svector<node>& ancestors) {
-        ancestors.reset();
-        while (m_pred[start] != -1) {                        
-            ancestors.push_back(start);
-            start = m_pred[start];
-        }
-    }
-        
-    /**
-        \brief add entering_edge, remove leaving_edge from spanning tree.
-
-        Old tree:                    New tree:
-                       root                      root
-                     /      \                  /      \
-                    x       y                  x       y
-                  /   \    / \               /   \    / \
-                       u     s                   u      s
-                       |    /                          /
-                       v    w                    v    w
-                      / \    \                  / \    \
-                         z   p                     z   p
-                          \                         \ /
-                           q                         q
-        */
-    template<typename Ext>
-    void thread_spanning_tree<Ext>::update(node p, node q, node u, node v) {
-        bool q_upwards = false;
-
-        // v is parent of u so T_u does not contain root node
-        if (m_pred[u] == v) {
-            std::swap(u, v);
-        }  
-        SASSERT(m_pred[v] == u);
-
-        if (is_ancestor_of(v, p)) {
-            std::swap(p, q);
-            q_upwards = true;
-        }
-        SASSERT(is_ancestor_of(v, q));
-
-        TRACE("network_flow", { 
-            tout << "update_spanning_tree: (" << p << ", " << q << ") enters, (";
-            tout << u << ", " << v << ") leaves\n";
-        });
-
-        // Update m_pred (for nodes in the stem from q to v)
-        // Note: m_pred[v] == u
-        // Initialize m_upwards[q] = q_upwards
-
-        bool prev_upwards = q_upwards;   
-        node old_pred = m_pred[q];        
-        if (q != v) {
-            for (node n = q; n != u; ) {
-                SASSERT(old_pred != u || n == v); // the last processed node is v
-                TRACE("network_flow", {
-                    tout << pp_vector("Predecessors", m_pred, true);
-                });
-                SASSERT(-1 != m_pred[old_pred]);
-                int next_old_pred = m_pred[old_pred];  
-                swap_order(n, old_pred);
-                std::swap(m_upwards[n], prev_upwards);
-                prev_upwards = !prev_upwards; // flip previous version of upwards.
-                n = old_pred;
-                old_pred = next_old_pred;
-            }     
-        }
-        m_pred[q] = p;
-
-        // m_thread were updated.
-        // update the depth.
-
-        fix_depth(q, get_final(q));
-
-        TRACE("network_flow", {
-            tout << pp_vector("Predecessors", m_pred, true) << pp_vector("Threads", m_thread); 
-            tout << pp_vector("Depths", m_depth) << pp_vector("Upwards", m_upwards);
-            });
-    }
-
-    /**
-        \brief Check invariants of main data-structures.
-
-        Spanning tree of m_graph + root is represented using:
-        
-        svector<edge_state> m_states;      edge_id |-> edge_state
-        svector<bool> m_upwards;           node |-> bool
-        svector<node> m_pred;              node |-> node
-        svector<int>  m_depth;             node |-> int
-        svector<node> m_thread;            node |-> node
-
-        Tree is determined by m_pred:
-        - m_pred[root] == -1
-        - m_pred[n] = m != n     for each node n, acyclic until reaching root.
-        - m_depth[m_pred[n]] + 1 == m_depth[n] for each n != root        
-
-        m_thread is a linked list traversing all nodes.
-        Furthermore, the nodes linked in m_thread follows a 
-        depth-first traversal order.
-
-        m_upwards direction of edge from i to m_pred[i] m_graph
-                
-    */
-    template<typename Ext>
-    bool thread_spanning_tree<Ext>::check_well_formed() {
-        node root = m_pred.size()-1;
-
-        // Check that m_thread traverses each node.
-        // This gets checked using union-find as well.
-        svector<bool> found(m_thread.size(), false);
-        found[root] = true;
-        for (node x = m_thread[root]; x != root; x = m_thread[x]) {
-            SASSERT(x != m_thread[x]);
-            found[x] = true;
-        }
-        for (unsigned i = 0; i < found.size(); ++i) {
-            SASSERT(found[i]);
-        }
-
-        // m_pred is acyclic, and points to root.
-        SASSERT(m_pred[root] == -1);
-        SASSERT(m_depth[root] == 0);
-        for (node i = 0; i < root; ++i) {
-            SASSERT(m_depth[m_pred[i]] < m_depth[i]);
-        }
-
-        // m_depth[x] denotes distance from x to the root node
-        for (node x = m_thread[root]; x != root; x = m_thread[x]) {
-            SASSERT(m_depth[x] > 0);
-            SASSERT(m_depth[x] == m_depth[m_pred[x]] + 1);
-        }
-
-        // m_thread forms a spanning tree over [0..root]
-        // Union-find structure
-        svector<int> roots(m_pred.size(), -1);
-                      
-        for (node x = m_thread[root]; x != root; x = m_thread[x]) {
-            node y = m_pred[x];
-            // We are now going to check the edge between x and y
-            SASSERT(find(roots, x) != find(roots, y));
-            merge(roots, x, y);
-        }
-
-        // All nodes belong to the same spanning tree
-        for (unsigned i = 0; i < roots.size(); ++i) {            
-            SASSERT(roots[i] + roots.size() == 0 || roots[i] >= 0);            
-        }        
-
-        return true;
-    }
-
-    template<typename Ext>
-    bool thread_spanning_tree<Ext>::get_arc_direction(node start) const {
-        return m_upwards[start];
-    }
-
-    template<typename Ext>
-    typename thread_spanning_tree<Ext>::node thread_spanning_tree<Ext>::get_parent(node start) {
-        return m_pred[start];
-    }
 }
 
 #endif