Chris@16: // Copyright (C) 2004-2006 The Trustees of Indiana University.
Chris@16: 
Chris@16: // Use, modification and distribution is subject to the Boost Software
Chris@16: // License, Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
Chris@16: // http://www.boost.org/LICENSE_1_0.txt)
Chris@16: 
Chris@16: //  Authors: Nick Edmonds
Chris@16: //           Douglas Gregor
Chris@16: //           Andrew Lumsdaine
Chris@16: #ifndef BOOST_GRAPH_PARALLEL_CC_HPP
Chris@16: #define BOOST_GRAPH_PARALLEL_CC_HPP
Chris@16: 
Chris@16: #ifndef BOOST_GRAPH_USE_MPI
Chris@16: #error "Parallel BGL files should not be included unless <boost/graph/use_mpi.hpp> has been included"
Chris@16: #endif
Chris@16: 
Chris@16: #include <boost/detail/is_sorted.hpp>
Chris@16: #include <boost/assert.hpp>
Chris@16: #include <boost/property_map/property_map.hpp>
Chris@16: #include <boost/property_map/parallel/caching_property_map.hpp>
Chris@16: #include <boost/graph/parallel/algorithm.hpp>
Chris@16: #include <boost/pending/indirect_cmp.hpp>
Chris@16: #include <boost/graph/graph_traits.hpp>
Chris@16: #include <boost/graph/overloading.hpp>
Chris@16: #include <boost/graph/distributed/concepts.hpp>
Chris@16: #include <boost/graph/parallel/properties.hpp>
Chris@16: #include <boost/graph/distributed/local_subgraph.hpp>
Chris@16: #include <boost/graph/connected_components.hpp>
Chris@16: #include <boost/graph/named_function_params.hpp>
Chris@16: #include <boost/graph/parallel/process_group.hpp>
Chris@16: #include <boost/optional.hpp>
Chris@16: #include <functional>
Chris@16: #include <algorithm>
Chris@16: #include <vector>
Chris@16: #include <list>
Chris@16: #include <boost/graph/parallel/container_traits.hpp>
Chris@16: #include <boost/graph/iteration_macros.hpp>
Chris@16: 
Chris@16: #define PBGL_IN_PLACE_MERGE /* In place merge instead of sorting */
Chris@16: //#define PBGL_SORT_ASSERT    /* Assert sorted for in place merge */
Chris@16: 
Chris@16: /* Explicit sychronization in pointer doubling step? */
Chris@16: #define PBGL_EXPLICIT_SYNCH
Chris@16: //#define PBGL_CONSTRUCT_METAGRAPH
Chris@16: #ifdef PBGL_CONSTRUCT_METAGRAPH
Chris@16: #  define MAX_VERTICES_IN_METAGRAPH 10000
Chris@16: #endif
Chris@16: 
Chris@16: namespace boost { namespace graph { namespace distributed {
Chris@16:   namespace cc_detail {
Chris@16:     enum connected_components_message { 
Chris@16:       edges_msg, req_parents_msg, parents_msg, root_adj_msg
Chris@16:     };
Chris@16: 
Chris@16:     template <typename Vertex>
Chris@16:     struct metaVertex {
Chris@16:       metaVertex() {}
Chris@16:       metaVertex(const Vertex& v) : name(v) {}
Chris@16: 
Chris@16:       template<typename Archiver>
Chris@16:       void serialize(Archiver& ar, const unsigned int /*version*/)
Chris@16:       {
Chris@16:         ar & name;
Chris@16:       }
Chris@16: 
Chris@16:       Vertex name;
Chris@16:     };
Chris@16: 
Chris@16: #ifdef PBGL_CONSTRUCT_METAGRAPH
Chris@16:     // Build meta-graph on result of local connected components
Chris@16:     template <typename Graph, typename ParentMap, typename RootIterator,
Chris@16:               typename AdjacencyMap>
Chris@16:     void
Chris@16:     build_local_metagraph(const Graph& g, ParentMap p, RootIterator r,
Chris@16:                           RootIterator r_end, AdjacencyMap& adj)
Chris@16:     {
Chris@16:       // TODO: Static assert that AdjacencyMap::value_type is std::vector<vertex_descriptor>
Chris@16: 
Chris@16:       typedef typename boost::graph::parallel::process_group_type<Graph>::type
Chris@16:         process_group_type;
Chris@16:       typedef typename process_group_type::process_id_type process_id_type;
Chris@16: 
Chris@16:       typedef typename graph_traits<Graph>::vertex_descriptor vertex_descriptor;
Chris@16: 
Chris@16:       BOOST_STATIC_ASSERT((is_same<typename AdjacencyMap::mapped_type,
Chris@16:                                      std::vector<vertex_descriptor> >::value));
Chris@16: 
Chris@16:       using boost::graph::parallel::process_group;
Chris@16: 
Chris@16:       process_group_type pg = process_group(g);
Chris@16:       process_id_type id = process_id(pg);
Chris@16:       
Chris@16:       if (id != 0) {
Chris@16: 
Chris@16:         // Send component roots and their associated edges to P0
Chris@16:         for ( ; r != r_end; ++r ) {
Chris@16:           std::vector<vertex_descriptor> adjs(1, *r); // Root
Chris@16:           adjs.reserve(adjs.size() + adj[*r].size());
Chris@16:           for (typename std::vector<vertex_descriptor>::iterator iter = adj[*r].begin();
Chris@16:                iter != adj[*r].end(); ++iter)
Chris@16:             adjs.push_back(get(p, *iter)); // Adjacencies
Chris@16: 
Chris@16:           send(pg, 0, root_adj_msg, adjs); 
Chris@16:         }
Chris@16:       }
Chris@16:       
Chris@16:       synchronize(pg);
Chris@16: 
Chris@16:       if (id == 0) {
Chris@16:         typedef metaVertex<vertex_descriptor> VertexProperties;
Chris@16: 
Chris@16:         typedef boost::adjacency_list<vecS, vecS, undirectedS, 
Chris@16:           VertexProperties> metaGraph;
Chris@16:         typedef typename graph_traits<metaGraph>::vertex_descriptor 
Chris@16:           meta_vertex_descriptor;
Chris@16: 
Chris@16:         std::map<vertex_descriptor, meta_vertex_descriptor> vertex_map;
Chris@16:         std::vector<std::pair<vertex_descriptor, vertex_descriptor> > edges;
Chris@16: 
Chris@16:         // Receive remote roots and edges
Chris@16:         while (optional<std::pair<process_id_type, int> > m = probe(pg)) {
Chris@16:           BOOST_ASSERT(m->second == root_adj_msg);
Chris@16: 
Chris@16:           std::vector<vertex_descriptor> adjs;
Chris@16:           receive(pg, m->first, m->second, adjs);
Chris@16: 
Chris@16:           vertex_map[adjs[0]] = graph_traits<metaGraph>::null_vertex();
Chris@16:           for (typename std::vector<vertex_descriptor>::iterator iter 
Chris@16:                  = ++adjs.begin(); iter != adjs.end(); ++iter)
Chris@16:             edges.push_back(std::make_pair(adjs[0], *iter));
Chris@16:         }
Chris@16: 
Chris@16:         // Add local roots and edges
Chris@16:         for ( ; r != r_end; ++r ) {
Chris@16:           vertex_map[*r] = graph_traits<metaGraph>::null_vertex();
Chris@16:           edges.reserve(edges.size() + adj[*r].size());
Chris@16:           for (typename std::vector<vertex_descriptor>::iterator iter = adj[*r].begin();
Chris@16:                iter != adj[*r].end(); ++iter)
Chris@16:             edges.push_back(std::make_pair(*r, get(p, *iter)));
Chris@16:         } 
Chris@16: 
Chris@16:         // Build local meta-graph
Chris@16:         metaGraph mg;
Chris@16: 
Chris@16:         // Add vertices with property to map back to distributed graph vertex
Chris@16:         for (typename std::map<vertex_descriptor, meta_vertex_descriptor>::iterator
Chris@16:                iter = vertex_map.begin(); iter != vertex_map.end(); ++iter)
Chris@16:           vertex_map[iter->first] 
Chris@16:             = add_vertex(metaVertex<vertex_descriptor>(iter->first), mg);
Chris@16: 
Chris@16:         // Build meta-vertex map
Chris@16:         typename property_map<metaGraph, vertex_descriptor VertexProperties::*>::type 
Chris@16:           metaVertexMap = get(&VertexProperties::name, mg);
Chris@16: 
Chris@16:         typename std::vector<std::pair<vertex_descriptor, vertex_descriptor> >
Chris@16:           ::iterator edge_iter = edges.begin();
Chris@16:         for ( ; edge_iter != edges.end(); ++edge_iter)
Chris@16:           add_edge(vertex_map[edge_iter->first], vertex_map[edge_iter->second], mg);
Chris@16:         
Chris@16:         edges.clear();
Chris@16:   
Chris@16:         // Call connected_components on it
Chris@16:         typedef typename property_map<metaGraph, vertex_index_t>::type 
Chris@16:           meta_index_map_type;
Chris@16:         meta_index_map_type meta_index = get(vertex_index, mg);
Chris@16: 
Chris@16:         std::vector<std::size_t> mg_component_vec(num_vertices(mg));
Chris@16:         typedef iterator_property_map<std::vector<std::size_t>::iterator,
Chris@16:                                       meta_index_map_type>
Chris@16:         meta_components_map_type;
Chris@16:         meta_components_map_type mg_component(mg_component_vec.begin(),
Chris@16:                                               meta_index);
Chris@16:         std::size_t num_comp = connected_components(mg, mg_component);
Chris@16: 
Chris@16:         // Update Parent pointers
Chris@16:         std::vector<meta_vertex_descriptor> roots(num_comp, graph_traits<metaGraph>::null_vertex());
Chris@16: 
Chris@16:         BGL_FORALL_VERTICES_T(v, mg, metaGraph) {
Chris@16:           size_t component = get(mg_component, v);
Chris@16:           if (roots[component] == graph_traits<metaGraph>::null_vertex() ||
Chris@16:               get(meta_index, v) < get(meta_index, roots[component])) 
Chris@16:             roots[component] = v;
Chris@16:         }
Chris@16: 
Chris@16:         // Set all the local parent pointers
Chris@16:         BGL_FORALL_VERTICES_T(v, mg, metaGraph) {
Chris@16:           // Problem in value being put (3rd parameter)
Chris@16:           put(p, get(metaVertexMap, v), get(metaVertexMap, roots[get(mg_component, v)]));
Chris@16:         }
Chris@16:       }
Chris@16: 
Chris@16:       synchronize(p);
Chris@16:     }
Chris@16: #endif
Chris@16: 
Chris@16:     /* Function object used to remove internal vertices and vertices >
Chris@16:        the current vertex from the adjacent vertex lists at each
Chris@16:        root */
Chris@16:     template <typename Vertex, typename ParentMap>
Chris@16:     class cull_adjacency_list
Chris@16:     {
Chris@16:     public:
Chris@16:       cull_adjacency_list(const Vertex v, const ParentMap p) : v(v), p(p) {}
Chris@16:       bool operator() (const Vertex x) { return (get(p, x) == v || x == v); }
Chris@16: 
Chris@16:     private:
Chris@16:       const Vertex    v;
Chris@16:       const ParentMap p;
Chris@16:     };
Chris@16: 
Chris@16:     /* Comparison operator used to choose targets for hooking s.t. vertices 
Chris@16:        that are hooked to are evenly distributed across processors */
Chris@16:     template <typename OwnerMap, typename LocalMap>
Chris@16:     class hashed_vertex_compare
Chris@16:     {
Chris@16:     public:
Chris@16:       hashed_vertex_compare (const OwnerMap& o, const LocalMap& l)
Chris@16:         : owner(o), local(l) { }
Chris@16: 
Chris@16:       template <typename Vertex>
Chris@16:       bool operator() (const Vertex x, const Vertex y) 
Chris@16:       { 
Chris@16:         if (get(local, x) < get(local, y))
Chris@16:           return true;
Chris@16:         else if (get(local, x) == get(local, y))
Chris@16:           return (get(owner, x) < get(owner, y));
Chris@16:         return false;
Chris@16:       }
Chris@16: 
Chris@16:     private:
Chris@16:       OwnerMap   owner;
Chris@16:       LocalMap   local;
Chris@16:     };
Chris@16: 
Chris@16: #ifdef PBGL_EXPLICIT_SYNCH
Chris@16:     template <typename Graph, typename ParentMap, typename VertexList>
Chris@16:     void
Chris@16:     request_parent_map_entries(const Graph& g, ParentMap p,
Chris@16:                                std::vector<VertexList>& parent_requests)
Chris@16:     {
Chris@16:       typedef typename boost::graph::parallel::process_group_type<Graph>
Chris@16:         ::type process_group_type;
Chris@16:       typedef typename process_group_type::process_id_type process_id_type;
Chris@16: 
Chris@16:       typedef typename graph_traits<Graph>::vertex_descriptor
Chris@16:         vertex_descriptor;
Chris@16: 
Chris@16:       process_group_type pg = process_group(g);
Chris@16:       
Chris@16:       /*
Chris@16:         This should probably be send_oob_with_reply, especially when Dave 
Chris@16:         finishes prefetch-batching
Chris@16:       */
Chris@16: 
Chris@16:       // Send root requests
Chris@16:       for (process_id_type i = 0; i < num_processes(pg); ++i) {
Chris@16:         if (!parent_requests[i].empty()) {
Chris@16:           std::vector<vertex_descriptor> reqs(parent_requests[i].begin(),
Chris@16:                                               parent_requests[i].end());
Chris@16:           send(pg, i, req_parents_msg, reqs);
Chris@16:         }
Chris@16:       }
Chris@16:       
Chris@16:       synchronize(pg);
Chris@16:       
Chris@16:       // Receive root requests and reply to them
Chris@16:       while (optional<std::pair<process_id_type, int> > m = probe(pg)) {
Chris@16:         std::vector<vertex_descriptor> requests;
Chris@16:         receive(pg, m->first, m->second, requests);
Chris@16:         for (std::size_t i = 0; i < requests.size(); ++i)
Chris@16:           requests[i] = get(p, requests[i]);
Chris@16:         send(pg, m->first, parents_msg, requests);
Chris@16:       }
Chris@16:       
Chris@16:       synchronize(pg);
Chris@16:       
Chris@16:       // Receive requested parents
Chris@16:       std::vector<vertex_descriptor> responses;
Chris@16:       for (process_id_type i = 0; i < num_processes(pg); ++i) {
Chris@16:         if (!parent_requests[i].empty()) {
Chris@16:           receive(pg, i, parents_msg, responses);
Chris@16:           std::size_t parent_idx = 0;
Chris@16:           for (typename VertexList::iterator v = parent_requests[i].begin();
Chris@16:                v != parent_requests[i].end(); ++v, ++parent_idx)
Chris@16:             put(p, *v, responses[parent_idx]);
Chris@16:         }
Chris@16:       }
Chris@16:     }
Chris@16: #endif
Chris@16:     
Chris@16:     template<typename DistributedGraph, typename ParentMap>
Chris@16:     void
Chris@16:     parallel_connected_components(DistributedGraph& g, ParentMap p)
Chris@16:     {
Chris@16:       using boost::connected_components;
Chris@16: 
Chris@16:       typedef typename graph_traits<DistributedGraph>::adjacency_iterator
Chris@16:         adjacency_iterator;
Chris@16:       typedef typename graph_traits<DistributedGraph>::vertex_descriptor
Chris@16:         vertex_descriptor;
Chris@16: 
Chris@16:       typedef typename boost::graph::parallel::process_group_type<DistributedGraph>
Chris@16:         ::type process_group_type;
Chris@16:       typedef typename process_group_type::process_id_type process_id_type;
Chris@16: 
Chris@16:       using boost::graph::parallel::process_group;
Chris@16: 
Chris@16:       process_group_type pg = process_group(g);
Chris@16:       process_id_type id = process_id(pg);
Chris@16: 
Chris@16:       // TODO (NGE): Should old_roots, roots, and completed_roots be std::list
Chris@16:       adjacency_iterator av1, av2;
Chris@16:       std::vector<vertex_descriptor> old_roots;
Chris@16:       typename std::vector<vertex_descriptor>::iterator liter;
Chris@16:       typename std::vector<vertex_descriptor>::iterator aliter;
Chris@16:       typename std::map<vertex_descriptor,
Chris@16:                         std::vector<vertex_descriptor> > adj;
Chris@16: 
Chris@16:       typedef typename property_map<DistributedGraph, vertex_owner_t>::const_type
Chris@16:         OwnerMap;
Chris@16:       OwnerMap owner = get(vertex_owner, g);
Chris@16:       typedef typename property_map<DistributedGraph, vertex_local_t>::const_type
Chris@16:         LocalMap;
Chris@16:       LocalMap local = get(vertex_local, g);
Chris@16: 
Chris@16:       // We need to hold on to all of the parent pointers
Chris@16:       p.set_max_ghost_cells(0);
Chris@16: 
Chris@16:       //
Chris@16:       // STAGE 1 : Compute local components
Chris@16:       //
Chris@16:       local_subgraph<const DistributedGraph> ls(g);
Chris@16:       typedef typename property_map<local_subgraph<const DistributedGraph>,
Chris@16:                                     vertex_index_t>::type local_index_map_type;
Chris@16:       local_index_map_type local_index = get(vertex_index, ls);
Chris@16: 
Chris@16:       // Compute local connected components
Chris@16:       std::vector<std::size_t> ls_components_vec(num_vertices(ls));
Chris@16:       typedef iterator_property_map<std::vector<std::size_t>::iterator,
Chris@16:                                     local_index_map_type>
Chris@16:         ls_components_map_type;
Chris@16:       ls_components_map_type ls_component(ls_components_vec.begin(),
Chris@16:                                           local_index);
Chris@16:       std::size_t num_comp = connected_components(ls, ls_component);
Chris@16: 
Chris@16:       std::vector<vertex_descriptor> 
Chris@16:         roots(num_comp, graph_traits<DistributedGraph>::null_vertex());
Chris@16: 
Chris@16:       BGL_FORALL_VERTICES_T(v, g, DistributedGraph) {
Chris@16:         size_t component = get(ls_component, v);
Chris@16:         if (roots[component] == graph_traits<DistributedGraph>::null_vertex() ||
Chris@16:             get(local_index, v) < get(local_index, roots[component])) 
Chris@16:           roots[component] = v;
Chris@16:       }
Chris@16: 
Chris@16:       // Set all the local parent pointers
Chris@16:       BGL_FORALL_VERTICES_T(v, g, DistributedGraph) {
Chris@16:         put(p, v, roots[get(ls_component, v)]);
Chris@16:       }
Chris@16: 
Chris@16:       if (num_processes(pg) == 1) return;
Chris@16: 
Chris@16:       // Build adjacency list for all roots
Chris@16:       BGL_FORALL_VERTICES_T(v, g, DistributedGraph) {
Chris@16:         std::vector<vertex_descriptor>& my_adj = adj[get(p, v)];
Chris@16:         for (boost::tie(av1, av2) = adjacent_vertices(v, g);
Chris@16:              av1 != av2; ++av1) {
Chris@16:           if (get(owner, *av1) != id) my_adj.push_back(*av1);
Chris@16:         }
Chris@16:       }
Chris@16: 
Chris@16:       // For all vertices adjacent to a local vertex get p(v)
Chris@16:       for ( liter = roots.begin(); liter != roots.end(); ++liter ) {
Chris@16:         std::vector<vertex_descriptor>& my_adj = adj[*liter];
Chris@16:         for ( aliter = my_adj.begin(); aliter != my_adj.end(); ++aliter )
Chris@16:           request(p, *aliter);
Chris@16:       }
Chris@16:       synchronize(p);
Chris@16: 
Chris@16:       // Update adjacency list at root to make sure all adjacent
Chris@16:       // vertices are roots of remote components
Chris@16:       for ( liter = roots.begin(); liter != roots.end(); ++liter )
Chris@16:         {
Chris@16:           std::vector<vertex_descriptor>& my_adj = adj[*liter];
Chris@16:           for ( aliter = my_adj.begin(); aliter != my_adj.end(); ++aliter )
Chris@16:             *aliter = get(p, *aliter);
Chris@16: 
Chris@16:           my_adj.erase
Chris@16:             (std::remove_if(my_adj.begin(), my_adj.end(),
Chris@16:                        cull_adjacency_list<vertex_descriptor, 
Chris@16:                                            ParentMap>(*liter, p) ),
Chris@16:              my_adj.end());
Chris@16:           // This sort needs to be here to make sure the initial
Chris@16:           // adjacency list is sorted
Chris@16:           std::sort(my_adj.begin(), my_adj.end(), std::less<vertex_descriptor>());
Chris@16:           my_adj.erase(std::unique(my_adj.begin(), my_adj.end()), my_adj.end());
Chris@16:         }
Chris@16: 
Chris@16:       // Get p(v) for the new adjacent roots
Chris@16:       p.clear();
Chris@16:       for ( liter = roots.begin(); liter != roots.end(); ++liter ) {
Chris@16:         std::vector<vertex_descriptor>& my_adj = adj[*liter];
Chris@16:         for ( aliter = my_adj.begin(); aliter != my_adj.end(); ++aliter )
Chris@16:           request(p, *aliter);
Chris@16:       }
Chris@16: #ifdef PBGL_EXPLICIT_SYNCH
Chris@16:       synchronize(p);
Chris@16: #endif
Chris@16: 
Chris@16:       // Lastly, remove roots with no adjacent vertices, this is
Chris@16:       // unnecessary but will speed up sparse graphs
Chris@16:       for ( liter = roots.begin(); liter != roots.end(); /*in loop*/)
Chris@16:         {
Chris@16:           if ( adj[*liter].empty() )
Chris@16:             liter = roots.erase(liter);
Chris@16:           else
Chris@16:             ++liter;
Chris@16:         }
Chris@16: 
Chris@16: #ifdef PBGL_CONSTRUCT_METAGRAPH
Chris@16:       /* TODO: If the number of roots is sufficiently small, we can 
Chris@16:                use a 'problem folding' approach like we do in MST
Chris@16:                to gather all the roots and their adjacencies on one proc
Chris@16:                and solve for the connected components of the meta-graph */
Chris@16:       using boost::parallel::all_reduce;
Chris@16:       std::size_t num_roots = all_reduce(pg, roots.size(), std::plus<std::size_t>());
Chris@16:       if (num_roots < MAX_VERTICES_IN_METAGRAPH) {
Chris@16:         build_local_metagraph(g, p, roots.begin(), roots.end(), adj);
Chris@16:         
Chris@16:         // For each vertex in g, p(v) = p(p(v)), assign parent of leaf
Chris@16:         // vertices from first step to final parent
Chris@16:         BGL_FORALL_VERTICES_T(v, g, DistributedGraph) {
Chris@16:           put(p, v, get(p, get(p, v)));
Chris@16:         }
Chris@16:         
Chris@16:         synchronize(p);
Chris@16:         
Chris@16:         return;
Chris@16:       }
Chris@16: #endif
Chris@16: 
Chris@16:       //
Chris@16:       // Parallel Phase
Chris@16:       //
Chris@16: 
Chris@16:       std::vector<vertex_descriptor> completed_roots;
Chris@16:       hashed_vertex_compare<OwnerMap, LocalMap> v_compare(owner, local);
Chris@16:       bool any_hooked;
Chris@16:       vertex_descriptor new_root;
Chris@16: 
Chris@16:       std::size_t steps = 0;
Chris@16: 
Chris@16:       do {
Chris@16:         ++steps;
Chris@16: 
Chris@16:         // Pull in new parents for hooking phase
Chris@16:         synchronize(p);
Chris@16: 
Chris@16:         //
Chris@16:         // Hooking
Chris@16:         //
Chris@16:         bool hooked = false;
Chris@16:         completed_roots.clear();
Chris@16:         for ( liter = roots.begin(); liter != roots.end(); )
Chris@16:           {
Chris@16:             new_root = graph_traits<DistributedGraph>::null_vertex();
Chris@16:             std::vector<vertex_descriptor>& my_adj = adj[*liter];
Chris@16:             for ( aliter = my_adj.begin(); aliter != my_adj.end(); ++aliter )
Chris@16:               // try to hook to better adjacent vertex
Chris@16:               if ( v_compare( get(p, *aliter), *liter ) )
Chris@16:                 new_root = get(p, *aliter);
Chris@16: 
Chris@16:             if ( new_root != graph_traits<DistributedGraph>::null_vertex() )
Chris@16:               {
Chris@16:                 hooked = true;
Chris@16:                 put(p, *liter, new_root);
Chris@16:                 old_roots.push_back(*liter);
Chris@16:                 completed_roots.push_back(*liter);
Chris@16:                 liter = roots.erase(liter);
Chris@16:               }
Chris@16:             else
Chris@16:               ++liter;
Chris@16:           }
Chris@16: 
Chris@16:         //
Chris@16:         // Pointer jumping, perform until new roots determined
Chris@16:         //
Chris@16: 
Chris@16:         // TODO: Implement cycle reduction rules to reduce this from
Chris@16:         // O(n) to O(log n) [n = cycle length]
Chris@16:         bool all_done;
Chris@16:         std::size_t parent_root_count;
Chris@16: 
Chris@16:         std::size_t double_steps = 0;
Chris@16: 
Chris@16:         do {
Chris@16:           ++double_steps;
Chris@16: #ifndef PBGL_EXPLICIT_SYNCH
Chris@16:           // Get p(p(v)) for all old roots, and p(v) for all current roots
Chris@16:           for ( liter = old_roots.begin(); liter != old_roots.end(); ++liter )
Chris@16:             request(p, get(p, *liter));
Chris@16: 
Chris@16:           synchronize(p);
Chris@16: #else
Chris@16:           // Build root requests
Chris@16:           typedef std::set<vertex_descriptor> VertexSet;
Chris@16:           std::vector<VertexSet> parent_requests(num_processes(pg));
Chris@16:           for ( liter = old_roots.begin(); liter != old_roots.end(); ++liter )
Chris@16:             {
Chris@16:               vertex_descriptor p1 = *liter;
Chris@16:               if (get(owner, p1) != id) parent_requests[get(owner, p1)].insert(p1);
Chris@16:               vertex_descriptor p2 = get(p, p1);
Chris@16:               if (get(owner, p2) != id) parent_requests[get(owner, p2)].insert(p2);
Chris@16:             }
Chris@16: 
Chris@16:           request_parent_map_entries(g, p, parent_requests);
Chris@16: #endif
Chris@16:           // Perform a pointer jumping step on all old roots
Chris@16:           for ( liter = old_roots.begin(); liter != old_roots.end(); ++liter )
Chris@16:               put(p, *liter, get(p, get(p, *liter)));
Chris@16: 
Chris@16:           // make sure the parent of all old roots is itself a root
Chris@16:           parent_root_count = 0;
Chris@16:           for ( liter = old_roots.begin(); liter != old_roots.end(); ++liter )
Chris@16:             if ( get(p, *liter) == get(p, get(p, *liter)) )
Chris@16:               parent_root_count++;
Chris@16: 
Chris@16:           bool done = parent_root_count == old_roots.size();
Chris@16: 
Chris@16:           all_reduce(pg, &done, &done+1, &all_done,
Chris@16:                      std::logical_and<bool>());
Chris@16:         } while ( !all_done );
Chris@16: #ifdef PARALLEL_BGL_DEBUG
Chris@16:         if (id == 0) std::cerr << double_steps << " doubling steps.\n";
Chris@16: #endif
Chris@16:         //
Chris@16:         // Add adjacent vertices of just completed roots to adjacent
Chris@16:         // vertex list at new parent
Chris@16:         //
Chris@16:         typename std::vector<vertex_descriptor> outgoing_edges;
Chris@16:         for ( liter = completed_roots.begin(); liter != completed_roots.end();
Chris@16:               ++liter )
Chris@16:           {
Chris@16:             vertex_descriptor new_parent = get(p, *liter);
Chris@16: 
Chris@16:             if ( get(owner, new_parent) == id )
Chris@16:               {
Chris@16:                 std::vector<vertex_descriptor>& my_adj = adj[new_parent];
Chris@16:                 my_adj.reserve(my_adj.size() + adj[*liter].size());
Chris@16:                 my_adj.insert( my_adj.end(),
Chris@16:                                adj[*liter].begin(), adj[*liter].end() );
Chris@16: #ifdef PBGL_IN_PLACE_MERGE
Chris@16: #ifdef PBGL_SORT_ASSERT
Chris@16:                 BOOST_ASSERT(::boost::detail::is_sorted(my_adj.begin(),
Chris@16:                                                   my_adj.end() - adj[*liter].size(),
Chris@16:                                                   std::less<vertex_descriptor>()));
Chris@16:                 BOOST_ASSERT(::boost::detail::is_sorted(my_adj.end() - adj[*liter].size(),
Chris@16:                                                   my_adj.end(),
Chris@16:                                                   std::less<vertex_descriptor>()));
Chris@16: #endif
Chris@16:                 std::inplace_merge(my_adj.begin(),
Chris@16:                                    my_adj.end() - adj[*liter].size(),
Chris@16:                                    my_adj.end(),
Chris@16:                                    std::less<vertex_descriptor>());
Chris@16: #endif
Chris@16: 
Chris@16: 
Chris@16:               }
Chris@16:             else if ( adj[*liter].begin() != adj[*liter].end() )
Chris@16:               {
Chris@16:                 outgoing_edges.clear();
Chris@16:                 outgoing_edges.reserve(adj[*liter].size() + 1);
Chris@16:                 // First element is the destination of the adjacency list
Chris@16:                 outgoing_edges.push_back(new_parent);
Chris@16:                 outgoing_edges.insert(outgoing_edges.end(),
Chris@16:                                       adj[*liter].begin(), adj[*liter].end() );
Chris@16:                 send(pg, get(owner, new_parent), edges_msg, outgoing_edges);
Chris@16:                 adj[*liter].clear();
Chris@16:               }
Chris@16:           }
Chris@16:         synchronize(pg);
Chris@16: 
Chris@16:         // Receive edges sent by remote nodes and add them to the
Chris@16:         // indicated vertex's adjacency list
Chris@16:         while (optional<std::pair<process_id_type, int> > m
Chris@16:                = probe(pg))
Chris@16:           {
Chris@16:             std::vector<vertex_descriptor> incoming_edges;
Chris@16:             receive(pg, m->first, edges_msg, incoming_edges);
Chris@16:             typename std::vector<vertex_descriptor>::iterator aviter
Chris@16:               = incoming_edges.begin();
Chris@16:             ++aviter;
Chris@16: 
Chris@16:             std::vector<vertex_descriptor>& my_adj = adj[incoming_edges[0]];
Chris@16: 
Chris@16:             my_adj.reserve(my_adj.size() + incoming_edges.size() - 1);
Chris@16:             my_adj.insert( my_adj.end(), aviter, incoming_edges.end() );
Chris@16: 
Chris@16: #ifdef PBGL_IN_PLACE_MERGE
Chris@16:             std::size_t num_incoming_edges = incoming_edges.size();
Chris@16: #ifdef PBGL_SORT_ASSERT
Chris@16:             BOOST_ASSERT(::boost::detail::is_sorted(my_adj.begin(),
Chris@16:                                               my_adj.end() - (num_incoming_edges-1),
Chris@16:                                               std::less<vertex_descriptor>()));
Chris@16:             BOOST_ASSERT(::boost::detail::is_sorted(my_adj.end() - (num_incoming_edges-1),
Chris@16:                                               my_adj.end(),
Chris@16:                                               std::less<vertex_descriptor>()));
Chris@16: #endif
Chris@16:             std::inplace_merge(my_adj.begin(),
Chris@16:                                my_adj.end() - (num_incoming_edges - 1),
Chris@16:                                my_adj.end(),
Chris@16:                                std::less<vertex_descriptor>());
Chris@16: #endif
Chris@16: 
Chris@16:           }
Chris@16: 
Chris@16: 
Chris@16:         // Remove any adjacent vertices that are in the same component
Chris@16:         // as a root from that root's list
Chris@16:         for ( liter = roots.begin(); liter != roots.end(); ++liter )
Chris@16:           {
Chris@16:             // We can probably get away without sorting and removing
Chris@16:             // duplicates Though sorting *may* cause root
Chris@16:             // determination to occur faster by choosing the root with
Chris@16:             // the most potential to hook to at each step
Chris@16:             std::vector<vertex_descriptor>& my_adj = adj[*liter];
Chris@16:             my_adj.erase
Chris@16:               (std::remove_if(my_adj.begin(), my_adj.end(),
Chris@16:                          cull_adjacency_list<vertex_descriptor,
Chris@16:                                              ParentMap>(*liter, p) ),
Chris@16:                my_adj.end());
Chris@16: #ifndef PBGL_IN_PLACE_MERGE
Chris@16:             std::sort(my_adj.begin(), my_adj.end(),
Chris@16:                  std::less<vertex_descriptor>() );
Chris@16: #endif
Chris@16:             my_adj.erase(std::unique(my_adj.begin(), my_adj.end()), my_adj.end());
Chris@16:           }
Chris@16: 
Chris@16:         // Reduce result of empty root list test
Chris@16:         all_reduce(pg, &hooked, &hooked+1, &any_hooked,
Chris@16:                    std::logical_or<bool>());
Chris@16:       } while ( any_hooked );
Chris@16: #ifdef PARALLEL_BGL_DEBUG
Chris@16:       if (id == 0) std::cerr << steps << " iterations.\n";
Chris@16: #endif
Chris@16:       //
Chris@16:       // Finalize
Chris@16:       //
Chris@16: 
Chris@16:       // For each vertex in g, p(v) = p(p(v)), assign parent of leaf
Chris@16:       // vertices from first step to final parent
Chris@16:       BGL_FORALL_VERTICES_T(v, g, DistributedGraph) {
Chris@16:         put(p, v, get(p, get(p, v)));
Chris@16:       }
Chris@16:       
Chris@16:       synchronize(p);
Chris@16:     }
Chris@16: 
Chris@16:   } // end namespace cc_detail
Chris@16: 
Chris@16:   template<typename Graph, typename ParentMap, typename ComponentMap>
Chris@16:   typename property_traits<ComponentMap>::value_type
Chris@16:   number_components_from_parents(const Graph& g, ParentMap p, ComponentMap c)
Chris@16:   {
Chris@16:     typedef typename graph_traits<Graph>::vertex_descriptor
Chris@16:       vertex_descriptor;
Chris@16:     typedef typename boost::graph::parallel::process_group_type<Graph>::type
Chris@16:       process_group_type;
Chris@16:     typedef typename property_traits<ComponentMap>::value_type
Chris@16:       ComponentMapType;
Chris@16: 
Chris@16:     process_group_type pg = process_group(g);
Chris@16: 
Chris@16:     /* Build list of roots */
Chris@16:     std::vector<vertex_descriptor> my_roots, all_roots;
Chris@16: 
Chris@16:     BGL_FORALL_VERTICES_T(v, g, Graph) {
Chris@16:       if( std::find( my_roots.begin(), my_roots.end(), get(p, v) )
Chris@16:           == my_roots.end() )
Chris@16:         my_roots.push_back( get(p, v) );
Chris@16:     }
Chris@16: 
Chris@16:     all_gather(pg, my_roots.begin(), my_roots.end(), all_roots);
Chris@16: 
Chris@16:     /* Number components */
Chris@16:     std::map<vertex_descriptor, ComponentMapType> comp_numbers;
Chris@16:     ComponentMapType c_num = 0;
Chris@16: 
Chris@16:     // Compute component numbers
Chris@16:     for (std::size_t i = 0; i < all_roots.size(); i++ )
Chris@16:       if ( comp_numbers.count(all_roots[i]) == 0 )
Chris@16:         comp_numbers[all_roots[i]] = c_num++;
Chris@16: 
Chris@16:     // Broadcast component numbers
Chris@16:     BGL_FORALL_VERTICES_T(v, g, Graph) {
Chris@16:       put( c, v, comp_numbers[get(p, v)] );
Chris@16:     }
Chris@16: 
Chris@16:     // Broadcast number of components
Chris@16:     if (process_id(pg) == 0) {
Chris@16:       typedef typename process_group_type::process_size_type
Chris@16:         process_size_type;
Chris@16:       for (process_size_type dest = 1, n = num_processes(pg);
Chris@16:            dest != n; ++dest)
Chris@16:         send(pg, dest, 0, c_num);
Chris@16:     }
Chris@16:     synchronize(pg);
Chris@16: 
Chris@16:     if (process_id(pg) != 0) receive(pg, 0, 0, c_num);
Chris@16: 
Chris@16:     synchronize(c);
Chris@16: 
Chris@16:     return c_num;
Chris@16:   }
Chris@16: 
Chris@16:   template<typename Graph, typename ParentMap>
Chris@16:   int
Chris@16:   number_components_from_parents(const Graph& g, ParentMap p, 
Chris@16:                                  dummy_property_map)
Chris@16:   {
Chris@16:     using boost::parallel::all_reduce;
Chris@16: 
Chris@16:     // Count local roots.
Chris@16:     int num_roots = 0;
Chris@16:     BGL_FORALL_VERTICES_T(v, g, Graph)
Chris@16:       if (get(p, v) == v) ++num_roots;
Chris@16:     return all_reduce(g.process_group(), num_roots, std::plus<int>());
Chris@16:   }
Chris@16: 
Chris@16:   template<typename Graph, typename ComponentMap, typename ParentMap>
Chris@16:   typename property_traits<ComponentMap>::value_type
Chris@16:   connected_components
Chris@16:     (const Graph& g, ComponentMap c, ParentMap p
Chris@16:      BOOST_GRAPH_ENABLE_IF_MODELS_PARM(Graph, distributed_graph_tag))
Chris@16:   {
Chris@16:     cc_detail::parallel_connected_components(g, p);
Chris@16:     return number_components_from_parents(g, p, c);
Chris@16:   }
Chris@16: 
Chris@16:   /* Construct ParentMap by default */
Chris@16:   template<typename Graph, typename ComponentMap>
Chris@16:   typename property_traits<ComponentMap>::value_type
Chris@16:   connected_components
Chris@16:     ( const Graph& g, ComponentMap c
Chris@16:       BOOST_GRAPH_ENABLE_IF_MODELS_PARM(Graph, distributed_graph_tag) )
Chris@16:   {
Chris@16:     typedef typename graph_traits<Graph>::vertex_descriptor vertex_descriptor;
Chris@16: 
Chris@16:     std::vector<vertex_descriptor> x(num_vertices(g));
Chris@16: 
Chris@16:     return connected_components
Chris@16:              (g, c,
Chris@16:               make_iterator_property_map(x.begin(), get(vertex_index, g)));
Chris@16:   }
Chris@16: } // end namespace distributed
Chris@16: 
Chris@16: using distributed::connected_components;
Chris@16: } // end namespace graph
Chris@16: 
Chris@16: using graph::distributed::connected_components;
Chris@16: } // end namespace boost
Chris@16: 
Chris@16: #endif // BOOST_GRAPH_PARALLEL_CC_HPP