boost/graph/minimum_degree_ordering.hpp
//-*-c++-*- //======================================================================= // Copyright 1997-2001 University of Notre Dame. // Authors: Lie-Quan Lee, Jeremy Siek // // Distributed under the Boost Software License, Version 1.0. (See // accompanying file LICENSE_1_0.txt or copy at // http://www.boost.org/LICENSE_1_0.txt) //======================================================================= // #ifndef BOOST_GRAPH_MINIMUM_DEGREE_ORDERING_HPP #define BOOST_GRAPH_MINIMUM_DEGREE_ORDERING_HPP #include <vector> #include <boost/assert.hpp> #include <boost/config.hpp> #include <boost/pending/bucket_sorter.hpp> #include <boost/detail/numeric_traits.hpp> // for integer_traits #include <boost/graph/graph_traits.hpp> #include <boost/property_map/property_map.hpp> namespace boost { namespace detail { // // Given a set of n integers (where the integer values range from // zero to n-1), we want to keep track of a collection of stacks // of integers. It so happens that an integer will appear in at // most one stack at a time, so the stacks form disjoint sets. // Because of these restrictions, we can use one big array to // store all the stacks, intertwined with one another. // No allocation/deallocation happens in the push()/pop() methods // so this is faster than using std::stack's. // template < class SignedInteger > class Stacks { typedef SignedInteger value_type; typedef typename std::vector< value_type >::size_type size_type; public: Stacks(size_type n) : data(n) {} //: stack class stack { typedef typename std::vector< value_type >::iterator Iterator; public: stack(Iterator _data, const value_type& head) : data(_data), current(head) { } // did not use default argument here to avoid internal compiler // error in g++. stack(Iterator _data) : data(_data), current(-(std::numeric_limits< value_type >::max)()) { } void pop() { BOOST_ASSERT(!empty()); current = data[current]; } void push(value_type v) { data[v] = current; current = v; } bool empty() { return current == -(std::numeric_limits< value_type >::max)(); } value_type& top() { return current; } private: Iterator data; value_type current; }; // To return a stack object stack make_stack() { return stack(data.begin()); } protected: std::vector< value_type > data; }; // marker class, a generalization of coloring. // // This class is to provide a generalization of coloring which has // complexity of amortized constant time to set all vertices' color // back to be untagged. It implemented by increasing a tag. // // The colors are: // not tagged // tagged // multiple_tagged // done // template < class SignedInteger, class Vertex, class VertexIndexMap > class Marker { typedef SignedInteger value_type; typedef typename std::vector< value_type >::size_type size_type; static value_type done() { return (std::numeric_limits< value_type >::max)() / 2; } public: Marker(size_type _num, VertexIndexMap index_map) : tag(1 - (std::numeric_limits< value_type >::max)()) , data(_num, -(std::numeric_limits< value_type >::max)()) , id(index_map) { } void mark_done(Vertex node) { data[get(id, node)] = done(); } bool is_done(Vertex node) { return data[get(id, node)] == done(); } void mark_tagged(Vertex node) { data[get(id, node)] = tag; } void mark_multiple_tagged(Vertex node) { data[get(id, node)] = multiple_tag; } bool is_tagged(Vertex node) const { return data[get(id, node)] >= tag; } bool is_not_tagged(Vertex node) const { return data[get(id, node)] < tag; } bool is_multiple_tagged(Vertex node) const { return data[get(id, node)] >= multiple_tag; } void increment_tag() { const size_type num = data.size(); ++tag; if (tag >= done()) { tag = 1 - (std::numeric_limits< value_type >::max)(); for (size_type i = 0; i < num; ++i) if (data[i] < done()) data[i] = -(std::numeric_limits< value_type >::max)(); } } void set_multiple_tag(value_type mdeg0) { const size_type num = data.size(); multiple_tag = tag + mdeg0; if (multiple_tag >= done()) { tag = 1 - (std::numeric_limits< value_type >::max)(); for (size_type i = 0; i < num; i++) if (data[i] < done()) data[i] = -(std::numeric_limits< value_type >::max)(); multiple_tag = tag + mdeg0; } } void set_tag_as_multiple_tag() { tag = multiple_tag; } protected: value_type tag; value_type multiple_tag; std::vector< value_type > data; VertexIndexMap id; }; template < class Iterator, class SignedInteger, class Vertex, class VertexIndexMap, int offset = 1 > class Numbering { typedef SignedInteger number_type; number_type num; // start from 1 instead of zero Iterator data; number_type max_num; VertexIndexMap id; public: Numbering(Iterator _data, number_type _max_num, VertexIndexMap id) : num(1), data(_data), max_num(_max_num), id(id) { } void operator()(Vertex node) { data[get(id, node)] = -num; } bool all_done(number_type i = 0) const { return num + i > max_num; } void increment(number_type i = 1) { num += i; } bool is_numbered(Vertex node) const { return data[get(id, node)] < 0; } void indistinguishable(Vertex i, Vertex j) { data[get(id, i)] = -(get(id, j) + offset); } }; template < class SignedInteger, class Vertex, class VertexIndexMap > class degreelists_marker { public: typedef SignedInteger value_type; typedef typename std::vector< value_type >::size_type size_type; degreelists_marker(size_type n, VertexIndexMap id) : marks(n, 0), id(id) { } void mark_need_update(Vertex i) { marks[get(id, i)] = 1; } bool need_update(Vertex i) { return marks[get(id, i)] == 1; } bool outmatched_or_done(Vertex i) { return marks[get(id, i)] == -1; } void mark(Vertex i) { marks[get(id, i)] = -1; } void unmark(Vertex i) { marks[get(id, i)] = 0; } private: std::vector< value_type > marks; VertexIndexMap id; }; // Helper function object for edge removal template < class Graph, class MarkerP, class NumberD, class Stack, class VertexIndexMap > class predicateRemoveEdge1 { typedef typename graph_traits< Graph >::vertex_descriptor vertex_t; typedef typename graph_traits< Graph >::edge_descriptor edge_t; public: predicateRemoveEdge1(Graph& _g, MarkerP& _marker, NumberD _numbering, Stack& n_e, VertexIndexMap id) : g(&_g) , marker(&_marker) , numbering(_numbering) , neighbor_elements(&n_e) , id(id) { } bool operator()(edge_t e) { vertex_t dist = target(e, *g); if (marker->is_tagged(dist)) return true; marker->mark_tagged(dist); if (numbering.is_numbered(dist)) { neighbor_elements->push(get(id, dist)); return true; } return false; } private: Graph* g; MarkerP* marker; NumberD numbering; Stack* neighbor_elements; VertexIndexMap id; }; // Helper function object for edge removal template < class Graph, class MarkerP > class predicate_remove_tagged_edges { typedef typename graph_traits< Graph >::vertex_descriptor vertex_t; typedef typename graph_traits< Graph >::edge_descriptor edge_t; public: predicate_remove_tagged_edges(Graph& _g, MarkerP& _marker) : g(&_g), marker(&_marker) { } bool operator()(edge_t e) { vertex_t dist = target(e, *g); if (marker->is_tagged(dist)) return true; return false; } private: Graph* g; MarkerP* marker; }; template < class Graph, class DegreeMap, class InversePermutationMap, class PermutationMap, class SuperNodeMap, class VertexIndexMap > class mmd_impl { // Typedefs typedef graph_traits< Graph > Traits; typedef typename Traits::vertices_size_type size_type; typedef typename detail::integer_traits< size_type >::difference_type diff_t; typedef typename Traits::vertex_descriptor vertex_t; typedef typename Traits::adjacency_iterator adj_iter; typedef iterator_property_map< vertex_t*, identity_property_map, vertex_t, vertex_t& > IndexVertexMap; typedef detail::Stacks< diff_t > Workspace; typedef bucket_sorter< size_type, vertex_t, DegreeMap, VertexIndexMap > DegreeLists; typedef Numbering< InversePermutationMap, diff_t, vertex_t, VertexIndexMap > NumberingD; typedef degreelists_marker< diff_t, vertex_t, VertexIndexMap > DegreeListsMarker; typedef Marker< diff_t, vertex_t, VertexIndexMap > MarkerP; // Data Members bool has_no_edges; // input parameters Graph& G; int delta; DegreeMap degree; InversePermutationMap inverse_perm; PermutationMap perm; SuperNodeMap supernode_size; VertexIndexMap vertex_index_map; // internal data-structures std::vector< vertex_t > index_vertex_vec; size_type n; IndexVertexMap index_vertex_map; DegreeLists degreelists; NumberingD numbering; DegreeListsMarker degree_lists_marker; MarkerP marker; Workspace work_space; public: mmd_impl(Graph& g, size_type n_, int delta, DegreeMap degree, InversePermutationMap inverse_perm, PermutationMap perm, SuperNodeMap supernode_size, VertexIndexMap id) : has_no_edges(true) , G(g) , delta(delta) , degree(degree) , inverse_perm(inverse_perm) , perm(perm) , supernode_size(supernode_size) , vertex_index_map(id) , index_vertex_vec(n_) , n(n_) , degreelists(n_ + 1, n_, degree, id) , numbering(inverse_perm, n_, vertex_index_map) , degree_lists_marker(n_, vertex_index_map) , marker(n_, vertex_index_map) , work_space(n_) { typename graph_traits< Graph >::vertex_iterator v, vend; size_type vid = 0; for (boost::tie(v, vend) = vertices(G); v != vend; ++v, ++vid) index_vertex_vec[vid] = *v; index_vertex_map = IndexVertexMap(&index_vertex_vec[0]); // Initialize degreelists. Degreelists organizes the nodes // according to their degree. for (boost::tie(v, vend) = vertices(G); v != vend; ++v) { typename Traits::degree_size_type d = out_degree(*v, G); put(degree, *v, d); if (0 < d) has_no_edges = false; degreelists.push(*v); } } void do_mmd() { // Eliminate the isolated nodes -- these are simply the nodes // with no neighbors, which are accessible as a list (really, a // stack) at location 0. Since these don't affect any other // nodes, we can eliminate them without doing degree updates. typename DegreeLists::stack list_isolated = degreelists[0]; while (!list_isolated.empty()) { vertex_t node = list_isolated.top(); marker.mark_done(node); numbering(node); numbering.increment(); list_isolated.pop(); } if (has_no_edges) { return; } size_type min_degree = 1; typename DegreeLists::stack list_min_degree = degreelists[min_degree]; while (list_min_degree.empty()) { ++min_degree; list_min_degree = degreelists[min_degree]; } // check if the whole eliminating process is done while (!numbering.all_done()) { size_type min_degree_limit = min_degree + delta; // WARNING typename Workspace::stack llist = work_space.make_stack(); // multiple elimination while (delta >= 0) { // Find the next non-empty degree for (list_min_degree = degreelists[min_degree]; list_min_degree.empty() && min_degree <= min_degree_limit; ++min_degree, list_min_degree = degreelists[min_degree]) ; if (min_degree > min_degree_limit) break; const vertex_t node = list_min_degree.top(); const size_type node_id = get(vertex_index_map, node); list_min_degree.pop(); numbering(node); // check if node is the last one if (numbering.all_done(supernode_size[node])) { numbering.increment(supernode_size[node]); break; } marker.increment_tag(); marker.mark_tagged(node); this->eliminate(node); numbering.increment(supernode_size[node]); llist.push(node_id); } // multiple elimination if (numbering.all_done()) break; this->update(llist, min_degree); } } // do_mmd() void eliminate(vertex_t node) { typename Workspace::stack element_neighbor = work_space.make_stack(); // Create two function objects for edge removal typedef typename Workspace::stack WorkStack; predicateRemoveEdge1< Graph, MarkerP, NumberingD, WorkStack, VertexIndexMap > p(G, marker, numbering, element_neighbor, vertex_index_map); predicate_remove_tagged_edges< Graph, MarkerP > p2(G, marker); // Reconstruct the adjacent node list, push element neighbor in a // List. remove_out_edge_if(node, p, G); // during removal element neighbors are collected. while (!element_neighbor.empty()) { // element absorb size_type e_id = element_neighbor.top(); vertex_t element = get(index_vertex_map, e_id); adj_iter i, i_end; for (boost::tie(i, i_end) = adjacent_vertices(element, G); i != i_end; ++i) { vertex_t i_node = *i; if (!marker.is_tagged(i_node) && !numbering.is_numbered(i_node)) { marker.mark_tagged(i_node); add_edge(node, i_node, G); } } element_neighbor.pop(); } adj_iter v, ve; for (boost::tie(v, ve) = adjacent_vertices(node, G); v != ve; ++v) { vertex_t v_node = *v; if (!degree_lists_marker.need_update(v_node) && !degree_lists_marker.outmatched_or_done(v_node)) { degreelists.remove(v_node); } // update out edges of v_node remove_out_edge_if(v_node, p2, G); if (out_degree(v_node, G) == 0) { // indistinguishable nodes supernode_size[node] += supernode_size[v_node]; supernode_size[v_node] = 0; numbering.indistinguishable(v_node, node); marker.mark_done(v_node); degree_lists_marker.mark(v_node); } else { // not indistinguishable nodes add_edge(v_node, node, G); degree_lists_marker.mark_need_update(v_node); } } } // eliminate() template < class Stack > void update(Stack llist, size_type& min_degree) { size_type min_degree0 = min_degree + delta + 1; while (!llist.empty()) { size_type deg, deg0 = 0; marker.set_multiple_tag(min_degree0); typename Workspace::stack q2list = work_space.make_stack(); typename Workspace::stack qxlist = work_space.make_stack(); vertex_t current = get(index_vertex_map, llist.top()); adj_iter i, ie; for (boost::tie(i, ie) = adjacent_vertices(current, G); i != ie; ++i) { vertex_t i_node = *i; const size_type i_id = get(vertex_index_map, i_node); if (supernode_size[i_node] != 0) { deg0 += supernode_size[i_node]; marker.mark_multiple_tagged(i_node); if (degree_lists_marker.need_update(i_node)) { if (out_degree(i_node, G) == 2) q2list.push(i_id); else qxlist.push(i_id); } } } while (!q2list.empty()) { const size_type u_id = q2list.top(); vertex_t u_node = get(index_vertex_map, u_id); // if u_id is outmatched by others, no need to update degree if (degree_lists_marker.outmatched_or_done(u_node)) { q2list.pop(); continue; } marker.increment_tag(); deg = deg0; adj_iter nu = adjacent_vertices(u_node, G).first; vertex_t neighbor = *nu; if (neighbor == u_node) { ++nu; neighbor = *nu; } if (numbering.is_numbered(neighbor)) { adj_iter i, ie; for (boost::tie(i, ie) = adjacent_vertices(neighbor, G); i != ie; ++i) { const vertex_t i_node = *i; if (i_node == u_node || supernode_size[i_node] == 0) continue; if (marker.is_tagged(i_node)) { if (degree_lists_marker.need_update(i_node)) { if (out_degree(i_node, G) == 2) { // is indistinguishable supernode_size[u_node] += supernode_size[i_node]; supernode_size[i_node] = 0; numbering.indistinguishable( i_node, u_node); marker.mark_done(i_node); degree_lists_marker.mark(i_node); } else // is outmatched degree_lists_marker.mark(i_node); } } else { marker.mark_tagged(i_node); deg += supernode_size[i_node]; } } } else deg += supernode_size[neighbor]; deg -= supernode_size[u_node]; degree[u_node] = deg; // update degree degreelists[deg].push(u_node); // u_id has been pushed back into degreelists degree_lists_marker.unmark(u_node); if (min_degree > deg) min_degree = deg; q2list.pop(); } // while (!q2list.empty()) while (!qxlist.empty()) { const size_type u_id = qxlist.top(); const vertex_t u_node = get(index_vertex_map, u_id); // if u_id is outmatched by others, no need to update degree if (degree_lists_marker.outmatched_or_done(u_node)) { qxlist.pop(); continue; } marker.increment_tag(); deg = deg0; adj_iter i, ie; for (boost::tie(i, ie) = adjacent_vertices(u_node, G); i != ie; ++i) { vertex_t i_node = *i; if (marker.is_tagged(i_node)) continue; marker.mark_tagged(i_node); if (numbering.is_numbered(i_node)) { adj_iter j, je; for (boost::tie(j, je) = adjacent_vertices(i_node, G); j != je; ++j) { const vertex_t j_node = *j; if (marker.is_not_tagged(j_node)) { marker.mark_tagged(j_node); deg += supernode_size[j_node]; } } } else deg += supernode_size[i_node]; } // for adjacent vertices of u_node deg -= supernode_size[u_node]; degree[u_node] = deg; degreelists[deg].push(u_node); // u_id has been pushed back into degreelists degree_lists_marker.unmark(u_node); if (min_degree > deg) min_degree = deg; qxlist.pop(); } // while (!qxlist.empty()) { marker.set_tag_as_multiple_tag(); llist.pop(); } // while (! llist.empty()) } // update() void build_permutation(InversePermutationMap next, PermutationMap prev) { // collect the permutation info size_type i; for (i = 0; i < n; ++i) { diff_t size = supernode_size[get(index_vertex_map, i)]; if (size <= 0) { prev[i] = next[i]; supernode_size[get(index_vertex_map, i)] = next[i] + 1; // record the supernode info } else prev[i] = -next[i]; } for (i = 1; i < n + 1; ++i) { if (prev[i - 1] > 0) continue; diff_t parent = i; while (prev[parent - 1] < 0) { parent = -prev[parent - 1]; } diff_t root = parent; diff_t num = prev[root - 1] + 1; next[i - 1] = -num; prev[root - 1] = num; parent = i; diff_t next_node = -prev[parent - 1]; while (next_node > 0) { prev[parent - 1] = -root; parent = next_node; next_node = -prev[parent - 1]; } } for (i = 0; i < n; i++) { diff_t num = -next[i] - 1; next[i] = num; prev[num] = i; } } // build_permutation() }; } // namespace detail // MMD algorithm // // The implementation presently includes the enhancements for mass // elimination, incomplete degree update, multiple elimination, and // external degree. // // Important Note: This implementation requires the BGL graph to be // directed. Therefore, nonzero entry (i, j) in a symmetrical matrix // A coresponds to two directed edges (i->j and j->i). // // see Alan George and Joseph W. H. Liu, The Evolution of the Minimum // Degree Ordering Algorithm, SIAM Review, 31, 1989, Page 1-19 template < class Graph, class DegreeMap, class InversePermutationMap, class PermutationMap, class SuperNodeMap, class VertexIndexMap > void minimum_degree_ordering(Graph& G, DegreeMap degree, InversePermutationMap inverse_perm, PermutationMap perm, SuperNodeMap supernode_size, int delta, VertexIndexMap vertex_index_map) { detail::mmd_impl< Graph, DegreeMap, InversePermutationMap, PermutationMap, SuperNodeMap, VertexIndexMap > impl(G, num_vertices(G), delta, degree, inverse_perm, perm, supernode_size, vertex_index_map); impl.do_mmd(); impl.build_permutation(inverse_perm, perm); } } // namespace boost #endif // BOOST_GRAPH_MINIMUM_DEGREE_ORDERING_HPP