1 // Copyright 2005 The Trustees of Indiana University.
3 // Use, modification and distribution is subject to the Boost Software
4 // License, Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
5 // http://www.boost.org/LICENSE_1_0.txt)
7 // Authors: Douglas Gregor
10 // An implementation of Walter Hohberg's distributed biconnected
11 // components algorithm, from:
13 // Walter Hohberg. How to Find Biconnected Components in Distributed
14 // Networks. J. Parallel Distrib. Comput., 9(4):374-386, 1990.
16 #ifndef BOOST_GRAPH_DISTRIBUTED_HOHBERG_BICONNECTED_COMPONENTS_HPP
17 #define BOOST_GRAPH_DISTRIBUTED_HOHBERG_BICONNECTED_COMPONENTS_HPP
19 #ifndef BOOST_GRAPH_USE_MPI
20 #error "Parallel BGL files should not be included unless <boost/graph/use_mpi.hpp> has been included"
23 /* You can define PBGL_HOHBERG_DEBUG to an integer value (1, 2, or 3)
24 * to enable debugging information. 1 includes only the phases of the
25 * algorithm and messages as their are received. 2 and 3 add
26 * additional levels of detail about internal data structures related
27 * to the algorithm itself.
29 * #define PBGL_HOHBERG_DEBUG 1
32 #include <boost/graph/graph_traits.hpp>
33 #include <boost/graph/parallel/container_traits.hpp>
34 #include <boost/graph/parallel/process_group.hpp>
35 #include <boost/static_assert.hpp>
36 #include <boost/mpi/operations.hpp>
37 #include <boost/type_traits/is_convertible.hpp>
38 #include <boost/graph/graph_concepts.hpp>
39 #include <boost/graph/iteration_macros.hpp>
40 #include <boost/optional.hpp>
41 #include <utility> // for std::pair
42 #include <boost/assert.hpp>
43 #include <algorithm> // for std::find, std::mismatch
45 #include <boost/graph/parallel/algorithm.hpp>
46 #include <boost/graph/distributed/connected_components.hpp>
47 #include <boost/concept/assert.hpp>
49 namespace boost { namespace graph { namespace distributed {
51 namespace hohberg_detail {
53 /* A header for the PATH message, stating which edge the message
54 is coming on and how many vertices will be following. The data
55 structure communicated will be a path_header. */
57 /* A message containing the vertices that make up a path. It will
58 always follow a msg_path_header and will contain vertex
61 /* A header for the TREE message, stating the value of gamma and
62 the number of vertices to come in the following
65 /* A message containing the vertices that make up the set of
66 vertices in the same bicomponent as the sender. It will always
67 follow a msg_tree_header and will contain vertex descriptors,
70 /* Provides a name for the biconnected component of the edge. */
74 // Payload for a msg_path_header message.
75 template<typename EdgeDescriptor>
78 // The edge over which the path is being communicated
81 // The length of the path, i.e., the number of vertex descriptors
82 // that will be coming in the next msg_path_vertices message.
83 std::size_t path_length;
85 template<typename Archiver>
86 void serialize(Archiver& ar, const unsigned int /*version*/)
88 ar & edge & path_length;
92 // Payload for a msg_tree_header message.
93 template<typename Vertex, typename Edge>
96 // The edge over which the tree is being communicated
99 // Gamma, which is the eta of the sender.
102 // The length of the list of vertices in the bicomponent, i.e.,
103 // the number of vertex descriptors that will be coming in the
104 // next msg_tree_vertices message.
105 std::size_t bicomp_length;
107 template<typename Archiver>
108 void serialize(Archiver& ar, const unsigned int /*version*/)
110 ar & edge & gamma & bicomp_length;
114 // Payload for the msg_name message.
115 template<typename EdgeDescriptor>
118 // The edge being communicated and named.
121 // The 0-based name of the component
124 template<typename Archiver>
125 void serialize(Archiver& ar, const unsigned int /*version*/)
131 /* Computes the branch point between two paths. The branch point is
132 the last position at which both paths are equivalent, beyond
133 which the paths diverge. Both paths must have length > 0 and the
134 initial elements of the paths must be equal. This is guaranteed
135 in Hohberg's algorithm because all paths start at the
136 leader. Returns the value at the branch point. */
138 T branch_point(const std::vector<T>& p1, const std::vector<T>& p2)
140 BOOST_ASSERT(!p1.empty());
141 BOOST_ASSERT(!p2.empty());
142 BOOST_ASSERT(p1.front() == p2.front());
144 typedef typename std::vector<T>::const_iterator iterator;
146 iterator mismatch_pos;
147 if (p1.size() <= p2.size())
148 mismatch_pos = std::mismatch(p1.begin(), p1.end(), p2.begin()).first;
150 mismatch_pos = std::mismatch(p2.begin(), p2.end(), p1.begin()).first;
152 return *mismatch_pos;
155 /* Computes the infimum of vertices a and b in the given path. The
156 infimum is the largest element that is on the paths from a to the
157 root and from b to the root. */
159 T infimum(const std::vector<T>& parent_path, T a, T b)
163 typedef typename std::vector<T>::const_iterator iterator;
164 iterator first = parent_path.begin(), last = parent_path.end();
166 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 2
167 std::cerr << "infimum(";
168 for (iterator i = first; i != last; ++i) {
169 if (i != first) std::cerr << ' ';
170 std::cerr << local(*i) << '@' << owner(*i);
172 std::cerr << ", " << local(a) << '@' << owner(a) << ", "
173 << local(b) << '@' << owner(b) << ") = ";
177 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 2
178 std::cerr << local(a) << '@' << owner(a) << std::endl;
183 // Try to find a or b, whichever is closest to the end
186 // If we match b, swap the two so we'll be looking for b later.
187 if (*last == b) { swap(a,b); break; }
190 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 2
191 std::cerr << local(*first) << '@' << owner(*first) << std::endl;
198 // Try to find b (which may originally have been a)
201 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 2
202 std::cerr << local(*first) << '@' << owner(*first) << std::endl;
209 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 2
210 std::cerr << local(*last) << '@' << owner(*last) << std::endl;
212 // We've found b; it's the infimum.
215 } // end namespace hohberg_detail
217 /* A message that will be stored for each edge by Hohberg's algorithm. */
218 template<typename Graph>
219 struct hohberg_message
221 typedef typename graph_traits<Graph>::vertex_descriptor Vertex;
222 typedef typename graph_traits<Graph>::edge_descriptor Edge;
224 // Assign from a path message
225 void assign(const std::vector<Vertex>& path)
227 gamma = graph_traits<Graph>::null_vertex();
228 path_or_bicomp = path;
231 // Assign from a tree message
232 void assign(Vertex gamma, const std::vector<Vertex>& in_same_bicomponent)
235 path_or_bicomp = in_same_bicomponent;
238 bool is_path() const { return gamma == graph_traits<Graph>::null_vertex(); }
239 bool is_tree() const { return gamma != graph_traits<Graph>::null_vertex(); }
241 /// The "gamma" of a tree message, or null_vertex() for a path message
244 // Either the path for a path message or the in_same_bicomponent
245 std::vector<Vertex> path_or_bicomp;
249 /* An abstraction of a vertex processor in Hohberg's algorithm. The
250 hohberg_vertex_processor class is responsible for processing
251 messages transmitted to it via its edges. */
252 template<typename Graph>
253 class hohberg_vertex_processor
255 typedef typename graph_traits<Graph>::vertex_descriptor Vertex;
256 typedef typename graph_traits<Graph>::edge_descriptor Edge;
257 typedef typename graph_traits<Graph>::degree_size_type degree_size_type;
258 typedef typename graph_traits<Graph>::edges_size_type edges_size_type;
259 typedef typename boost::graph::parallel::process_group_type<Graph>::type
261 typedef std::vector<Vertex> path_t;
262 typedef typename path_t::iterator path_iterator;
265 hohberg_vertex_processor()
267 parent(graph_traits<Graph>::null_vertex()),
268 eta(graph_traits<Graph>::null_vertex())
272 // Called to initialize a leader in the algorithm, which involves
273 // sending out the initial path messages and being ready to receive
275 void initialize_leader(Vertex alpha, const Graph& g);
277 /// Handle a path message on edge e. The path will be destroyed by
280 operator()(Edge e, path_t& path, const Graph& g);
282 /// Handle a tree message on edge e. in_same_bicomponent will be
283 /// destroyed by this operation.
285 operator()(Edge e, Vertex gamma, path_t& in_same_bicomponent,
288 // Handle a name message.
289 void operator()(Edge e, edges_size_type name, const Graph& g);
291 // Retrieve the phase
292 unsigned char get_phase() const { return phase; }
294 // Start the naming phase. The current phase must be 3 (echo), and
295 // the offset contains the offset at which this processor should
296 // begin when labelling its bicomponents. The offset is just a
297 // parallel prefix sum of the number of bicomponents in each
298 // processor that precedes it (globally).
300 start_naming_phase(Vertex alpha, const Graph& g, edges_size_type offset);
302 /* Determine the number of bicomponents that we will be naming
305 edges_size_type num_starting_bicomponents(Vertex alpha, const Graph& g);
307 // Fill in the edge component map with biconnected component
309 template<typename ComponentMap>
310 void fill_edge_map(Vertex alpha, const Graph& g, ComponentMap& component);
313 /* Start the echo phase (phase 3) where we propagate information up
315 void echo_phase(Vertex alpha, const Graph& g);
317 /* Retrieve the index of edge in the out-edges list of target(e, g). */
318 std::size_t get_edge_index(Edge e, const Graph& g);
320 /* Retrieve the index of the edge incidence on v in the out-edges
322 std::size_t get_incident_edge_index(Vertex u, Vertex v, const Graph& g);
324 /* Keeps track of which phase of the algorithm we are in. There are
325 * four phases plus the "finished" phase:
327 * 1) Building the spanning tree
328 * 2) Discovering cycles
329 * 3) Echoing back up the spanning tree
330 * 4) Labelling biconnected components
335 /* The parent of this vertex in the spanning tree. This value will
336 be graph_traits<Graph>::null_vertex() for the leader. */
339 /* The farthest ancestor up the tree that resides in the same
340 biconnected component as we do. This information is approximate:
341 we might not know about the actual farthest ancestor, but this is
342 the farthest one we've seen so far. */
345 /* The number of edges that have not yet transmitted any messages to
346 us. This starts at the degree of the vertex and decreases as we
347 receive messages. When this counter hits zero, we're done with
348 the second phase of the algorithm. In Hohberg's paper, the actual
349 remaining edge set E is stored with termination when all edges
350 have been removed from E, but we only need to detect termination
351 so the set E need not be explicitly represented. */
352 degree_size_type num_edges_not_transmitted;
354 /* The path from the root of the spanning tree to this vertex. This
355 vector will always store only the parts of the path leading up to
356 this vertex, and not the vertex itself. Thus, it will be empty
358 std::vector<Vertex> path_from_root;
360 /* Structure containing all of the extra data we need to keep around
361 PER EDGE. This information can not be contained within a property
362 map, because it can't be shared among vertices without breaking
363 the algorithm. Decreasing the size of this structure will drastically */
366 hohberg_message<Graph> msg;
367 std::vector<Vertex> M;
369 degree_size_type partition;
372 /* Data for each edge in the graph. This structure will be indexed
373 by the position of the edge in the out_edges() list. */
374 std::vector<per_edge_data> edge_data;
376 /* The mapping from local partition numbers (0..n-1) to global
377 partition numbers. */
378 std::vector<edges_size_type> local_to_global_partitions;
380 friend class boost::serialization::access;
382 // We cannot actually serialize a vertex processor, nor would we
383 // want to. However, the fact that we're putting instances into a
384 // distributed_property_map means that we need to have a serialize()
385 // function available.
386 template<typename Archiver>
387 void serialize(Archiver&, const unsigned int /*version*/)
393 template<typename Graph>
395 hohberg_vertex_processor<Graph>::initialize_leader(Vertex alpha,
398 using namespace hohberg_detail;
400 ProcessGroup pg = process_group(g);
402 typename property_map<Graph, vertex_owner_t>::const_type
403 owner = get(vertex_owner, g);
405 path_header<Edge> header;
406 header.path_length = 1;
407 BGL_FORALL_OUTEDGES_T(alpha, e, g, Graph) {
409 send(pg, get(owner, target(e, g)), msg_path_header, header);
410 send(pg, get(owner, target(e, g)), msg_path_vertices, alpha);
413 num_edges_not_transmitted = degree(alpha, g);
414 edge_data.resize(num_edges_not_transmitted);
418 template<typename Graph>
420 hohberg_vertex_processor<Graph>::operator()(Edge e, path_t& path,
423 using namespace hohberg_detail;
425 typename property_map<Graph, vertex_owner_t>::const_type
426 owner = get(vertex_owner, g);
428 #ifdef PBGL_HOHBERG_DEBUG
429 // std::cerr << local(source(e, g)) << '@' << owner(source(e, g)) << " -> "
430 // << local(target(e, g)) << '@' << owner(target(e, g)) << ": path(";
431 // for (std::size_t i = 0; i < path.size(); ++i) {
432 // if (i > 0) std::cerr << ' ';
433 // std::cerr << local(path[i]) << '@' << owner(path[i]);
435 std::cerr << "), phase = " << (int)phase << std::endl;
438 // Get access to edge-specific data
439 if (edge_data.empty())
440 edge_data.resize(degree(target(e, g), g));
441 per_edge_data& edata = edge_data[get_edge_index(e, g)];
443 // Record the message. We'll need it in phase 3.
444 edata.msg.assign(path);
446 // Note: "alpha" refers to the vertex "processor" receiving the
448 Vertex alpha = target(e, g);
453 num_edges_not_transmitted = degree(alpha, g) - 1;
454 edata.is_tree_edge = true;
455 parent = path.back();
457 edata.M.clear(); edata.M.push_back(parent);
459 // Broadcast the path from the root to our potential children in
460 // the spanning tree.
461 path.push_back(alpha);
462 path_header<Edge> header;
463 header.path_length = path.size();
464 ProcessGroup pg = process_group(g);
465 BGL_FORALL_OUTEDGES_T(alpha, oe, g, Graph) {
466 // Skip the tree edge we just received
467 if (target(oe, g) != source(e, g)) {
469 send(pg, get(owner, target(oe, g)), msg_path_header, header);
470 send(pg, get(owner, target(oe, g)), msg_path_vertices, &path[0],
476 // Swap the old path in, to save some extra copying. Nobody
477 path_from_root.swap(path);
479 // Once we have received our place in the spanning tree, move on
483 // If we only had only edge, skip to phase 3.
484 if (num_edges_not_transmitted == 0)
485 echo_phase(alpha, g);
491 --num_edges_not_transmitted;
492 edata.is_tree_edge = false;
494 // Determine if alpha (our vertex) is in the path
495 path_iterator pos = std::find(path.begin(), path.end(), alpha);
496 if (pos != path.end()) {
497 // Case A: There is a cycle alpha beta ... gamma alpha
498 // M(e) <- {beta, gammar}
501 // If pos == path.end(), we have a self-loop
502 if (pos != path.end()) {
504 edata.M.push_back(*pos);
507 // If pos == path.end(), we have a self-loop or beta == gamma
508 // (parallel edge). Otherwise, add gamma.
509 if (pos != path.end()) edata.M.push_back(path.back());
511 // Case B: There is a cycle but we haven't seen alpha yet.
512 // M(e) = {parent, path.back()}
514 edata.M.push_back(path.back());
515 if (parent != path.back()) edata.M.push_back(parent);
517 // eta = inf(eta, bra(pi_t, pi))
518 eta = infimum(path_from_root, eta, branch_point(path_from_root, path));
520 if (num_edges_not_transmitted == 0)
521 echo_phase(alpha, g);
526 // std::cerr << "Phase is " << int(phase) << "\n";
531 template<typename Graph>
533 hohberg_vertex_processor<Graph>::operator()(Edge e, Vertex gamma,
534 path_t& in_same_bicomponent,
537 using namespace hohberg_detail;
539 #ifdef PBGL_HOHBERG_DEBUG
540 std::cerr << local(source(e, g)) << '@' << owner(source(e, g)) << " -> "
541 << local(target(e, g)) << '@' << owner(target(e, g)) << ": tree("
542 << local(gamma) << '@' << owner(gamma) << ", ";
543 for (std::size_t i = 0; i < in_same_bicomponent.size(); ++i) {
544 if (i > 0) std::cerr << ' ';
545 std::cerr << local(in_same_bicomponent[i]) << '@'
546 << owner(in_same_bicomponent[i]);
548 std::cerr << ", " << local(source(e, g)) << '@' << owner(source(e, g))
549 << "), phase = " << (int)phase << std::endl;
552 // Get access to edge-specific data
553 per_edge_data& edata = edge_data[get_edge_index(e, g)];
555 // Record the message. We'll need it in phase 3.
556 edata.msg.assign(gamma, in_same_bicomponent);
558 // Note: "alpha" refers to the vertex "processor" receiving the
560 Vertex alpha = target(e, g);
561 Vertex beta = source(e, g);
565 --num_edges_not_transmitted;
566 edata.is_tree_edge = true;
568 if (gamma == alpha) {
570 edata.M.swap(in_same_bicomponent);
574 edata.M.push_back(parent);
575 if (beta != parent) edata.M.push_back(beta);
576 eta = infimum(path_from_root, eta, gamma);
578 if (num_edges_not_transmitted == 0)
579 echo_phase(alpha, g);
587 template<typename Graph>
589 hohberg_vertex_processor<Graph>::operator()(Edge e, edges_size_type name,
592 using namespace hohberg_detail;
594 #ifdef PBGL_HOHBERG_DEBUG
595 std::cerr << local(source(e, g)) << '@' << owner(source(e, g)) << " -> "
596 << local(target(e, g)) << '@' << owner(target(e, g)) << ": name("
597 << name << "), phase = " << (int)phase << std::endl;
600 BOOST_ASSERT(phase == 4);
602 typename property_map<Graph, vertex_owner_t>::const_type
603 owner = get(vertex_owner, g);
605 // Send name messages along the spanning tree edges that are in the
606 // same bicomponent as the edge to our parent.
607 ProcessGroup pg = process_group(g);
609 Vertex alpha = target(e, g);
612 BGL_FORALL_OUTEDGES_T(alpha, e, g, Graph) {
613 per_edge_data& edata = edge_data[idx++];
614 if (edata.is_tree_edge
615 && find(edata.M.begin(), edata.M.end(), parent) != edata.M.end()
616 && target(e, g) != parent) {
617 // Notify our children in the spanning tree of this name
618 name_header<Edge> header;
621 send(pg, get(owner, target(e, g)), msg_name, header);
622 } else if (target(e, g) == parent) {
623 // Map from local partition numbers to global bicomponent numbers
624 local_to_global_partitions[edata.partition] = name;
632 template<typename Graph>
633 typename hohberg_vertex_processor<Graph>::edges_size_type
634 hohberg_vertex_processor<Graph>::
635 num_starting_bicomponents(Vertex alpha, const Graph& g)
637 edges_size_type not_mapped = (std::numeric_limits<edges_size_type>::max)();
639 edges_size_type result = 0;
641 BGL_FORALL_OUTEDGES_T(alpha, e, g, Graph) {
642 per_edge_data& edata = edge_data[idx++];
643 if (edata.is_tree_edge
644 && find(edata.M.begin(), edata.M.end(), parent) == edata.M.end()) {
645 // Map from local partition numbers to global bicomponent numbers
646 if (local_to_global_partitions[edata.partition] == not_mapped)
647 local_to_global_partitions[edata.partition] = result++;
651 #ifdef PBGL_HOHBERG_DEBUG
652 std::cerr << local(alpha) << '@' << owner(alpha) << " has " << result
653 << " bicomponents originating at it." << std::endl;
659 template<typename Graph>
660 template<typename ComponentMap>
662 hohberg_vertex_processor<Graph>::
663 fill_edge_map(Vertex alpha, const Graph& g, ComponentMap& component)
666 BGL_FORALL_OUTEDGES_T(alpha, e, g, Graph) {
667 per_edge_data& edata = edge_data[idx++];
668 local_put(component, e, local_to_global_partitions[edata.partition]);
670 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 2
671 std::cerr << "component("
672 << local(source(e, g)) << '@' << owner(source(e, g)) << " -> "
673 << local(target(e, g)) << '@' << owner(target(e, g)) << ") = "
674 << local_to_global_partitions[edata.partition]
675 << " (partition = " << edata.partition << " of "
676 << local_to_global_partitions.size() << ")" << std::endl;
681 template<typename Graph>
683 hohberg_vertex_processor<Graph>::
684 start_naming_phase(Vertex alpha, const Graph& g, edges_size_type offset)
686 using namespace hohberg_detail;
688 BOOST_ASSERT(phase == 4);
690 typename property_map<Graph, vertex_owner_t>::const_type
691 owner = get(vertex_owner, g);
693 // Send name messages along the spanning tree edges of the
694 // components that we get to number.
695 ProcessGroup pg = process_group(g);
697 bool has_more_children_to_name = false;
699 // Map from local partition numbers to global bicomponent numbers
700 edges_size_type not_mapped = (std::numeric_limits<edges_size_type>::max)();
701 for (std::size_t i = 0; i < local_to_global_partitions.size(); ++i) {
702 if (local_to_global_partitions[i] != not_mapped)
703 local_to_global_partitions[i] += offset;
707 BGL_FORALL_OUTEDGES_T(alpha, e, g, Graph) {
708 per_edge_data& edata = edge_data[idx++];
709 if (edata.is_tree_edge
710 && find(edata.M.begin(), edata.M.end(), parent) == edata.M.end()) {
711 // Notify our children in the spanning tree of this new name
712 name_header<Edge> header;
714 header.name = local_to_global_partitions[edata.partition];
715 send(pg, get(owner, target(e, g)), msg_name, header);
716 } else if (edata.is_tree_edge) {
717 has_more_children_to_name = true;
719 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 2
720 std::cerr << "M[" << local(source(e, g)) << '@' << owner(source(e, g))
721 << " -> " << local(target(e, g)) << '@' << owner(target(e, g))
723 for (std::size_t i = 0; i < edata.M.size(); ++i) {
724 std::cerr << local(edata.M[i]) << '@' << owner(edata.M[i]) << ' ';
726 std::cerr << std::endl;
730 // See if we're done.
731 if (!has_more_children_to_name)
736 template<typename Graph>
738 hohberg_vertex_processor<Graph>::echo_phase(Vertex alpha, const Graph& g)
740 using namespace hohberg_detail;
742 typename property_map<Graph, vertex_owner_t>::const_type
743 owner = get(vertex_owner, g);
745 /* We're entering the echo phase. */
748 if (parent != graph_traits<Graph>::null_vertex()) {
751 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 1
752 std::cerr << local(alpha) << '@' << owner(alpha) << " echo: parent = "
753 << local(parent) << '@' << owner(parent) << ", eta = "
754 << local(eta) << '@' << owner(eta) << ", Gamma = ";
757 std::vector<Vertex> bicomp;
758 std::size_t e_index = 0;
759 BGL_FORALL_OUTEDGES_T(alpha, e, g, Graph) {
760 if (target(e, g) == parent && parent == eta) {
762 if (find(bicomp.begin(), bicomp.end(), alpha) == bicomp.end()) {
763 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 1
764 std::cerr << local(alpha) << '@' << owner(alpha) << ' ';
766 bicomp.push_back(alpha);
769 if (target(e, g) == parent) edge_to_parent = e;
771 per_edge_data& edata = edge_data[e_index];
773 if (edata.msg.is_path()) {
774 path_iterator pos = std::find(edata.msg.path_or_bicomp.begin(),
775 edata.msg.path_or_bicomp.end(),
777 if (pos != edata.msg.path_or_bicomp.end()) {
779 if (pos != edata.msg.path_or_bicomp.end()
780 && find(bicomp.begin(), bicomp.end(), *pos) == bicomp.end()) {
781 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 1
782 std::cerr << local(*pos) << '@' << owner(*pos) << ' ';
784 bicomp.push_back(*pos);
787 } else if (edata.msg.is_tree() && edata.msg.gamma == eta) {
788 for (path_iterator i = edata.msg.path_or_bicomp.begin();
789 i != edata.msg.path_or_bicomp.end(); ++i) {
790 if (find(bicomp.begin(), bicomp.end(), *i) == bicomp.end()) {
791 #if defined(PBGL_HOHBERG_DEBUG) && PBGL_HOHBERG_DEBUG > 1
792 std::cerr << local(*i) << '@' << owner(*i) << ' ';
794 bicomp.push_back(*i);
801 #ifdef PBGL_HOHBERG_DEBUG
802 std::cerr << std::endl;
805 // Send tree(eta, bicomp) to parent
806 tree_header<Vertex, Edge> header;
807 header.edge = edge_to_parent;
809 header.bicomp_length = bicomp.size();
810 ProcessGroup pg = process_group(g);
811 send(pg, get(owner, parent), msg_tree_header, header);
812 send(pg, get(owner, parent), msg_tree_vertices, &bicomp[0],
813 header.bicomp_length);
816 // Compute the partition of edges such that iff two edges e1 and e2
817 // are in different subsets then M(e1) is disjoint from M(e2).
819 // Start by putting each edge in a different partition
820 std::vector<degree_size_type> parent_vec(edge_data.size());
821 degree_size_type idx = 0;
822 for (idx = 0; idx < edge_data.size(); ++idx)
823 parent_vec[idx] = idx;
825 // Go through each edge e, performing a union() on the edges
826 // incident on all vertices in M[e].
828 BGL_FORALL_OUTEDGES_T(alpha, e, g, Graph) {
829 per_edge_data& edata = edge_data[idx++];
831 // Compute union of vertices in M
832 if (!edata.M.empty()) {
833 degree_size_type e1 = get_incident_edge_index(alpha, edata.M.front(), g);
834 while (parent_vec[e1] != e1) e1 = parent_vec[e1];
836 for (std::size_t i = 1; i < edata.M.size(); ++i) {
837 degree_size_type e2 = get_incident_edge_index(alpha, edata.M[i], g);
838 while (parent_vec[e2] != e2) e2 = parent_vec[e2];
844 edges_size_type not_mapped = (std::numeric_limits<edges_size_type>::max)();
846 // Determine the number of partitions
847 for (idx = 0; idx < parent_vec.size(); ++idx) {
848 if (parent_vec[idx] == idx) {
849 edge_data[idx].partition = local_to_global_partitions.size();
850 local_to_global_partitions.push_back(not_mapped);
854 // Assign partition numbers to each edge
855 for (idx = 0; idx < parent_vec.size(); ++idx) {
856 degree_size_type rep = parent_vec[idx];
857 while (rep != parent_vec[rep]) rep = parent_vec[rep];
858 edge_data[idx].partition = edge_data[rep].partition;
861 // Enter the naming phase (but don't send anything yet).
865 template<typename Graph>
867 hohberg_vertex_processor<Graph>::get_edge_index(Edge e, const Graph& g)
869 std::size_t result = 0;
870 BGL_FORALL_OUTEDGES_T(target(e, g), oe, g, Graph) {
871 if (source(e, g) == target(oe, g)) return result;
877 template<typename Graph>
879 hohberg_vertex_processor<Graph>::get_incident_edge_index(Vertex u, Vertex v,
882 std::size_t result = 0;
883 BGL_FORALL_OUTEDGES_T(u, e, g, Graph) {
884 if (target(e, g) == v) return result;
890 template<typename Graph, typename InputIterator, typename ComponentMap,
891 typename VertexProcessorMap>
892 typename graph_traits<Graph>::edges_size_type
893 hohberg_biconnected_components
895 ComponentMap component,
896 InputIterator first, InputIterator last,
897 VertexProcessorMap vertex_processor)
899 using namespace boost::graph::parallel;
900 using namespace hohberg_detail;
901 using boost::parallel::all_reduce;
903 typename property_map<Graph, vertex_owner_t>::const_type
904 owner = get(vertex_owner, g);
906 // The graph must be undirected
908 (is_convertible<typename graph_traits<Graph>::directed_category,
909 undirected_tag>::value));
911 // The graph must model Incidence Graph
912 BOOST_CONCEPT_ASSERT(( IncidenceGraphConcept<Graph> ));
914 typedef typename graph_traits<Graph>::edges_size_type edges_size_type;
915 typedef typename graph_traits<Graph>::vertex_descriptor vertex_descriptor;
916 typedef typename graph_traits<Graph>::edge_descriptor edge_descriptor;
918 // Retrieve the process group we will use for communication
919 typedef typename process_group_type<Graph>::type process_group_type;
920 process_group_type pg = process_group(g);
922 // Keeps track of the edges that we know to be tree edges.
923 std::vector<edge_descriptor> tree_edges;
925 // The leaders send out a path message to initiate the algorithm
926 while (first != last) {
927 vertex_descriptor leader = *first;
928 if (process_id(pg) == get(owner, leader))
929 vertex_processor[leader].initialize_leader(leader, g);
934 // Will hold the number of bicomponents in the graph.
935 edges_size_type num_bicomponents = 0;
937 // Keep track of the path length that we should expect, based on the
938 // level in the breadth-first search tree. At present, this is only
939 // used as a sanity check. TBD: This could be used to decrease the
940 // amount of communication required per-edge (by about 4 bytes).
941 std::size_t path_length = 1;
943 typedef std::vector<vertex_descriptor> path_t;
945 unsigned char minimum_phase = 5;
947 while (optional<std::pair<int, int> > msg = probe(pg)) {
948 switch (msg->second) {
949 case msg_path_header:
951 // Receive the path header
952 path_header<edge_descriptor> header;
953 receive(pg, msg->first, msg->second, header);
954 BOOST_ASSERT(path_length == header.path_length);
956 // Receive the path itself
957 path_t path(path_length);
958 receive(pg, msg->first, msg_path_vertices, &path[0], path_length);
960 edge_descriptor e = header.edge;
961 vertex_processor[target(e, g)](e, path, g);
965 case msg_path_vertices:
966 // Should be handled in msg_path_header case, unless we're going
971 case msg_tree_header:
973 // Receive the tree header
974 tree_header<vertex_descriptor, edge_descriptor> header;
975 receive(pg, msg->first, msg->second, header);
977 // Receive the tree itself
978 path_t in_same_bicomponent(header.bicomp_length);
979 receive(pg, msg->first, msg_tree_vertices, &in_same_bicomponent[0],
980 header.bicomp_length);
982 edge_descriptor e = header.edge;
983 vertex_processor[target(e, g)](e, header.gamma, in_same_bicomponent,
988 case msg_tree_vertices:
989 // Should be handled in msg_tree_header case, unless we're
996 name_header<edge_descriptor> header;
997 receive(pg, msg->first, msg->second, header);
998 edge_descriptor e = header.edge;
999 vertex_processor[target(e, g)](e, header.name, g);
1004 BOOST_ASSERT(false);
1009 // Compute minimum phase locally
1011 unsigned char maximum_phase = 1;
1012 BGL_FORALL_VERTICES_T(v, g, Graph) {
1013 minimum_phase = (std::min)(minimum_phase, vertex_processor[v].get_phase());
1014 maximum_phase = (std::max)(maximum_phase, vertex_processor[v].get_phase());
1017 #ifdef PBGL_HOHBERG_DEBUG
1018 if (process_id(pg) == 0)
1019 std::cerr << "<---------End of stage------------->" << std::endl;
1021 // Compute minimum phase globally
1022 minimum_phase = all_reduce(pg, minimum_phase, boost::mpi::minimum<char>());
1024 #ifdef PBGL_HOHBERG_DEBUG
1025 if (process_id(pg) == 0)
1026 std::cerr << "Minimum phase = " << (int)minimum_phase << std::endl;
1029 if (minimum_phase == 4
1030 && all_reduce(pg, maximum_phase, boost::mpi::maximum<char>()) == 4) {
1032 #ifdef PBGL_HOHBERG_DEBUG
1033 if (process_id(pg) == 0)
1034 std::cerr << "<---------Naming phase------------->" << std::endl;
1036 // Compute the biconnected component number offsets for each
1038 std::vector<edges_size_type> local_offsets;
1039 local_offsets.reserve(num_vertices(g));
1040 edges_size_type num_local_bicomponents = 0;
1041 BGL_FORALL_VERTICES_T(v, g, Graph) {
1042 local_offsets.push_back(num_local_bicomponents);
1043 num_local_bicomponents +=
1044 vertex_processor[v].num_starting_bicomponents(v, g);
1049 // Find our the number of bicomponent names that will originate
1050 // from each process. This tells us how many bicomponents are in
1051 // the entire graph and what our global offset is for computing
1052 // our own biconnected component names.
1053 std::vector<edges_size_type> all_bicomponents(num_processes(pg));
1054 all_gather(pg, &num_local_bicomponents, &num_local_bicomponents + 1,
1056 num_bicomponents = 0;
1057 edges_size_type my_global_offset = 0;
1058 for (std::size_t i = 0; i < all_bicomponents.size(); ++i) {
1059 if (i == (std::size_t)process_id(pg))
1060 my_global_offset = num_bicomponents;
1061 num_bicomponents += all_bicomponents[i];
1064 std::size_t index = 0;
1065 BGL_FORALL_VERTICES_T(v, g, Graph) {
1066 edges_size_type offset = my_global_offset + local_offsets[index++];
1067 vertex_processor[v].start_naming_phase(v, g, offset);
1072 } while (minimum_phase < 5);
1074 // Number the edges appropriately.
1075 BGL_FORALL_VERTICES_T(v, g, Graph)
1076 vertex_processor[v].fill_edge_map(v, g, component);
1078 return num_bicomponents;
1081 template<typename Graph, typename ComponentMap, typename InputIterator>
1082 typename graph_traits<Graph>::edges_size_type
1083 hohberg_biconnected_components
1084 (const Graph& g, ComponentMap component,
1085 InputIterator first, InputIterator last)
1088 std::vector<hohberg_vertex_processor<Graph> >
1089 vertex_processors(num_vertices(g));
1090 return hohberg_biconnected_components
1091 (g, component, first, last,
1092 make_iterator_property_map(vertex_processors.begin(),
1093 get(vertex_index, g)));
1096 template<typename Graph, typename ComponentMap, typename ParentMap>
1097 typename graph_traits<Graph>::edges_size_type
1098 hohberg_biconnected_components(const Graph& g, ComponentMap component,
1101 // We need the connected components of the graph, but we don't care
1102 // about component numbers.
1103 connected_components(g, dummy_property_map(), parent);
1105 // Each root in the parent map is a leader
1106 typedef typename graph_traits<Graph>::vertex_descriptor vertex_descriptor;
1107 std::vector<vertex_descriptor> leaders;
1108 BGL_FORALL_VERTICES_T(v, g, Graph)
1109 if (get(parent, v) == v) leaders.push_back(v);
1111 return hohberg_biconnected_components(g, component,
1112 leaders.begin(), leaders.end());
1115 template<typename Graph, typename ComponentMap>
1116 typename graph_traits<Graph>::edges_size_type
1117 hohberg_biconnected_components(const Graph& g, ComponentMap component)
1119 typedef typename graph_traits<Graph>::vertex_descriptor vertex_descriptor;
1120 std::vector<vertex_descriptor> parents(num_vertices(g));
1121 return hohberg_biconnected_components
1122 (g, component, make_iterator_property_map(parents.begin(),
1123 get(vertex_index, g)));
1126 } } } // end namespace boost::graph::distributed
1128 #endif // BOOST_GRAPH_DISTRIBUTED_HOHBERG_BICONNECTED_COMPONENTS_HPP