• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 //=======================================================================
2 // Copyright (c) Aaron Windsor 2007
3 //
4 // Distributed under the Boost Software License, Version 1.0. (See
5 // accompanying file LICENSE_1_0.txt or copy at
6 // http://www.boost.org/LICENSE_1_0.txt)
7 //=======================================================================
8 #ifndef __BOYER_MYRVOLD_IMPL_HPP__
9 #define __BOYER_MYRVOLD_IMPL_HPP__
10 
11 #include <vector>
12 #include <list>
13 #include <boost/next_prior.hpp>
14 #include <boost/config.hpp> //for std::min macros
15 #include <boost/shared_ptr.hpp>
16 #include <boost/tuple/tuple.hpp>
17 #include <boost/property_map/property_map.hpp>
18 #include <boost/graph/graph_traits.hpp>
19 #include <boost/graph/depth_first_search.hpp>
20 #include <boost/graph/planar_detail/face_handles.hpp>
21 #include <boost/graph/planar_detail/face_iterators.hpp>
22 #include <boost/graph/planar_detail/bucket_sort.hpp>
23 
24 namespace boost
25 {
26 namespace detail
27 {
28     enum bm_case_t
29     {
30         BM_NO_CASE_CHOSEN,
31         BM_CASE_A,
32         BM_CASE_B,
33         BM_CASE_C,
34         BM_CASE_D,
35         BM_CASE_E
36     };
37 }
38 
39 template < typename LowPointMap, typename DFSParentMap, typename DFSNumberMap,
40     typename LeastAncestorMap, typename DFSParentEdgeMap, typename SizeType >
41 struct planar_dfs_visitor : public dfs_visitor<>
42 {
planar_dfs_visitorboost::planar_dfs_visitor43     planar_dfs_visitor(LowPointMap lpm, DFSParentMap dfs_p, DFSNumberMap dfs_n,
44         LeastAncestorMap lam, DFSParentEdgeMap dfs_edge)
45     : low(lpm)
46     , parent(dfs_p)
47     , df_number(dfs_n)
48     , least_ancestor(lam)
49     , df_edge(dfs_edge)
50     , count(0)
51     {
52     }
53 
54     template < typename Vertex, typename Graph >
start_vertexboost::planar_dfs_visitor55     void start_vertex(const Vertex& u, Graph&)
56     {
57         put(parent, u, u);
58         put(least_ancestor, u, count);
59     }
60 
61     template < typename Vertex, typename Graph >
discover_vertexboost::planar_dfs_visitor62     void discover_vertex(const Vertex& u, Graph&)
63     {
64         put(low, u, count);
65         put(df_number, u, count);
66         ++count;
67     }
68 
69     template < typename Edge, typename Graph >
tree_edgeboost::planar_dfs_visitor70     void tree_edge(const Edge& e, Graph& g)
71     {
72         typedef typename graph_traits< Graph >::vertex_descriptor vertex_t;
73         vertex_t s(source(e, g));
74         vertex_t t(target(e, g));
75 
76         put(parent, t, s);
77         put(df_edge, t, e);
78         put(least_ancestor, t, get(df_number, s));
79     }
80 
81     template < typename Edge, typename Graph >
back_edgeboost::planar_dfs_visitor82     void back_edge(const Edge& e, Graph& g)
83     {
84         typedef typename graph_traits< Graph >::vertex_descriptor vertex_t;
85         typedef typename graph_traits< Graph >::vertices_size_type v_size_t;
86 
87         vertex_t s(source(e, g));
88         vertex_t t(target(e, g));
89         BOOST_USING_STD_MIN();
90 
91         if (t != get(parent, s))
92         {
93             v_size_t s_low_df_number = get(low, s);
94             v_size_t t_df_number = get(df_number, t);
95             v_size_t s_least_ancestor_df_number = get(least_ancestor, s);
96 
97             put(low, s,
98                 min BOOST_PREVENT_MACRO_SUBSTITUTION(
99                     s_low_df_number, t_df_number));
100 
101             put(least_ancestor, s,
102                 min BOOST_PREVENT_MACRO_SUBSTITUTION(
103                     s_least_ancestor_df_number, t_df_number));
104         }
105     }
106 
107     template < typename Vertex, typename Graph >
finish_vertexboost::planar_dfs_visitor108     void finish_vertex(const Vertex& u, Graph&)
109     {
110         typedef typename graph_traits< Graph >::vertices_size_type v_size_t;
111 
112         Vertex u_parent = get(parent, u);
113         v_size_t u_parent_lowpoint = get(low, u_parent);
114         v_size_t u_lowpoint = get(low, u);
115         BOOST_USING_STD_MIN();
116 
117         if (u_parent != u)
118         {
119             put(low, u_parent,
120                 min BOOST_PREVENT_MACRO_SUBSTITUTION(
121                     u_lowpoint, u_parent_lowpoint));
122         }
123     }
124 
125     LowPointMap low;
126     DFSParentMap parent;
127     DFSNumberMap df_number;
128     LeastAncestorMap least_ancestor;
129     DFSParentEdgeMap df_edge;
130     SizeType count;
131 };
132 
133 template < typename Graph, typename VertexIndexMap,
134     typename StoreOldHandlesPolicy = graph::detail::store_old_handles,
135     typename StoreEmbeddingPolicy = graph::detail::recursive_lazy_list >
136 class boyer_myrvold_impl
137 {
138 
139     typedef typename graph_traits< Graph >::vertices_size_type v_size_t;
140     typedef typename graph_traits< Graph >::vertex_descriptor vertex_t;
141     typedef typename graph_traits< Graph >::edge_descriptor edge_t;
142     typedef typename graph_traits< Graph >::vertex_iterator vertex_iterator_t;
143     typedef typename graph_traits< Graph >::edge_iterator edge_iterator_t;
144     typedef
145         typename graph_traits< Graph >::out_edge_iterator out_edge_iterator_t;
146     typedef graph::detail::face_handle< Graph, StoreOldHandlesPolicy,
147         StoreEmbeddingPolicy >
148         face_handle_t;
149     typedef std::vector< vertex_t > vertex_vector_t;
150     typedef std::vector< edge_t > edge_vector_t;
151     typedef std::list< vertex_t > vertex_list_t;
152     typedef std::list< face_handle_t > face_handle_list_t;
153     typedef boost::shared_ptr< face_handle_list_t > face_handle_list_ptr_t;
154     typedef boost::shared_ptr< vertex_list_t > vertex_list_ptr_t;
155     typedef boost::tuple< vertex_t, bool, bool > merge_stack_frame_t;
156     typedef std::vector< merge_stack_frame_t > merge_stack_t;
157 
158     template < typename T > struct map_vertex_to_
159     {
160         typedef iterator_property_map< typename std::vector< T >::iterator,
161             VertexIndexMap >
162             type;
163     };
164 
165     typedef typename map_vertex_to_< v_size_t >::type vertex_to_v_size_map_t;
166     typedef typename map_vertex_to_< vertex_t >::type vertex_to_vertex_map_t;
167     typedef typename map_vertex_to_< edge_t >::type vertex_to_edge_map_t;
168     typedef typename map_vertex_to_< vertex_list_ptr_t >::type
169         vertex_to_vertex_list_ptr_map_t;
170     typedef typename map_vertex_to_< edge_vector_t >::type
171         vertex_to_edge_vector_map_t;
172     typedef typename map_vertex_to_< bool >::type vertex_to_bool_map_t;
173     typedef typename map_vertex_to_< face_handle_t >::type
174         vertex_to_face_handle_map_t;
175     typedef typename map_vertex_to_< face_handle_list_ptr_t >::type
176         vertex_to_face_handle_list_ptr_map_t;
177     typedef typename map_vertex_to_< typename vertex_list_t::iterator >::type
178         vertex_to_separated_node_map_t;
179 
180     template < typename BicompSideToTraverse = single_side,
181         typename VisitorType = lead_visitor, typename Time = current_iteration >
182     struct face_vertex_iterator
183     {
184         typedef face_iterator< Graph, vertex_to_face_handle_map_t, vertex_t,
185             BicompSideToTraverse, VisitorType, Time >
186             type;
187     };
188 
189     template < typename BicompSideToTraverse = single_side,
190         typename Time = current_iteration >
191     struct face_edge_iterator
192     {
193         typedef face_iterator< Graph, vertex_to_face_handle_map_t, edge_t,
194             BicompSideToTraverse, lead_visitor, Time >
195             type;
196     };
197 
198 public:
boyer_myrvold_impl(const Graph & arg_g,VertexIndexMap arg_vm)199     boyer_myrvold_impl(const Graph& arg_g, VertexIndexMap arg_vm)
200     : g(arg_g)
201     , vm(arg_vm)
202     ,
203 
204         low_point_vector(num_vertices(g))
205     , dfs_parent_vector(num_vertices(g))
206     , dfs_number_vector(num_vertices(g))
207     , least_ancestor_vector(num_vertices(g))
208     , pertinent_roots_vector(num_vertices(g))
209     , backedge_flag_vector(num_vertices(g), num_vertices(g) + 1)
210     , visited_vector(num_vertices(g), num_vertices(g) + 1)
211     , face_handles_vector(num_vertices(g))
212     , dfs_child_handles_vector(num_vertices(g))
213     , separated_dfs_child_list_vector(num_vertices(g))
214     , separated_node_in_parent_list_vector(num_vertices(g))
215     , canonical_dfs_child_vector(num_vertices(g))
216     , flipped_vector(num_vertices(g), false)
217     , backedges_vector(num_vertices(g))
218     , dfs_parent_edge_vector(num_vertices(g))
219     ,
220 
221         vertices_by_dfs_num(num_vertices(g))
222     ,
223 
224         low_point(low_point_vector.begin(), vm)
225     , dfs_parent(dfs_parent_vector.begin(), vm)
226     , dfs_number(dfs_number_vector.begin(), vm)
227     , least_ancestor(least_ancestor_vector.begin(), vm)
228     , pertinent_roots(pertinent_roots_vector.begin(), vm)
229     , backedge_flag(backedge_flag_vector.begin(), vm)
230     , visited(visited_vector.begin(), vm)
231     , face_handles(face_handles_vector.begin(), vm)
232     , dfs_child_handles(dfs_child_handles_vector.begin(), vm)
233     , separated_dfs_child_list(separated_dfs_child_list_vector.begin(), vm)
234     , separated_node_in_parent_list(
235           separated_node_in_parent_list_vector.begin(), vm)
236     , canonical_dfs_child(canonical_dfs_child_vector.begin(), vm)
237     , flipped(flipped_vector.begin(), vm)
238     , backedges(backedges_vector.begin(), vm)
239     , dfs_parent_edge(dfs_parent_edge_vector.begin(), vm)
240 
241     {
242 
243         planar_dfs_visitor< vertex_to_v_size_map_t, vertex_to_vertex_map_t,
244             vertex_to_v_size_map_t, vertex_to_v_size_map_t,
245             vertex_to_edge_map_t, v_size_t >
246             vis(low_point, dfs_parent, dfs_number, least_ancestor,
247                 dfs_parent_edge);
248 
249         // Perform a depth-first search to find each vertex's low point, least
250         // ancestor, and dfs tree information
251         depth_first_search(g, visitor(vis).vertex_index_map(vm));
252 
253         // Sort vertices by their lowpoint - need this later in the constructor
254         vertex_vector_t vertices_by_lowpoint(num_vertices(g));
255         std::copy(vertices(g).first, vertices(g).second,
256             vertices_by_lowpoint.begin());
257         bucket_sort(vertices_by_lowpoint.begin(), vertices_by_lowpoint.end(),
258             low_point, num_vertices(g));
259 
260         // Sort vertices by their dfs number - need this to iterate by reverse
261         // DFS number in the main loop.
262         std::copy(
263             vertices(g).first, vertices(g).second, vertices_by_dfs_num.begin());
264         bucket_sort(vertices_by_dfs_num.begin(), vertices_by_dfs_num.end(),
265             dfs_number, num_vertices(g));
266 
267         // Initialize face handles. A face handle is an abstraction that serves
268         // two uses in our implementation - it allows us to efficiently move
269         // along the outer face of embedded bicomps in a partially embedded
270         // graph, and it provides storage for the planar embedding. Face
271         // handles are implemented by a sequence of edges and are associated
272         // with a particular vertex - the sequence of edges represents the
273         // current embedding of edges around that vertex, and the first and
274         // last edges in the sequence represent the pair of edges on the outer
275         // face that are adjacent to the associated vertex. This lets us embed
276         // edges in the graph by just pushing them on the front or back of the
277         // sequence of edges held by the face handles.
278         //
279         // Our algorithm starts with a DFS tree of edges (where every vertex is
280         // an articulation point and every edge is a singleton bicomp) and
281         // repeatedly merges bicomps by embedding additional edges. Note that
282         // any bicomp at any point in the algorithm can be associated with a
283         // unique edge connecting the vertex of that bicomp with the lowest DFS
284         // number (which we refer to as the "root" of the bicomp) with its DFS
285         // child in the bicomp: the existence of two such edges would contradict
286         // the properties of a DFS tree. We refer to the DFS child of the root
287         // of a bicomp as the "canonical DFS child" of the bicomp. Note that a
288         // vertex can be the root of more than one bicomp.
289         //
290         // We move around the external faces of a bicomp using a few property
291         // maps, which we'll initialize presently:
292         //
293         // - face_handles: maps a vertex to a face handle that can be used to
294         //   move "up" a bicomp. For a vertex that isn't an articulation point,
295         //   this holds the face handles that can be used to move around that
296         //   vertex's unique bicomp. For a vertex that is an articulation point,
297         //   this holds the face handles associated with the unique bicomp that
298         //   the vertex is NOT the root of. These handles can therefore be used
299         //   to move from any point on the outer face of the tree of bicomps
300         //   around the current outer face towards the root of the DFS tree.
301         //
302         // - dfs_child_handles: these are used to hold face handles for
303         //   vertices that are articulation points - dfs_child_handles[v] holds
304         //   the face handles corresponding to vertex u in the bicomp with root
305         //   u and canonical DFS child v.
306         //
307         // - canonical_dfs_child: this property map allows one to determine the
308         //   canonical DFS child of a bicomp while traversing the outer face.
309         //   This property map is only valid when applied to one of the two
310         //   vertices adjacent to the root of the bicomp on the outer face. To
311         //   be more precise, if v is the canonical DFS child of a bicomp,
312         //   canonical_dfs_child[dfs_child_handles[v].first_vertex()] == v and
313         //   canonical_dfs_child[dfs_child_handles[v].second_vertex()] == v.
314         //
315         // - pertinent_roots: given a vertex v, pertinent_roots[v] contains a
316         //   list of face handles pointing to the top of bicomps that need to
317         //   be visited by the current walkdown traversal (since they lead to
318         //   backedges that need to be embedded). These lists are populated by
319         //   the walkup and consumed by the walkdown.
320 
321         vertex_iterator_t vi, vi_end;
322         for (boost::tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi)
323         {
324             vertex_t v(*vi);
325             vertex_t parent = dfs_parent[v];
326 
327             if (parent != v)
328             {
329                 edge_t parent_edge = dfs_parent_edge[v];
330                 add_to_embedded_edges(parent_edge, StoreOldHandlesPolicy());
331                 face_handles[v] = face_handle_t(v, parent_edge, g);
332                 dfs_child_handles[v] = face_handle_t(parent, parent_edge, g);
333             }
334             else
335             {
336                 face_handles[v] = face_handle_t(v);
337                 dfs_child_handles[v] = face_handle_t(parent);
338             }
339 
340             canonical_dfs_child[v] = v;
341             pertinent_roots[v] = face_handle_list_ptr_t(new face_handle_list_t);
342             separated_dfs_child_list[v] = vertex_list_ptr_t(new vertex_list_t);
343         }
344 
345         // We need to create a list of not-yet-merged depth-first children for
346         // each vertex that will be updated as bicomps get merged. We sort each
347         // list by ascending lowpoint, which allows the externally_active
348         // function to run in constant time, and we keep a pointer to each
349         // vertex's representation in its parent's list, which allows merging
350         // in constant time.
351 
352         for (typename vertex_vector_t::iterator itr
353              = vertices_by_lowpoint.begin();
354              itr != vertices_by_lowpoint.end(); ++itr)
355         {
356             vertex_t v(*itr);
357             vertex_t parent(dfs_parent[v]);
358             if (v != parent)
359             {
360                 separated_node_in_parent_list[v]
361                     = separated_dfs_child_list[parent]->insert(
362                         separated_dfs_child_list[parent]->end(), v);
363             }
364         }
365 
366         // The merge stack holds path information during a walkdown iteration
367         merge_stack.reserve(num_vertices(g));
368     }
369 
is_planar()370     bool is_planar()
371     {
372 
373         // This is the main algorithm: starting with a DFS tree of embedded
374         // edges (which, since it's a tree, is planar), iterate through all
375         // vertices by reverse DFS number, attempting to embed all backedges
376         // connecting the current vertex to vertices with higher DFS numbers.
377         //
378         // The walkup is a procedure that examines all such backedges and sets
379         // up the required data structures so that they can be searched by the
380         // walkdown in linear time. The walkdown does the actual work of
381         // embedding edges and flipping bicomps, and can identify when it has
382         // come across a kuratowski subgraph.
383         //
384         // store_old_face_handles caches face handles from the previous
385         // iteration - this is used only for the kuratowski subgraph isolation,
386         // and is therefore dispatched based on the StoreOldHandlesPolicy.
387         //
388         // clean_up_embedding does some clean-up and fills in values that have
389         // to be computed lazily during the actual execution of the algorithm
390         // (for instance, whether or not a bicomp is flipped in the final
391         // embedding). It's dispatched on the the StoreEmbeddingPolicy, since
392         // it's not needed if an embedding isn't desired.
393 
394         typename vertex_vector_t::reverse_iterator vi, vi_end;
395 
396         vi_end = vertices_by_dfs_num.rend();
397         for (vi = vertices_by_dfs_num.rbegin(); vi != vi_end; ++vi)
398         {
399 
400             store_old_face_handles(StoreOldHandlesPolicy());
401 
402             vertex_t v(*vi);
403 
404             walkup(v);
405 
406             if (!walkdown(v))
407                 return false;
408         }
409 
410         clean_up_embedding(StoreEmbeddingPolicy());
411 
412         return true;
413     }
414 
415 private:
walkup(vertex_t v)416     void walkup(vertex_t v)
417     {
418 
419         // The point of the walkup is to follow all backedges from v to
420         // vertices with higher DFS numbers, and update pertinent_roots
421         // for the bicomp roots on the path from backedge endpoints up
422         // to v. This will set the stage for the walkdown to efficiently
423         // traverse the graph of bicomps down from v.
424 
425         typedef
426             typename face_vertex_iterator< both_sides >::type walkup_iterator_t;
427 
428         out_edge_iterator_t oi, oi_end;
429         for (boost::tie(oi, oi_end) = out_edges(v, g); oi != oi_end; ++oi)
430         {
431             edge_t e(*oi);
432             vertex_t e_source(source(e, g));
433             vertex_t e_target(target(e, g));
434 
435             if (e_source == e_target)
436             {
437                 self_loops.push_back(e);
438                 continue;
439             }
440 
441             vertex_t w(e_source == v ? e_target : e_source);
442 
443             // continue if not a back edge or already embedded
444             if (dfs_number[w] < dfs_number[v] || e == dfs_parent_edge[w])
445                 continue;
446 
447             backedges[w].push_back(e);
448 
449             v_size_t timestamp = dfs_number[v];
450             backedge_flag[w] = timestamp;
451 
452             walkup_iterator_t walkup_itr(w, face_handles);
453             walkup_iterator_t walkup_end;
454             vertex_t lead_vertex = w;
455 
456             while (true)
457             {
458 
459                 // Move to the root of the current bicomp or the first visited
460                 // vertex on the bicomp by going up each side in parallel
461 
462                 while (walkup_itr != walkup_end
463                     && visited[*walkup_itr] != timestamp)
464                 {
465                     lead_vertex = *walkup_itr;
466                     visited[lead_vertex] = timestamp;
467                     ++walkup_itr;
468                 }
469 
470                 // If we've found the root of a bicomp through a path we haven't
471                 // seen before, update pertinent_roots with a handle to the
472                 // current bicomp. Otherwise, we've just seen a path we've been
473                 // up before, so break out of the main while loop.
474 
475                 if (walkup_itr == walkup_end)
476                 {
477                     vertex_t dfs_child = canonical_dfs_child[lead_vertex];
478                     vertex_t parent = dfs_parent[dfs_child];
479 
480                     visited[dfs_child_handles[dfs_child].first_vertex()]
481                         = timestamp;
482                     visited[dfs_child_handles[dfs_child].second_vertex()]
483                         = timestamp;
484 
485                     if (low_point[dfs_child] < dfs_number[v]
486                         || least_ancestor[dfs_child] < dfs_number[v])
487                     {
488                         pertinent_roots[parent]->push_back(
489                             dfs_child_handles[dfs_child]);
490                     }
491                     else
492                     {
493                         pertinent_roots[parent]->push_front(
494                             dfs_child_handles[dfs_child]);
495                     }
496 
497                     if (parent != v && visited[parent] != timestamp)
498                     {
499                         walkup_itr = walkup_iterator_t(parent, face_handles);
500                         lead_vertex = parent;
501                     }
502                     else
503                         break;
504                 }
505                 else
506                     break;
507             }
508         }
509     }
510 
walkdown(vertex_t v)511     bool walkdown(vertex_t v)
512     {
513         // This procedure is where all of the action is - pertinent_roots
514         // has already been set up by the walkup, so we just need to move
515         // down bicomps from v until we find vertices that have been
516         // labeled as backedge endpoints. Once we find such a vertex, we
517         // embed the corresponding edge and glue together the bicomps on
518         // the path connecting the two vertices in the edge. This may
519         // involve flipping bicomps along the way.
520 
521         vertex_t w; // the other endpoint of the edge we're embedding
522 
523         while (!pertinent_roots[v]->empty())
524         {
525 
526             face_handle_t root_face_handle = pertinent_roots[v]->front();
527             face_handle_t curr_face_handle = root_face_handle;
528             pertinent_roots[v]->pop_front();
529 
530             merge_stack.clear();
531 
532             while (true)
533             {
534 
535                 typename face_vertex_iterator<>::type first_face_itr,
536                     second_face_itr, face_end;
537                 vertex_t first_side_vertex
538                     = graph_traits< Graph >::null_vertex();
539                 vertex_t second_side_vertex
540                     = graph_traits< Graph >::null_vertex();
541                 vertex_t first_tail, second_tail;
542 
543                 first_tail = second_tail = curr_face_handle.get_anchor();
544                 first_face_itr = typename face_vertex_iterator<>::type(
545                     curr_face_handle, face_handles, first_side());
546                 second_face_itr = typename face_vertex_iterator<>::type(
547                     curr_face_handle, face_handles, second_side());
548 
549                 for (; first_face_itr != face_end; ++first_face_itr)
550                 {
551                     vertex_t face_vertex(*first_face_itr);
552                     if (pertinent(face_vertex, v)
553                         || externally_active(face_vertex, v))
554                     {
555                         first_side_vertex = face_vertex;
556                         second_side_vertex = face_vertex;
557                         break;
558                     }
559                     first_tail = face_vertex;
560                 }
561 
562                 if (first_side_vertex == graph_traits< Graph >::null_vertex()
563                     || first_side_vertex == curr_face_handle.get_anchor())
564                     break;
565 
566                 for (; second_face_itr != face_end; ++second_face_itr)
567                 {
568                     vertex_t face_vertex(*second_face_itr);
569                     if (pertinent(face_vertex, v)
570                         || externally_active(face_vertex, v))
571                     {
572                         second_side_vertex = face_vertex;
573                         break;
574                     }
575                     second_tail = face_vertex;
576                 }
577 
578                 vertex_t chosen;
579                 bool chose_first_upper_path;
580                 if (internally_active(first_side_vertex, v))
581                 {
582                     chosen = first_side_vertex;
583                     chose_first_upper_path = true;
584                 }
585                 else if (internally_active(second_side_vertex, v))
586                 {
587                     chosen = second_side_vertex;
588                     chose_first_upper_path = false;
589                 }
590                 else if (pertinent(first_side_vertex, v))
591                 {
592                     chosen = first_side_vertex;
593                     chose_first_upper_path = true;
594                 }
595                 else if (pertinent(second_side_vertex, v))
596                 {
597                     chosen = second_side_vertex;
598                     chose_first_upper_path = false;
599                 }
600                 else
601                 {
602 
603                     // If there's a pertinent vertex on the lower face
604                     // between the first_face_itr and the second_face_itr,
605                     // this graph isn't planar.
606                     for (; *first_face_itr != second_side_vertex;
607                          ++first_face_itr)
608                     {
609                         vertex_t p(*first_face_itr);
610                         if (pertinent(p, v))
611                         {
612                             // Found a Kuratowski subgraph
613                             kuratowski_v = v;
614                             kuratowski_x = first_side_vertex;
615                             kuratowski_y = second_side_vertex;
616                             return false;
617                         }
618                     }
619 
620                     // Otherwise, the fact that we didn't find a pertinent
621                     // vertex on this face is fine - we should set the
622                     // short-circuit edges and break out of this loop to
623                     // start looking at a different pertinent root.
624 
625                     if (first_side_vertex == second_side_vertex)
626                     {
627                         if (first_tail != v)
628                         {
629                             vertex_t first
630                                 = face_handles[first_tail].first_vertex();
631                             vertex_t second
632                                 = face_handles[first_tail].second_vertex();
633                             boost::tie(first_side_vertex, first_tail)
634                                 = make_tuple(first_tail,
635                                     first == first_side_vertex ? second
636                                                                : first);
637                         }
638                         else if (second_tail != v)
639                         {
640                             vertex_t first
641                                 = face_handles[second_tail].first_vertex();
642                             vertex_t second
643                                 = face_handles[second_tail].second_vertex();
644                             boost::tie(second_side_vertex, second_tail)
645                                 = make_tuple(second_tail,
646                                     first == second_side_vertex ? second
647                                                                 : first);
648                         }
649                         else
650                             break;
651                     }
652 
653                     canonical_dfs_child[first_side_vertex]
654                         = canonical_dfs_child[root_face_handle.first_vertex()];
655                     canonical_dfs_child[second_side_vertex]
656                         = canonical_dfs_child[root_face_handle.second_vertex()];
657                     root_face_handle.set_first_vertex(first_side_vertex);
658                     root_face_handle.set_second_vertex(second_side_vertex);
659 
660                     if (face_handles[first_side_vertex].first_vertex()
661                         == first_tail)
662                         face_handles[first_side_vertex].set_first_vertex(v);
663                     else
664                         face_handles[first_side_vertex].set_second_vertex(v);
665 
666                     if (face_handles[second_side_vertex].first_vertex()
667                         == second_tail)
668                         face_handles[second_side_vertex].set_first_vertex(v);
669                     else
670                         face_handles[second_side_vertex].set_second_vertex(v);
671 
672                     break;
673                 }
674 
675                 // When we unwind the stack, we need to know which direction
676                 // we came down from on the top face handle
677 
678                 bool chose_first_lower_path
679                     = (chose_first_upper_path
680                           && face_handles[chosen].first_vertex() == first_tail)
681                     || (!chose_first_upper_path
682                         && face_handles[chosen].first_vertex() == second_tail);
683 
684                 // If there's a backedge at the chosen vertex, embed it now
685                 if (backedge_flag[chosen] == dfs_number[v])
686                 {
687                     w = chosen;
688 
689                     backedge_flag[chosen] = num_vertices(g) + 1;
690                     add_to_merge_points(chosen, StoreOldHandlesPolicy());
691 
692                     typename edge_vector_t::iterator ei, ei_end;
693                     ei_end = backedges[chosen].end();
694                     for (ei = backedges[chosen].begin(); ei != ei_end; ++ei)
695                     {
696                         edge_t e(*ei);
697                         add_to_embedded_edges(e, StoreOldHandlesPolicy());
698 
699                         if (chose_first_lower_path)
700                             face_handles[chosen].push_first(e, g);
701                         else
702                             face_handles[chosen].push_second(e, g);
703                     }
704                 }
705                 else
706                 {
707                     merge_stack.push_back(make_tuple(chosen,
708                         chose_first_upper_path, chose_first_lower_path));
709                     curr_face_handle = *pertinent_roots[chosen]->begin();
710                     continue;
711                 }
712 
713                 // Unwind the merge stack to the root, merging all bicomps
714 
715                 bool bottom_path_follows_first;
716                 bool top_path_follows_first;
717                 bool next_bottom_follows_first = chose_first_upper_path;
718 
719                 vertex_t merge_point = chosen;
720 
721                 while (!merge_stack.empty())
722                 {
723 
724                     bottom_path_follows_first = next_bottom_follows_first;
725                     boost::tie(merge_point, next_bottom_follows_first,
726                         top_path_follows_first)
727                         = merge_stack.back();
728                     merge_stack.pop_back();
729 
730                     face_handle_t top_handle(face_handles[merge_point]);
731                     face_handle_t bottom_handle(
732                         *pertinent_roots[merge_point]->begin());
733 
734                     vertex_t bottom_dfs_child = canonical_dfs_child
735                         [pertinent_roots[merge_point]->begin()->first_vertex()];
736 
737                     remove_vertex_from_separated_dfs_child_list(
738                         canonical_dfs_child[pertinent_roots[merge_point]
739                                                 ->begin()
740                                                 ->first_vertex()]);
741 
742                     pertinent_roots[merge_point]->pop_front();
743 
744                     add_to_merge_points(
745                         top_handle.get_anchor(), StoreOldHandlesPolicy());
746 
747                     if (top_path_follows_first && bottom_path_follows_first)
748                     {
749                         bottom_handle.flip();
750                         top_handle.glue_first_to_second(bottom_handle);
751                     }
752                     else if (!top_path_follows_first
753                         && bottom_path_follows_first)
754                     {
755                         flipped[bottom_dfs_child] = true;
756                         top_handle.glue_second_to_first(bottom_handle);
757                     }
758                     else if (top_path_follows_first
759                         && !bottom_path_follows_first)
760                     {
761                         flipped[bottom_dfs_child] = true;
762                         top_handle.glue_first_to_second(bottom_handle);
763                     }
764                     else //! top_path_follows_first &&
765                          //! !bottom_path_follows_first
766                     {
767                         bottom_handle.flip();
768                         top_handle.glue_second_to_first(bottom_handle);
769                     }
770                 }
771 
772                 // Finally, embed all edges (v,w) at their upper end points
773                 canonical_dfs_child[w]
774                     = canonical_dfs_child[root_face_handle.first_vertex()];
775 
776                 add_to_merge_points(
777                     root_face_handle.get_anchor(), StoreOldHandlesPolicy());
778 
779                 typename edge_vector_t::iterator ei, ei_end;
780                 ei_end = backedges[chosen].end();
781                 for (ei = backedges[chosen].begin(); ei != ei_end; ++ei)
782                 {
783                     if (next_bottom_follows_first)
784                         root_face_handle.push_first(*ei, g);
785                     else
786                         root_face_handle.push_second(*ei, g);
787                 }
788 
789                 backedges[chosen].clear();
790                 curr_face_handle = root_face_handle;
791 
792             } // while(true)
793 
794         } // while(!pertinent_roots[v]->empty())
795 
796         return true;
797     }
798 
store_old_face_handles(graph::detail::no_old_handles)799     void store_old_face_handles(graph::detail::no_old_handles) {}
800 
store_old_face_handles(graph::detail::store_old_handles)801     void store_old_face_handles(graph::detail::store_old_handles)
802     {
803         for (typename std::vector< vertex_t >::iterator mp_itr
804              = current_merge_points.begin();
805              mp_itr != current_merge_points.end(); ++mp_itr)
806         {
807             face_handles[*mp_itr].store_old_face_handles();
808         }
809         current_merge_points.clear();
810     }
811 
add_to_merge_points(vertex_t,graph::detail::no_old_handles)812     void add_to_merge_points(vertex_t, graph::detail::no_old_handles) {}
813 
add_to_merge_points(vertex_t v,graph::detail::store_old_handles)814     void add_to_merge_points(vertex_t v, graph::detail::store_old_handles)
815     {
816         current_merge_points.push_back(v);
817     }
818 
add_to_embedded_edges(edge_t,graph::detail::no_old_handles)819     void add_to_embedded_edges(edge_t, graph::detail::no_old_handles) {}
820 
add_to_embedded_edges(edge_t e,graph::detail::store_old_handles)821     void add_to_embedded_edges(edge_t e, graph::detail::store_old_handles)
822     {
823         embedded_edges.push_back(e);
824     }
825 
clean_up_embedding(graph::detail::no_embedding)826     void clean_up_embedding(graph::detail::no_embedding) {}
827 
clean_up_embedding(graph::detail::store_embedding)828     void clean_up_embedding(graph::detail::store_embedding)
829     {
830 
831         // If the graph isn't biconnected, we'll still have entries
832         // in the separated_dfs_child_list for some vertices. Since
833         // these represent articulation points, we can obtain a
834         // planar embedding no matter what order we embed them in.
835 
836         vertex_iterator_t xi, xi_end;
837         for (boost::tie(xi, xi_end) = vertices(g); xi != xi_end; ++xi)
838         {
839             if (!separated_dfs_child_list[*xi]->empty())
840             {
841                 typename vertex_list_t::iterator yi, yi_end;
842                 yi_end = separated_dfs_child_list[*xi]->end();
843                 for (yi = separated_dfs_child_list[*xi]->begin(); yi != yi_end;
844                      ++yi)
845                 {
846                     dfs_child_handles[*yi].flip();
847                     face_handles[*xi].glue_first_to_second(
848                         dfs_child_handles[*yi]);
849                 }
850             }
851         }
852 
853         // Up until this point, we've flipped bicomps lazily by setting
854         // flipped[v] to true if the bicomp rooted at v was flipped (the
855         // lazy aspect of this flip is that all descendents of that vertex
856         // need to have their orientations reversed as well). Now, we
857         // traverse the DFS tree by DFS number and perform the actual
858         // flipping as needed
859 
860         typedef typename vertex_vector_t::iterator vertex_vector_itr_t;
861         vertex_vector_itr_t vi_end = vertices_by_dfs_num.end();
862         for (vertex_vector_itr_t vi = vertices_by_dfs_num.begin(); vi != vi_end;
863              ++vi)
864         {
865             vertex_t v(*vi);
866             bool v_flipped = flipped[v];
867             bool p_flipped = flipped[dfs_parent[v]];
868             if (v_flipped && !p_flipped)
869             {
870                 face_handles[v].flip();
871             }
872             else if (p_flipped && !v_flipped)
873             {
874                 face_handles[v].flip();
875                 flipped[v] = true;
876             }
877             else
878             {
879                 flipped[v] = false;
880             }
881         }
882 
883         // If there are any self-loops in the graph, they were flagged
884         // during the walkup, and we should add them to the embedding now.
885         // Adding a self loop anywhere in the embedding could never
886         // invalidate the embedding, but they would complicate the traversal
887         // if they were added during the walkup/walkdown.
888 
889         typename edge_vector_t::iterator ei, ei_end;
890         ei_end = self_loops.end();
891         for (ei = self_loops.begin(); ei != ei_end; ++ei)
892         {
893             edge_t e(*ei);
894             face_handles[source(e, g)].push_second(e, g);
895         }
896     }
897 
pertinent(vertex_t w,vertex_t v)898     bool pertinent(vertex_t w, vertex_t v)
899     {
900         // w is pertinent with respect to v if there is a backedge (v,w) or if
901         // w is the root of a bicomp that contains a pertinent vertex.
902 
903         return backedge_flag[w] == dfs_number[v]
904             || !pertinent_roots[w]->empty();
905     }
906 
externally_active(vertex_t w,vertex_t v)907     bool externally_active(vertex_t w, vertex_t v)
908     {
909         // Let a be any proper depth-first search ancestor of v. w is externally
910         // active with respect to v if there exists a backedge (a,w) or a
911         // backedge (a,w_0) for some w_0 in a descendent bicomp of w.
912 
913         v_size_t dfs_number_of_v = dfs_number[v];
914         return (least_ancestor[w] < dfs_number_of_v)
915             || (!separated_dfs_child_list[w]->empty()
916                 && low_point[separated_dfs_child_list[w]->front()]
917                     < dfs_number_of_v);
918     }
919 
internally_active(vertex_t w,vertex_t v)920     bool internally_active(vertex_t w, vertex_t v)
921     {
922         return pertinent(w, v) && !externally_active(w, v);
923     }
924 
remove_vertex_from_separated_dfs_child_list(vertex_t v)925     void remove_vertex_from_separated_dfs_child_list(vertex_t v)
926     {
927         typename vertex_list_t::iterator to_delete
928             = separated_node_in_parent_list[v];
929         garbage.splice(garbage.end(), *separated_dfs_child_list[dfs_parent[v]],
930             to_delete, boost::next(to_delete));
931     }
932 
933     // End of the implementation of the basic Boyer-Myrvold Algorithm. The rest
934     // of the code below implements the isolation of a Kuratowski subgraph in
935     // the case that the input graph is not planar. This is by far the most
936     // complicated part of the implementation.
937 
938 public:
939     template < typename EdgeToBoolPropertyMap, typename EdgeContainer >
kuratowski_walkup(vertex_t v,EdgeToBoolPropertyMap forbidden_edge,EdgeToBoolPropertyMap goal_edge,EdgeToBoolPropertyMap is_embedded,EdgeContainer & path_edges)940     vertex_t kuratowski_walkup(vertex_t v, EdgeToBoolPropertyMap forbidden_edge,
941         EdgeToBoolPropertyMap goal_edge, EdgeToBoolPropertyMap is_embedded,
942         EdgeContainer& path_edges)
943     {
944         vertex_t current_endpoint;
945         bool seen_goal_edge = false;
946         out_edge_iterator_t oi, oi_end;
947 
948         for (boost::tie(oi, oi_end) = out_edges(v, g); oi != oi_end; ++oi)
949             forbidden_edge[*oi] = true;
950 
951         for (boost::tie(oi, oi_end) = out_edges(v, g); oi != oi_end; ++oi)
952         {
953             path_edges.clear();
954 
955             edge_t e(*oi);
956             current_endpoint
957                 = target(*oi, g) == v ? source(*oi, g) : target(*oi, g);
958 
959             if (dfs_number[current_endpoint] < dfs_number[v] || is_embedded[e]
960                 || v == current_endpoint // self-loop
961             )
962             {
963                 // Not a backedge
964                 continue;
965             }
966 
967             path_edges.push_back(e);
968             if (goal_edge[e])
969             {
970                 return current_endpoint;
971             }
972 
973             typedef typename face_edge_iterator<>::type walkup_itr_t;
974 
975             walkup_itr_t walkup_itr(
976                 current_endpoint, face_handles, first_side());
977             walkup_itr_t walkup_end;
978 
979             seen_goal_edge = false;
980 
981             while (true)
982             {
983 
984                 if (walkup_itr != walkup_end && forbidden_edge[*walkup_itr])
985                     break;
986 
987                 while (walkup_itr != walkup_end && !goal_edge[*walkup_itr]
988                     && !forbidden_edge[*walkup_itr])
989                 {
990                     edge_t f(*walkup_itr);
991                     forbidden_edge[f] = true;
992                     path_edges.push_back(f);
993                     current_endpoint = source(f, g) == current_endpoint
994                         ? target(f, g)
995                         : source(f, g);
996                     ++walkup_itr;
997                 }
998 
999                 if (walkup_itr != walkup_end && goal_edge[*walkup_itr])
1000                 {
1001                     path_edges.push_back(*walkup_itr);
1002                     seen_goal_edge = true;
1003                     break;
1004                 }
1005 
1006                 walkup_itr = walkup_itr_t(
1007                     current_endpoint, face_handles, first_side());
1008             }
1009 
1010             if (seen_goal_edge)
1011                 break;
1012         }
1013 
1014         if (seen_goal_edge)
1015             return current_endpoint;
1016         else
1017             return graph_traits< Graph >::null_vertex();
1018     }
1019 
1020     template < typename OutputIterator, typename EdgeIndexMap >
extract_kuratowski_subgraph(OutputIterator o_itr,EdgeIndexMap em)1021     void extract_kuratowski_subgraph(OutputIterator o_itr, EdgeIndexMap em)
1022     {
1023 
1024         // If the main algorithm has failed to embed one of the back-edges from
1025         // a vertex v, we can use the current state of the algorithm to isolate
1026         // a Kuratowksi subgraph. The isolation process breaks down into five
1027         // cases, A - E. The general configuration of all five cases is shown in
1028         //                  figure 1. There is a vertex v from which the planar
1029         //         v        embedding process could not proceed. This means that
1030         //         |        there exists some bicomp containing three vertices
1031         //       -----      x,y, and z as shown such that x and y are externally
1032         //      |     |     active with respect to v (which means that there are
1033         //      x     y     two vertices x_0 and y_0 such that (1) both x_0 and
1034         //      |     |     y_0 are proper depth-first search ancestors of v and
1035         //       --z--      (2) there are two disjoint paths, one connecting x
1036         //                  and x_0 and one connecting y and y_0, both
1037         //                  consisting
1038         //       fig. 1     entirely of unembedded edges). Furthermore, there
1039         //                  exists a vertex z_0 such that z is a depth-first
1040         // search ancestor of z_0 and (v,z_0) is an unembedded back-edge from v.
1041         // x,y and z all exist on the same bicomp, which consists entirely of
1042         // embedded edges. The five subcases break down as follows, and are
1043         // handled by the algorithm logically in the order A-E: First, if v is
1044         // not on the same bicomp as x,y, and z, a K_3_3 can be isolated - this
1045         // is case A. So, we'll assume that v is on the same bicomp as x,y, and
1046         // z. If z_0 is on a different bicomp than x,y, and z, a K_3_3 can also
1047         // be isolated - this is a case B - so we'll assume from now on that v
1048         // is on the same bicomp as x, y, and z=z_0. In this case, one can use
1049         // properties of the Boyer-Myrvold algorithm to show the existence of an
1050         // "x-y path" connecting some vertex on the "left side" of the x,y,z
1051         // bicomp with some vertex on the "right side" of the bicomp (where the
1052         // left and right are split by a line drawn through v and z.If either of
1053         // the endpoints of the x-y path is above x or y on the bicomp, a K_3_3
1054         // can be isolated - this is a case C. Otherwise, both endpoints are at
1055         // or below x and y on the bicomp. If there is a vertex alpha on the x-y
1056         // path such that alpha is not x or y and there's a path from alpha to v
1057         // that's disjoint from any of the edges on the bicomp and the x-y path,
1058         // a K_3_3 can be isolated - this is a case D. Otherwise, properties of
1059         // the Boyer-Myrvold algorithm can be used to show that another vertex
1060         // w exists on the lower half of the bicomp such that w is externally
1061         // active with respect to v. w can then be used to isolate a K_5 - this
1062         // is the configuration of case E.
1063 
1064         vertex_iterator_t vi, vi_end;
1065         edge_iterator_t ei, ei_end;
1066         out_edge_iterator_t oei, oei_end;
1067         typename std::vector< edge_t >::iterator xi, xi_end;
1068 
1069         // Clear the short-circuit edges - these are needed for the planar
1070         // testing/embedding algorithm to run in linear time, but they'll
1071         // complicate the kuratowski subgraph isolation
1072         for (boost::tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi)
1073         {
1074             face_handles[*vi].reset_vertex_cache();
1075             dfs_child_handles[*vi].reset_vertex_cache();
1076         }
1077 
1078         vertex_t v = kuratowski_v;
1079         vertex_t x = kuratowski_x;
1080         vertex_t y = kuratowski_y;
1081 
1082         typedef iterator_property_map< typename std::vector< bool >::iterator,
1083             EdgeIndexMap >
1084             edge_to_bool_map_t;
1085 
1086         std::vector< bool > is_in_subgraph_vector(num_edges(g), false);
1087         edge_to_bool_map_t is_in_subgraph(is_in_subgraph_vector.begin(), em);
1088 
1089         std::vector< bool > is_embedded_vector(num_edges(g), false);
1090         edge_to_bool_map_t is_embedded(is_embedded_vector.begin(), em);
1091 
1092         typename std::vector< edge_t >::iterator embedded_itr, embedded_end;
1093         embedded_end = embedded_edges.end();
1094         for (embedded_itr = embedded_edges.begin();
1095              embedded_itr != embedded_end; ++embedded_itr)
1096             is_embedded[*embedded_itr] = true;
1097 
1098         // upper_face_vertex is true for x,y, and all vertices above x and y in
1099         // the bicomp
1100         std::vector< bool > upper_face_vertex_vector(num_vertices(g), false);
1101         vertex_to_bool_map_t upper_face_vertex(
1102             upper_face_vertex_vector.begin(), vm);
1103 
1104         std::vector< bool > lower_face_vertex_vector(num_vertices(g), false);
1105         vertex_to_bool_map_t lower_face_vertex(
1106             lower_face_vertex_vector.begin(), vm);
1107 
1108         // These next few variable declarations are all things that we need
1109         // to find.
1110         vertex_t z = graph_traits< Graph >::null_vertex();
1111         vertex_t bicomp_root;
1112         vertex_t w = graph_traits< Graph >::null_vertex();
1113         face_handle_t w_handle;
1114         face_handle_t v_dfchild_handle;
1115         vertex_t first_x_y_path_endpoint = graph_traits< Graph >::null_vertex();
1116         vertex_t second_x_y_path_endpoint
1117             = graph_traits< Graph >::null_vertex();
1118         vertex_t w_ancestor = v;
1119 
1120         detail::bm_case_t chosen_case = detail::BM_NO_CASE_CHOSEN;
1121 
1122         std::vector< edge_t > x_external_path;
1123         std::vector< edge_t > y_external_path;
1124         std::vector< edge_t > case_d_edges;
1125 
1126         std::vector< edge_t > z_v_path;
1127         std::vector< edge_t > w_path;
1128 
1129         // first, use a walkup to find a path from V that starts with a
1130         // backedge from V, then goes up until it hits either X or Y
1131         //(but doesn't find X or Y as the root of a bicomp)
1132 
1133         typename face_vertex_iterator<>::type x_upper_itr(
1134             x, face_handles, first_side());
1135         typename face_vertex_iterator<>::type x_lower_itr(
1136             x, face_handles, second_side());
1137         typename face_vertex_iterator<>::type face_itr, face_end;
1138 
1139         // Don't know which path from x is the upper or lower path -
1140         // we'll find out here
1141         for (face_itr = x_upper_itr; face_itr != face_end; ++face_itr)
1142         {
1143             if (*face_itr == y)
1144             {
1145                 std::swap(x_upper_itr, x_lower_itr);
1146                 break;
1147             }
1148         }
1149 
1150         upper_face_vertex[x] = true;
1151 
1152         vertex_t current_vertex = x;
1153         vertex_t previous_vertex;
1154         for (face_itr = x_upper_itr; face_itr != face_end; ++face_itr)
1155         {
1156             previous_vertex = current_vertex;
1157             current_vertex = *face_itr;
1158             upper_face_vertex[current_vertex] = true;
1159         }
1160 
1161         v_dfchild_handle
1162             = dfs_child_handles[canonical_dfs_child[previous_vertex]];
1163 
1164         for (face_itr = x_lower_itr; *face_itr != y; ++face_itr)
1165         {
1166             vertex_t current_vertex(*face_itr);
1167             lower_face_vertex[current_vertex] = true;
1168 
1169             typename face_handle_list_t::iterator roots_itr, roots_end;
1170 
1171             if (w == graph_traits< Graph >::null_vertex()) // haven't found a w
1172                                                            // yet
1173             {
1174                 roots_end = pertinent_roots[current_vertex]->end();
1175                 for (roots_itr = pertinent_roots[current_vertex]->begin();
1176                      roots_itr != roots_end; ++roots_itr)
1177                 {
1178                     if (low_point
1179                             [canonical_dfs_child[roots_itr->first_vertex()]]
1180                         < dfs_number[v])
1181                     {
1182                         w = current_vertex;
1183                         w_handle = *roots_itr;
1184                         break;
1185                     }
1186                 }
1187             }
1188         }
1189 
1190         for (; face_itr != face_end; ++face_itr)
1191         {
1192             vertex_t current_vertex(*face_itr);
1193             upper_face_vertex[current_vertex] = true;
1194             bicomp_root = current_vertex;
1195         }
1196 
1197         typedef typename face_edge_iterator<>::type walkup_itr_t;
1198 
1199         std::vector< bool > outer_face_edge_vector(num_edges(g), false);
1200         edge_to_bool_map_t outer_face_edge(outer_face_edge_vector.begin(), em);
1201 
1202         walkup_itr_t walkup_end;
1203         for (walkup_itr_t walkup_itr(x, face_handles, first_side());
1204              walkup_itr != walkup_end; ++walkup_itr)
1205         {
1206             outer_face_edge[*walkup_itr] = true;
1207             is_in_subgraph[*walkup_itr] = true;
1208         }
1209 
1210         for (walkup_itr_t walkup_itr(x, face_handles, second_side());
1211              walkup_itr != walkup_end; ++walkup_itr)
1212         {
1213             outer_face_edge[*walkup_itr] = true;
1214             is_in_subgraph[*walkup_itr] = true;
1215         }
1216 
1217         std::vector< bool > forbidden_edge_vector(num_edges(g), false);
1218         edge_to_bool_map_t forbidden_edge(forbidden_edge_vector.begin(), em);
1219 
1220         std::vector< bool > goal_edge_vector(num_edges(g), false);
1221         edge_to_bool_map_t goal_edge(goal_edge_vector.begin(), em);
1222 
1223         // Find external path to x and to y
1224 
1225         for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1226         {
1227             edge_t e(*ei);
1228             goal_edge[e] = !outer_face_edge[e]
1229                 && (source(e, g) == x || target(e, g) == x);
1230             forbidden_edge[*ei] = outer_face_edge[*ei];
1231         }
1232 
1233         vertex_t x_ancestor = v;
1234         vertex_t x_endpoint = graph_traits< Graph >::null_vertex();
1235 
1236         while (x_endpoint == graph_traits< Graph >::null_vertex())
1237         {
1238             x_ancestor = dfs_parent[x_ancestor];
1239             x_endpoint = kuratowski_walkup(x_ancestor, forbidden_edge,
1240                 goal_edge, is_embedded, x_external_path);
1241         }
1242 
1243         for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1244         {
1245             edge_t e(*ei);
1246             goal_edge[e] = !outer_face_edge[e]
1247                 && (source(e, g) == y || target(e, g) == y);
1248             forbidden_edge[*ei] = outer_face_edge[*ei];
1249         }
1250 
1251         vertex_t y_ancestor = v;
1252         vertex_t y_endpoint = graph_traits< Graph >::null_vertex();
1253 
1254         while (y_endpoint == graph_traits< Graph >::null_vertex())
1255         {
1256             y_ancestor = dfs_parent[y_ancestor];
1257             y_endpoint = kuratowski_walkup(y_ancestor, forbidden_edge,
1258                 goal_edge, is_embedded, y_external_path);
1259         }
1260 
1261         vertex_t parent, child;
1262 
1263         // If v isn't on the same bicomp as x and y, it's a case A
1264         if (bicomp_root != v)
1265         {
1266             chosen_case = detail::BM_CASE_A;
1267 
1268             for (boost::tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi)
1269                 if (lower_face_vertex[*vi])
1270                     for (boost::tie(oei, oei_end) = out_edges(*vi, g);
1271                          oei != oei_end; ++oei)
1272                         if (!outer_face_edge[*oei])
1273                             goal_edge[*oei] = true;
1274 
1275             for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1276                 forbidden_edge[*ei] = outer_face_edge[*ei];
1277 
1278             z = kuratowski_walkup(
1279                 v, forbidden_edge, goal_edge, is_embedded, z_v_path);
1280         }
1281         else if (w != graph_traits< Graph >::null_vertex())
1282         {
1283             chosen_case = detail::BM_CASE_B;
1284 
1285             for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1286             {
1287                 edge_t e(*ei);
1288                 goal_edge[e] = false;
1289                 forbidden_edge[e] = outer_face_edge[e];
1290             }
1291 
1292             goal_edge[w_handle.first_edge()] = true;
1293             goal_edge[w_handle.second_edge()] = true;
1294 
1295             z = kuratowski_walkup(
1296                 v, forbidden_edge, goal_edge, is_embedded, z_v_path);
1297 
1298             for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1299             {
1300                 forbidden_edge[*ei] = outer_face_edge[*ei];
1301             }
1302 
1303             typename std::vector< edge_t >::iterator pi, pi_end;
1304             pi_end = z_v_path.end();
1305             for (pi = z_v_path.begin(); pi != pi_end; ++pi)
1306             {
1307                 goal_edge[*pi] = true;
1308             }
1309 
1310             w_ancestor = v;
1311             vertex_t w_endpoint = graph_traits< Graph >::null_vertex();
1312 
1313             while (w_endpoint == graph_traits< Graph >::null_vertex())
1314             {
1315                 w_ancestor = dfs_parent[w_ancestor];
1316                 w_endpoint = kuratowski_walkup(
1317                     w_ancestor, forbidden_edge, goal_edge, is_embedded, w_path);
1318             }
1319 
1320             // We really want both the w walkup and the z walkup to finish on
1321             // exactly the same edge, but for convenience (since we don't have
1322             // control over which side of a bicomp a walkup moves up) we've
1323             // defined the walkup to either end at w_handle.first_edge() or
1324             // w_handle.second_edge(). If both walkups ended at different edges,
1325             // we'll do a little surgery on the w walkup path to make it follow
1326             // the other side of the final bicomp.
1327 
1328             if ((w_path.back() == w_handle.first_edge()
1329                     && z_v_path.back() == w_handle.second_edge())
1330                 || (w_path.back() == w_handle.second_edge()
1331                     && z_v_path.back() == w_handle.first_edge()))
1332             {
1333                 walkup_itr_t wi, wi_end;
1334                 edge_t final_edge = w_path.back();
1335                 vertex_t anchor = source(final_edge, g) == w_handle.get_anchor()
1336                     ? target(final_edge, g)
1337                     : source(final_edge, g);
1338                 if (face_handles[anchor].first_edge() == final_edge)
1339                     wi = walkup_itr_t(anchor, face_handles, second_side());
1340                 else
1341                     wi = walkup_itr_t(anchor, face_handles, first_side());
1342 
1343                 w_path.pop_back();
1344 
1345                 for (; wi != wi_end; ++wi)
1346                 {
1347                     edge_t e(*wi);
1348                     if (w_path.back() == e)
1349                         w_path.pop_back();
1350                     else
1351                         w_path.push_back(e);
1352                 }
1353             }
1354         }
1355         else
1356         {
1357 
1358             // We need to find a valid z, since the x-y path re-defines the
1359             // lower face, and the z we found earlier may now be on the upper
1360             // face.
1361 
1362             chosen_case = detail::BM_CASE_E;
1363 
1364             // The z we've used so far is just an externally active vertex on
1365             // the lower face path, but may not be the z we need for a case C,
1366             // D, or E subgraph. the z we need now is any externally active
1367             // vertex on the lower face path with both old_face_handles edges on
1368             // the outer face. Since we know an x-y path exists, such a z must
1369             // also exist.
1370 
1371             // TODO: find this z in the first place.
1372 
1373             // find the new z
1374 
1375             for (face_itr = x_lower_itr; *face_itr != y; ++face_itr)
1376             {
1377                 vertex_t possible_z(*face_itr);
1378                 if (pertinent(possible_z, v)
1379                     && outer_face_edge[face_handles[possible_z]
1380                                            .old_first_edge()]
1381                     && outer_face_edge[face_handles[possible_z]
1382                                            .old_second_edge()])
1383                 {
1384                     z = possible_z;
1385                     break;
1386                 }
1387             }
1388 
1389             // find x-y path, and a w if one exists.
1390 
1391             if (externally_active(z, v))
1392                 w = z;
1393 
1394             typedef typename face_edge_iterator< single_side,
1395                 previous_iteration >::type old_face_iterator_t;
1396 
1397             old_face_iterator_t first_old_face_itr(
1398                 z, face_handles, first_side());
1399             old_face_iterator_t second_old_face_itr(
1400                 z, face_handles, second_side());
1401             old_face_iterator_t old_face_itr, old_face_end;
1402 
1403             std::vector< old_face_iterator_t > old_face_iterators;
1404             old_face_iterators.push_back(first_old_face_itr);
1405             old_face_iterators.push_back(second_old_face_itr);
1406 
1407             std::vector< bool > x_y_path_vertex_vector(num_vertices(g), false);
1408             vertex_to_bool_map_t x_y_path_vertex(
1409                 x_y_path_vertex_vector.begin(), vm);
1410 
1411             typename std::vector< old_face_iterator_t >::iterator of_itr,
1412                 of_itr_end;
1413             of_itr_end = old_face_iterators.end();
1414             for (of_itr = old_face_iterators.begin(); of_itr != of_itr_end;
1415                  ++of_itr)
1416             {
1417 
1418                 old_face_itr = *of_itr;
1419 
1420                 vertex_t previous_vertex;
1421                 bool seen_x_or_y = false;
1422                 vertex_t current_vertex = z;
1423                 for (; old_face_itr != old_face_end; ++old_face_itr)
1424                 {
1425                     edge_t e(*old_face_itr);
1426                     previous_vertex = current_vertex;
1427                     current_vertex = source(e, g) == current_vertex
1428                         ? target(e, g)
1429                         : source(e, g);
1430 
1431                     if (current_vertex == x || current_vertex == y)
1432                         seen_x_or_y = true;
1433 
1434                     if (w == graph_traits< Graph >::null_vertex()
1435                         && externally_active(current_vertex, v)
1436                         && outer_face_edge[e]
1437                         && outer_face_edge[*boost::next(old_face_itr)]
1438                         && !seen_x_or_y)
1439                     {
1440                         w = current_vertex;
1441                     }
1442 
1443                     if (!outer_face_edge[e])
1444                     {
1445                         if (!upper_face_vertex[current_vertex]
1446                             && !lower_face_vertex[current_vertex])
1447                         {
1448                             x_y_path_vertex[current_vertex] = true;
1449                         }
1450 
1451                         is_in_subgraph[e] = true;
1452                         if (upper_face_vertex[source(e, g)]
1453                             || lower_face_vertex[source(e, g)])
1454                         {
1455                             if (first_x_y_path_endpoint
1456                                 == graph_traits< Graph >::null_vertex())
1457                                 first_x_y_path_endpoint = source(e, g);
1458                             else
1459                                 second_x_y_path_endpoint = source(e, g);
1460                         }
1461                         if (upper_face_vertex[target(e, g)]
1462                             || lower_face_vertex[target(e, g)])
1463                         {
1464                             if (first_x_y_path_endpoint
1465                                 == graph_traits< Graph >::null_vertex())
1466                                 first_x_y_path_endpoint = target(e, g);
1467                             else
1468                                 second_x_y_path_endpoint = target(e, g);
1469                         }
1470                     }
1471                     else if (previous_vertex == x || previous_vertex == y)
1472                     {
1473                         chosen_case = detail::BM_CASE_C;
1474                     }
1475                 }
1476             }
1477 
1478             // Look for a case D - one of v's embedded edges will connect to the
1479             // x-y path along an inner face path.
1480 
1481             // First, get a list of all of v's embedded child edges
1482 
1483             out_edge_iterator_t v_edge_itr, v_edge_end;
1484             for (boost::tie(v_edge_itr, v_edge_end) = out_edges(v, g);
1485                  v_edge_itr != v_edge_end; ++v_edge_itr)
1486             {
1487                 edge_t embedded_edge(*v_edge_itr);
1488 
1489                 if (!is_embedded[embedded_edge]
1490                     || embedded_edge == dfs_parent_edge[v])
1491                     continue;
1492 
1493                 case_d_edges.push_back(embedded_edge);
1494 
1495                 vertex_t current_vertex = source(embedded_edge, g) == v
1496                     ? target(embedded_edge, g)
1497                     : source(embedded_edge, g);
1498 
1499                 typename face_edge_iterator<>::type internal_face_itr,
1500                     internal_face_end;
1501                 if (face_handles[current_vertex].first_vertex() == v)
1502                 {
1503                     internal_face_itr = typename face_edge_iterator<>::type(
1504                         current_vertex, face_handles, second_side());
1505                 }
1506                 else
1507                 {
1508                     internal_face_itr = typename face_edge_iterator<>::type(
1509                         current_vertex, face_handles, first_side());
1510                 }
1511 
1512                 while (internal_face_itr != internal_face_end
1513                     && !outer_face_edge[*internal_face_itr]
1514                     && !x_y_path_vertex[current_vertex])
1515                 {
1516                     edge_t e(*internal_face_itr);
1517                     case_d_edges.push_back(e);
1518                     current_vertex = source(e, g) == current_vertex
1519                         ? target(e, g)
1520                         : source(e, g);
1521                     ++internal_face_itr;
1522                 }
1523 
1524                 if (x_y_path_vertex[current_vertex])
1525                 {
1526                     chosen_case = detail::BM_CASE_D;
1527                     break;
1528                 }
1529                 else
1530                 {
1531                     case_d_edges.clear();
1532                 }
1533             }
1534         }
1535 
1536         if (chosen_case != detail::BM_CASE_B
1537             && chosen_case != detail::BM_CASE_A)
1538         {
1539 
1540             // Finding z and w.
1541 
1542             for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1543             {
1544                 edge_t e(*ei);
1545                 goal_edge[e] = !outer_face_edge[e]
1546                     && (source(e, g) == z || target(e, g) == z);
1547                 forbidden_edge[e] = outer_face_edge[e];
1548             }
1549 
1550             kuratowski_walkup(
1551                 v, forbidden_edge, goal_edge, is_embedded, z_v_path);
1552 
1553             if (chosen_case == detail::BM_CASE_E)
1554             {
1555 
1556                 for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1557                 {
1558                     forbidden_edge[*ei] = outer_face_edge[*ei];
1559                     goal_edge[*ei] = !outer_face_edge[*ei]
1560                         && (source(*ei, g) == w || target(*ei, g) == w);
1561                 }
1562 
1563                 for (boost::tie(oei, oei_end) = out_edges(w, g); oei != oei_end;
1564                      ++oei)
1565                 {
1566                     if (!outer_face_edge[*oei])
1567                         goal_edge[*oei] = true;
1568                 }
1569 
1570                 typename std::vector< edge_t >::iterator pi, pi_end;
1571                 pi_end = z_v_path.end();
1572                 for (pi = z_v_path.begin(); pi != pi_end; ++pi)
1573                 {
1574                     goal_edge[*pi] = true;
1575                 }
1576 
1577                 w_ancestor = v;
1578                 vertex_t w_endpoint = graph_traits< Graph >::null_vertex();
1579 
1580                 while (w_endpoint == graph_traits< Graph >::null_vertex())
1581                 {
1582                     w_ancestor = dfs_parent[w_ancestor];
1583                     w_endpoint = kuratowski_walkup(w_ancestor, forbidden_edge,
1584                         goal_edge, is_embedded, w_path);
1585                 }
1586             }
1587         }
1588 
1589         // We're done isolating the Kuratowski subgraph at this point -
1590         // but there's still some cleaning up to do.
1591 
1592         // Update is_in_subgraph with the paths we just found
1593 
1594         xi_end = x_external_path.end();
1595         for (xi = x_external_path.begin(); xi != xi_end; ++xi)
1596             is_in_subgraph[*xi] = true;
1597 
1598         xi_end = y_external_path.end();
1599         for (xi = y_external_path.begin(); xi != xi_end; ++xi)
1600             is_in_subgraph[*xi] = true;
1601 
1602         xi_end = z_v_path.end();
1603         for (xi = z_v_path.begin(); xi != xi_end; ++xi)
1604             is_in_subgraph[*xi] = true;
1605 
1606         xi_end = case_d_edges.end();
1607         for (xi = case_d_edges.begin(); xi != xi_end; ++xi)
1608             is_in_subgraph[*xi] = true;
1609 
1610         xi_end = w_path.end();
1611         for (xi = w_path.begin(); xi != xi_end; ++xi)
1612             is_in_subgraph[*xi] = true;
1613 
1614         child = bicomp_root;
1615         parent = dfs_parent[child];
1616         while (child != parent)
1617         {
1618             is_in_subgraph[dfs_parent_edge[child]] = true;
1619             boost::tie(parent, child)
1620                 = std::make_pair(dfs_parent[parent], parent);
1621         }
1622 
1623         // At this point, we've already isolated the Kuratowski subgraph and
1624         // collected all of the edges that compose it in the is_in_subgraph
1625         // property map. But we want the verification of such a subgraph to be
1626         // a deterministic process, and we can simplify the function
1627         // is_kuratowski_subgraph by cleaning up some edges here.
1628 
1629         if (chosen_case == detail::BM_CASE_B)
1630         {
1631             is_in_subgraph[dfs_parent_edge[v]] = false;
1632         }
1633         else if (chosen_case == detail::BM_CASE_C)
1634         {
1635             // In a case C subgraph, at least one of the x-y path endpoints
1636             // (call it alpha) is above either x or y on the outer face. The
1637             // other endpoint may be attached at x or y OR above OR below. In
1638             // any of these three cases, we can form a K_3_3 by removing the
1639             // edge attached to v on the outer face that is NOT on the path to
1640             // alpha.
1641 
1642             typename face_vertex_iterator< single_side, follow_visitor >::type
1643                 face_itr,
1644                 face_end;
1645             if (face_handles[v_dfchild_handle.first_vertex()].first_edge()
1646                 == v_dfchild_handle.first_edge())
1647             {
1648                 face_itr = typename face_vertex_iterator< single_side,
1649                     follow_visitor >::type(v_dfchild_handle.first_vertex(),
1650                     face_handles, second_side());
1651             }
1652             else
1653             {
1654                 face_itr = typename face_vertex_iterator< single_side,
1655                     follow_visitor >::type(v_dfchild_handle.first_vertex(),
1656                     face_handles, first_side());
1657             }
1658 
1659             for (; true; ++face_itr)
1660             {
1661                 vertex_t current_vertex(*face_itr);
1662                 if (current_vertex == x || current_vertex == y)
1663                 {
1664                     is_in_subgraph[v_dfchild_handle.first_edge()] = false;
1665                     break;
1666                 }
1667                 else if (current_vertex == first_x_y_path_endpoint
1668                     || current_vertex == second_x_y_path_endpoint)
1669                 {
1670                     is_in_subgraph[v_dfchild_handle.second_edge()] = false;
1671                     break;
1672                 }
1673             }
1674         }
1675         else if (chosen_case == detail::BM_CASE_D)
1676         {
1677             // Need to remove both of the edges adjacent to v on the outer face.
1678             // remove the connecting edges from v to bicomp, then
1679             // is_kuratowski_subgraph will shrink vertices of degree 1
1680             // automatically...
1681 
1682             is_in_subgraph[v_dfchild_handle.first_edge()] = false;
1683             is_in_subgraph[v_dfchild_handle.second_edge()] = false;
1684         }
1685         else if (chosen_case == detail::BM_CASE_E)
1686         {
1687             // Similarly to case C, if the endpoints of the x-y path are both
1688             // below x and y, we should remove an edge to allow the subgraph to
1689             // contract to a K_3_3.
1690 
1691             if ((first_x_y_path_endpoint != x && first_x_y_path_endpoint != y)
1692                 || (second_x_y_path_endpoint != x
1693                     && second_x_y_path_endpoint != y))
1694             {
1695                 is_in_subgraph[dfs_parent_edge[v]] = false;
1696 
1697                 vertex_t deletion_endpoint, other_endpoint;
1698                 if (lower_face_vertex[first_x_y_path_endpoint])
1699                 {
1700                     deletion_endpoint = second_x_y_path_endpoint;
1701                     other_endpoint = first_x_y_path_endpoint;
1702                 }
1703                 else
1704                 {
1705                     deletion_endpoint = first_x_y_path_endpoint;
1706                     other_endpoint = second_x_y_path_endpoint;
1707                 }
1708 
1709                 typename face_edge_iterator<>::type face_itr, face_end;
1710 
1711                 bool found_other_endpoint = false;
1712                 for (face_itr = typename face_edge_iterator<>::type(
1713                          deletion_endpoint, face_handles, first_side());
1714                      face_itr != face_end; ++face_itr)
1715                 {
1716                     edge_t e(*face_itr);
1717                     if (source(e, g) == other_endpoint
1718                         || target(e, g) == other_endpoint)
1719                     {
1720                         found_other_endpoint = true;
1721                         break;
1722                     }
1723                 }
1724 
1725                 if (found_other_endpoint)
1726                 {
1727                     is_in_subgraph[face_handles[deletion_endpoint].first_edge()]
1728                         = false;
1729                 }
1730                 else
1731                 {
1732                     is_in_subgraph[face_handles[deletion_endpoint]
1733                                        .second_edge()]
1734                         = false;
1735                 }
1736             }
1737         }
1738 
1739         for (boost::tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1740             if (is_in_subgraph[*ei])
1741                 *o_itr = *ei;
1742     }
1743 
1744     template < typename EdgePermutation >
make_edge_permutation(EdgePermutation perm)1745     void make_edge_permutation(EdgePermutation perm)
1746     {
1747         vertex_iterator_t vi, vi_end;
1748         for (boost::tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi)
1749         {
1750             vertex_t v(*vi);
1751             perm[v].clear();
1752             face_handles[v].get_list(std::back_inserter(perm[v]));
1753         }
1754     }
1755 
1756 private:
1757     const Graph& g;
1758     VertexIndexMap vm;
1759 
1760     vertex_t kuratowski_v;
1761     vertex_t kuratowski_x;
1762     vertex_t kuratowski_y;
1763 
1764     vertex_list_t garbage; // we delete items from linked lists by
1765                            // splicing them into garbage
1766 
1767     // only need these two for kuratowski subgraph isolation
1768     std::vector< vertex_t > current_merge_points;
1769     std::vector< edge_t > embedded_edges;
1770 
1771     // property map storage
1772     std::vector< v_size_t > low_point_vector;
1773     std::vector< vertex_t > dfs_parent_vector;
1774     std::vector< v_size_t > dfs_number_vector;
1775     std::vector< v_size_t > least_ancestor_vector;
1776     std::vector< face_handle_list_ptr_t > pertinent_roots_vector;
1777     std::vector< v_size_t > backedge_flag_vector;
1778     std::vector< v_size_t > visited_vector;
1779     std::vector< face_handle_t > face_handles_vector;
1780     std::vector< face_handle_t > dfs_child_handles_vector;
1781     std::vector< vertex_list_ptr_t > separated_dfs_child_list_vector;
1782     std::vector< typename vertex_list_t::iterator >
1783         separated_node_in_parent_list_vector;
1784     std::vector< vertex_t > canonical_dfs_child_vector;
1785     std::vector< bool > flipped_vector;
1786     std::vector< edge_vector_t > backedges_vector;
1787     edge_vector_t self_loops;
1788     std::vector< edge_t > dfs_parent_edge_vector;
1789     vertex_vector_t vertices_by_dfs_num;
1790 
1791     // property maps
1792     vertex_to_v_size_map_t low_point;
1793     vertex_to_vertex_map_t dfs_parent;
1794     vertex_to_v_size_map_t dfs_number;
1795     vertex_to_v_size_map_t least_ancestor;
1796     vertex_to_face_handle_list_ptr_map_t pertinent_roots;
1797     vertex_to_v_size_map_t backedge_flag;
1798     vertex_to_v_size_map_t visited;
1799     vertex_to_face_handle_map_t face_handles;
1800     vertex_to_face_handle_map_t dfs_child_handles;
1801     vertex_to_vertex_list_ptr_map_t separated_dfs_child_list;
1802     vertex_to_separated_node_map_t separated_node_in_parent_list;
1803     vertex_to_vertex_map_t canonical_dfs_child;
1804     vertex_to_bool_map_t flipped;
1805     vertex_to_edge_vector_map_t backedges;
1806     vertex_to_edge_map_t dfs_parent_edge; // only need for kuratowski
1807 
1808     merge_stack_t merge_stack;
1809 };
1810 
1811 } // namespace boost
1812 
1813 #endif //__BOYER_MYRVOLD_IMPL_HPP__
1814