lemon/planarity.h
changeset 797 30cb42e3e43a
child 798 58c330ad0b5c
equal deleted inserted replaced
-1:000000000000 0:c1e331e61456
       
     1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
       
     2  *
       
     3  * This file is a part of LEMON, a generic C++ optimization library.
       
     4  *
       
     5  * Copyright (C) 2003-2009
       
     6  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
       
     7  * (Egervary Research Group on Combinatorial Optimization, EGRES).
       
     8  *
       
     9  * Permission to use, modify and distribute this software is granted
       
    10  * provided that this copyright notice appears in all copies. For
       
    11  * precise terms see the accompanying LICENSE file.
       
    12  *
       
    13  * This software is provided "AS IS" with no warranty of any kind,
       
    14  * express or implied, and with no claim as to its suitability for any
       
    15  * purpose.
       
    16  *
       
    17  */
       
    18 
       
    19 #ifndef LEMON_PLANARITY_H
       
    20 #define LEMON_PLANARITY_H
       
    21 
       
    22 /// \ingroup planar
       
    23 /// \file
       
    24 /// \brief Planarity checking, embedding, drawing and coloring
       
    25 
       
    26 #include <vector>
       
    27 #include <list>
       
    28 
       
    29 #include <lemon/dfs.h>
       
    30 #include <lemon/bfs.h>
       
    31 #include <lemon/radix_sort.h>
       
    32 #include <lemon/maps.h>
       
    33 #include <lemon/path.h>
       
    34 #include <lemon/bucket_heap.h>
       
    35 #include <lemon/adaptors.h>
       
    36 #include <lemon/edge_set.h>
       
    37 #include <lemon/color.h>
       
    38 #include <lemon/dim2.h>
       
    39 
       
    40 namespace lemon {
       
    41 
       
    42   namespace _planarity_bits {
       
    43 
       
    44     template <typename Graph>
       
    45     struct PlanarityVisitor : DfsVisitor<Graph> {
       
    46 
       
    47       TEMPLATE_GRAPH_TYPEDEFS(Graph);
       
    48 
       
    49       typedef typename Graph::template NodeMap<Arc> PredMap;
       
    50 
       
    51       typedef typename Graph::template EdgeMap<bool> TreeMap;
       
    52 
       
    53       typedef typename Graph::template NodeMap<int> OrderMap;
       
    54       typedef std::vector<Node> OrderList;
       
    55 
       
    56       typedef typename Graph::template NodeMap<int> LowMap;
       
    57       typedef typename Graph::template NodeMap<int> AncestorMap;
       
    58 
       
    59       PlanarityVisitor(const Graph& graph,
       
    60                        PredMap& pred_map, TreeMap& tree_map,
       
    61                        OrderMap& order_map, OrderList& order_list,
       
    62                        AncestorMap& ancestor_map, LowMap& low_map)
       
    63         : _graph(graph), _pred_map(pred_map), _tree_map(tree_map),
       
    64           _order_map(order_map), _order_list(order_list),
       
    65           _ancestor_map(ancestor_map), _low_map(low_map) {}
       
    66 
       
    67       void reach(const Node& node) {
       
    68         _order_map[node] = _order_list.size();
       
    69         _low_map[node] = _order_list.size();
       
    70         _ancestor_map[node] = _order_list.size();
       
    71         _order_list.push_back(node);
       
    72       }
       
    73 
       
    74       void discover(const Arc& arc) {
       
    75         Node source = _graph.source(arc);
       
    76         Node target = _graph.target(arc);
       
    77 
       
    78         _tree_map[arc] = true;
       
    79         _pred_map[target] = arc;
       
    80       }
       
    81 
       
    82       void examine(const Arc& arc) {
       
    83         Node source = _graph.source(arc);
       
    84         Node target = _graph.target(arc);
       
    85 
       
    86         if (_order_map[target] < _order_map[source] && !_tree_map[arc]) {
       
    87           if (_low_map[source] > _order_map[target]) {
       
    88             _low_map[source] = _order_map[target];
       
    89           }
       
    90           if (_ancestor_map[source] > _order_map[target]) {
       
    91             _ancestor_map[source] = _order_map[target];
       
    92           }
       
    93         }
       
    94       }
       
    95 
       
    96       void backtrack(const Arc& arc) {
       
    97         Node source = _graph.source(arc);
       
    98         Node target = _graph.target(arc);
       
    99 
       
   100         if (_low_map[source] > _low_map[target]) {
       
   101           _low_map[source] = _low_map[target];
       
   102         }
       
   103       }
       
   104 
       
   105       const Graph& _graph;
       
   106       PredMap& _pred_map;
       
   107       TreeMap& _tree_map;
       
   108       OrderMap& _order_map;
       
   109       OrderList& _order_list;
       
   110       AncestorMap& _ancestor_map;
       
   111       LowMap& _low_map;
       
   112     };
       
   113 
       
   114     template <typename Graph, bool embedding = true>
       
   115     struct NodeDataNode {
       
   116       int prev, next;
       
   117       int visited;
       
   118       typename Graph::Arc first;
       
   119       bool inverted;
       
   120     };
       
   121 
       
   122     template <typename Graph>
       
   123     struct NodeDataNode<Graph, false> {
       
   124       int prev, next;
       
   125       int visited;
       
   126     };
       
   127 
       
   128     template <typename Graph>
       
   129     struct ChildListNode {
       
   130       typedef typename Graph::Node Node;
       
   131       Node first;
       
   132       Node prev, next;
       
   133     };
       
   134 
       
   135     template <typename Graph>
       
   136     struct ArcListNode {
       
   137       typename Graph::Arc prev, next;
       
   138     };
       
   139 
       
   140   }
       
   141 
       
   142   /// \ingroup planar
       
   143   ///
       
   144   /// \brief Planarity checking of an undirected simple graph
       
   145   ///
       
   146   /// This class implements the Boyer-Myrvold algorithm for planarity
       
   147   /// checking of an undirected graph. This class is a simplified
       
   148   /// version of the PlanarEmbedding algorithm class because neither
       
   149   /// the embedding nor the kuratowski subdivisons are not computed.
       
   150   template <typename Graph>
       
   151   class PlanarityChecking {
       
   152   private:
       
   153 
       
   154     TEMPLATE_GRAPH_TYPEDEFS(Graph);
       
   155 
       
   156     const Graph& _graph;
       
   157 
       
   158   private:
       
   159 
       
   160     typedef typename Graph::template NodeMap<Arc> PredMap;
       
   161 
       
   162     typedef typename Graph::template EdgeMap<bool> TreeMap;
       
   163 
       
   164     typedef typename Graph::template NodeMap<int> OrderMap;
       
   165     typedef std::vector<Node> OrderList;
       
   166 
       
   167     typedef typename Graph::template NodeMap<int> LowMap;
       
   168     typedef typename Graph::template NodeMap<int> AncestorMap;
       
   169 
       
   170     typedef _planarity_bits::NodeDataNode<Graph> NodeDataNode;
       
   171     typedef std::vector<NodeDataNode> NodeData;
       
   172 
       
   173     typedef _planarity_bits::ChildListNode<Graph> ChildListNode;
       
   174     typedef typename Graph::template NodeMap<ChildListNode> ChildLists;
       
   175 
       
   176     typedef typename Graph::template NodeMap<std::list<int> > MergeRoots;
       
   177 
       
   178     typedef typename Graph::template NodeMap<bool> EmbedArc;
       
   179 
       
   180   public:
       
   181 
       
   182     /// \brief Constructor
       
   183     ///
       
   184     /// \note The graph should be simple, i.e. parallel and loop arc
       
   185     /// free.
       
   186     PlanarityChecking(const Graph& graph) : _graph(graph) {}
       
   187 
       
   188     /// \brief Runs the algorithm.
       
   189     ///
       
   190     /// Runs the algorithm.
       
   191     /// \return %True when the graph is planar.
       
   192     bool run() {
       
   193       typedef _planarity_bits::PlanarityVisitor<Graph> Visitor;
       
   194 
       
   195       PredMap pred_map(_graph, INVALID);
       
   196       TreeMap tree_map(_graph, false);
       
   197 
       
   198       OrderMap order_map(_graph, -1);
       
   199       OrderList order_list;
       
   200 
       
   201       AncestorMap ancestor_map(_graph, -1);
       
   202       LowMap low_map(_graph, -1);
       
   203 
       
   204       Visitor visitor(_graph, pred_map, tree_map,
       
   205                       order_map, order_list, ancestor_map, low_map);
       
   206       DfsVisit<Graph, Visitor> visit(_graph, visitor);
       
   207       visit.run();
       
   208 
       
   209       ChildLists child_lists(_graph);
       
   210       createChildLists(tree_map, order_map, low_map, child_lists);
       
   211 
       
   212       NodeData node_data(2 * order_list.size());
       
   213 
       
   214       EmbedArc embed_arc(_graph, false);
       
   215 
       
   216       MergeRoots merge_roots(_graph);
       
   217 
       
   218       for (int i = order_list.size() - 1; i >= 0; --i) {
       
   219 
       
   220         Node node = order_list[i];
       
   221 
       
   222         Node source = node;
       
   223         for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
   224           Node target = _graph.target(e);
       
   225 
       
   226           if (order_map[source] < order_map[target] && tree_map[e]) {
       
   227             initFace(target, node_data, order_map, order_list);
       
   228           }
       
   229         }
       
   230 
       
   231         for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
   232           Node target = _graph.target(e);
       
   233 
       
   234           if (order_map[source] < order_map[target] && !tree_map[e]) {
       
   235             embed_arc[target] = true;
       
   236             walkUp(target, source, i, pred_map, low_map,
       
   237                    order_map, order_list, node_data, merge_roots);
       
   238           }
       
   239         }
       
   240 
       
   241         for (typename MergeRoots::Value::iterator it =
       
   242                merge_roots[node].begin(); it != merge_roots[node].end(); ++it) {
       
   243           int rn = *it;
       
   244           walkDown(rn, i, node_data, order_list, child_lists,
       
   245                    ancestor_map, low_map, embed_arc, merge_roots);
       
   246         }
       
   247         merge_roots[node].clear();
       
   248 
       
   249         for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
   250           Node target = _graph.target(e);
       
   251 
       
   252           if (order_map[source] < order_map[target] && !tree_map[e]) {
       
   253             if (embed_arc[target]) {
       
   254               return false;
       
   255             }
       
   256           }
       
   257         }
       
   258       }
       
   259 
       
   260       return true;
       
   261     }
       
   262 
       
   263   private:
       
   264 
       
   265     void createChildLists(const TreeMap& tree_map, const OrderMap& order_map,
       
   266                           const LowMap& low_map, ChildLists& child_lists) {
       
   267 
       
   268       for (NodeIt n(_graph); n != INVALID; ++n) {
       
   269         Node source = n;
       
   270 
       
   271         std::vector<Node> targets;
       
   272         for (OutArcIt e(_graph, n); e != INVALID; ++e) {
       
   273           Node target = _graph.target(e);
       
   274 
       
   275           if (order_map[source] < order_map[target] && tree_map[e]) {
       
   276             targets.push_back(target);
       
   277           }
       
   278         }
       
   279 
       
   280         if (targets.size() == 0) {
       
   281           child_lists[source].first = INVALID;
       
   282         } else if (targets.size() == 1) {
       
   283           child_lists[source].first = targets[0];
       
   284           child_lists[targets[0]].prev = INVALID;
       
   285           child_lists[targets[0]].next = INVALID;
       
   286         } else {
       
   287           radixSort(targets.begin(), targets.end(), mapToFunctor(low_map));
       
   288           for (int i = 1; i < int(targets.size()); ++i) {
       
   289             child_lists[targets[i]].prev = targets[i - 1];
       
   290             child_lists[targets[i - 1]].next = targets[i];
       
   291           }
       
   292           child_lists[targets.back()].next = INVALID;
       
   293           child_lists[targets.front()].prev = INVALID;
       
   294           child_lists[source].first = targets.front();
       
   295         }
       
   296       }
       
   297     }
       
   298 
       
   299     void walkUp(const Node& node, Node root, int rorder,
       
   300                 const PredMap& pred_map, const LowMap& low_map,
       
   301                 const OrderMap& order_map, const OrderList& order_list,
       
   302                 NodeData& node_data, MergeRoots& merge_roots) {
       
   303 
       
   304       int na, nb;
       
   305       bool da, db;
       
   306 
       
   307       na = nb = order_map[node];
       
   308       da = true; db = false;
       
   309 
       
   310       while (true) {
       
   311 
       
   312         if (node_data[na].visited == rorder) break;
       
   313         if (node_data[nb].visited == rorder) break;
       
   314 
       
   315         node_data[na].visited = rorder;
       
   316         node_data[nb].visited = rorder;
       
   317 
       
   318         int rn = -1;
       
   319 
       
   320         if (na >= int(order_list.size())) {
       
   321           rn = na;
       
   322         } else if (nb >= int(order_list.size())) {
       
   323           rn = nb;
       
   324         }
       
   325 
       
   326         if (rn == -1) {
       
   327           int nn;
       
   328 
       
   329           nn = da ? node_data[na].prev : node_data[na].next;
       
   330           da = node_data[nn].prev != na;
       
   331           na = nn;
       
   332 
       
   333           nn = db ? node_data[nb].prev : node_data[nb].next;
       
   334           db = node_data[nn].prev != nb;
       
   335           nb = nn;
       
   336 
       
   337         } else {
       
   338 
       
   339           Node rep = order_list[rn - order_list.size()];
       
   340           Node parent = _graph.source(pred_map[rep]);
       
   341 
       
   342           if (low_map[rep] < rorder) {
       
   343             merge_roots[parent].push_back(rn);
       
   344           } else {
       
   345             merge_roots[parent].push_front(rn);
       
   346           }
       
   347 
       
   348           if (parent != root) {
       
   349             na = nb = order_map[parent];
       
   350             da = true; db = false;
       
   351           } else {
       
   352             break;
       
   353           }
       
   354         }
       
   355       }
       
   356     }
       
   357 
       
   358     void walkDown(int rn, int rorder, NodeData& node_data,
       
   359                   OrderList& order_list, ChildLists& child_lists,
       
   360                   AncestorMap& ancestor_map, LowMap& low_map,
       
   361                   EmbedArc& embed_arc, MergeRoots& merge_roots) {
       
   362 
       
   363       std::vector<std::pair<int, bool> > merge_stack;
       
   364 
       
   365       for (int di = 0; di < 2; ++di) {
       
   366         bool rd = di == 0;
       
   367         int pn = rn;
       
   368         int n = rd ? node_data[rn].next : node_data[rn].prev;
       
   369 
       
   370         while (n != rn) {
       
   371 
       
   372           Node node = order_list[n];
       
   373 
       
   374           if (embed_arc[node]) {
       
   375 
       
   376             // Merging components on the critical path
       
   377             while (!merge_stack.empty()) {
       
   378 
       
   379               // Component root
       
   380               int cn = merge_stack.back().first;
       
   381               bool cd = merge_stack.back().second;
       
   382               merge_stack.pop_back();
       
   383 
       
   384               // Parent of component
       
   385               int dn = merge_stack.back().first;
       
   386               bool dd = merge_stack.back().second;
       
   387               merge_stack.pop_back();
       
   388 
       
   389               Node parent = order_list[dn];
       
   390 
       
   391               // Erasing from merge_roots
       
   392               merge_roots[parent].pop_front();
       
   393 
       
   394               Node child = order_list[cn - order_list.size()];
       
   395 
       
   396               // Erasing from child_lists
       
   397               if (child_lists[child].prev != INVALID) {
       
   398                 child_lists[child_lists[child].prev].next =
       
   399                   child_lists[child].next;
       
   400               } else {
       
   401                 child_lists[parent].first = child_lists[child].next;
       
   402               }
       
   403 
       
   404               if (child_lists[child].next != INVALID) {
       
   405                 child_lists[child_lists[child].next].prev =
       
   406                   child_lists[child].prev;
       
   407               }
       
   408 
       
   409               // Merging external faces
       
   410               {
       
   411                 int en = cn;
       
   412                 cn = cd ? node_data[cn].prev : node_data[cn].next;
       
   413                 cd = node_data[cn].next == en;
       
   414 
       
   415               }
       
   416 
       
   417               if (cd) node_data[cn].next = dn; else node_data[cn].prev = dn;
       
   418               if (dd) node_data[dn].prev = cn; else node_data[dn].next = cn;
       
   419 
       
   420             }
       
   421 
       
   422             bool d = pn == node_data[n].prev;
       
   423 
       
   424             if (node_data[n].prev == node_data[n].next &&
       
   425                 node_data[n].inverted) {
       
   426               d = !d;
       
   427             }
       
   428 
       
   429             // Embedding arc into external face
       
   430             if (rd) node_data[rn].next = n; else node_data[rn].prev = n;
       
   431             if (d) node_data[n].prev = rn; else node_data[n].next = rn;
       
   432             pn = rn;
       
   433 
       
   434             embed_arc[order_list[n]] = false;
       
   435           }
       
   436 
       
   437           if (!merge_roots[node].empty()) {
       
   438 
       
   439             bool d = pn == node_data[n].prev;
       
   440 
       
   441             merge_stack.push_back(std::make_pair(n, d));
       
   442 
       
   443             int rn = merge_roots[node].front();
       
   444 
       
   445             int xn = node_data[rn].next;
       
   446             Node xnode = order_list[xn];
       
   447 
       
   448             int yn = node_data[rn].prev;
       
   449             Node ynode = order_list[yn];
       
   450 
       
   451             bool rd;
       
   452             if (!external(xnode, rorder, child_lists, ancestor_map, low_map)) {
       
   453               rd = true;
       
   454             } else if (!external(ynode, rorder, child_lists,
       
   455                                  ancestor_map, low_map)) {
       
   456               rd = false;
       
   457             } else if (pertinent(xnode, embed_arc, merge_roots)) {
       
   458               rd = true;
       
   459             } else {
       
   460               rd = false;
       
   461             }
       
   462 
       
   463             merge_stack.push_back(std::make_pair(rn, rd));
       
   464 
       
   465             pn = rn;
       
   466             n = rd ? xn : yn;
       
   467 
       
   468           } else if (!external(node, rorder, child_lists,
       
   469                                ancestor_map, low_map)) {
       
   470             int nn = (node_data[n].next != pn ?
       
   471                       node_data[n].next : node_data[n].prev);
       
   472 
       
   473             bool nd = n == node_data[nn].prev;
       
   474 
       
   475             if (nd) node_data[nn].prev = pn;
       
   476             else node_data[nn].next = pn;
       
   477 
       
   478             if (n == node_data[pn].prev) node_data[pn].prev = nn;
       
   479             else node_data[pn].next = nn;
       
   480 
       
   481             node_data[nn].inverted =
       
   482               (node_data[nn].prev == node_data[nn].next && nd != rd);
       
   483 
       
   484             n = nn;
       
   485           }
       
   486           else break;
       
   487 
       
   488         }
       
   489 
       
   490         if (!merge_stack.empty() || n == rn) {
       
   491           break;
       
   492         }
       
   493       }
       
   494     }
       
   495 
       
   496     void initFace(const Node& node, NodeData& node_data,
       
   497                   const OrderMap& order_map, const OrderList& order_list) {
       
   498       int n = order_map[node];
       
   499       int rn = n + order_list.size();
       
   500 
       
   501       node_data[n].next = node_data[n].prev = rn;
       
   502       node_data[rn].next = node_data[rn].prev = n;
       
   503 
       
   504       node_data[n].visited = order_list.size();
       
   505       node_data[rn].visited = order_list.size();
       
   506 
       
   507     }
       
   508 
       
   509     bool external(const Node& node, int rorder,
       
   510                   ChildLists& child_lists, AncestorMap& ancestor_map,
       
   511                   LowMap& low_map) {
       
   512       Node child = child_lists[node].first;
       
   513 
       
   514       if (child != INVALID) {
       
   515         if (low_map[child] < rorder) return true;
       
   516       }
       
   517 
       
   518       if (ancestor_map[node] < rorder) return true;
       
   519 
       
   520       return false;
       
   521     }
       
   522 
       
   523     bool pertinent(const Node& node, const EmbedArc& embed_arc,
       
   524                    const MergeRoots& merge_roots) {
       
   525       return !merge_roots[node].empty() || embed_arc[node];
       
   526     }
       
   527 
       
   528   };
       
   529 
       
   530   /// \ingroup planar
       
   531   ///
       
   532   /// \brief Planar embedding of an undirected simple graph
       
   533   ///
       
   534   /// This class implements the Boyer-Myrvold algorithm for planar
       
   535   /// embedding of an undirected graph. The planar embedding is an
       
   536   /// ordering of the outgoing edges of the nodes, which is a possible
       
   537   /// configuration to draw the graph in the plane. If there is not
       
   538   /// such ordering then the graph contains a \f$ K_5 \f$ (full graph
       
   539   /// with 5 nodes) or a \f$ K_{3,3} \f$ (complete bipartite graph on
       
   540   /// 3 ANode and 3 BNode) subdivision.
       
   541   ///
       
   542   /// The current implementation calculates either an embedding or a
       
   543   /// Kuratowski subdivision. The running time of the algorithm is 
       
   544   /// \f$ O(n) \f$.
       
   545   template <typename Graph>
       
   546   class PlanarEmbedding {
       
   547   private:
       
   548 
       
   549     TEMPLATE_GRAPH_TYPEDEFS(Graph);
       
   550 
       
   551     const Graph& _graph;
       
   552     typename Graph::template ArcMap<Arc> _embedding;
       
   553 
       
   554     typename Graph::template EdgeMap<bool> _kuratowski;
       
   555 
       
   556   private:
       
   557 
       
   558     typedef typename Graph::template NodeMap<Arc> PredMap;
       
   559 
       
   560     typedef typename Graph::template EdgeMap<bool> TreeMap;
       
   561 
       
   562     typedef typename Graph::template NodeMap<int> OrderMap;
       
   563     typedef std::vector<Node> OrderList;
       
   564 
       
   565     typedef typename Graph::template NodeMap<int> LowMap;
       
   566     typedef typename Graph::template NodeMap<int> AncestorMap;
       
   567 
       
   568     typedef _planarity_bits::NodeDataNode<Graph> NodeDataNode;
       
   569     typedef std::vector<NodeDataNode> NodeData;
       
   570 
       
   571     typedef _planarity_bits::ChildListNode<Graph> ChildListNode;
       
   572     typedef typename Graph::template NodeMap<ChildListNode> ChildLists;
       
   573 
       
   574     typedef typename Graph::template NodeMap<std::list<int> > MergeRoots;
       
   575 
       
   576     typedef typename Graph::template NodeMap<Arc> EmbedArc;
       
   577 
       
   578     typedef _planarity_bits::ArcListNode<Graph> ArcListNode;
       
   579     typedef typename Graph::template ArcMap<ArcListNode> ArcLists;
       
   580 
       
   581     typedef typename Graph::template NodeMap<bool> FlipMap;
       
   582 
       
   583     typedef typename Graph::template NodeMap<int> TypeMap;
       
   584 
       
   585     enum IsolatorNodeType {
       
   586       HIGHX = 6, LOWX = 7,
       
   587       HIGHY = 8, LOWY = 9,
       
   588       ROOT = 10, PERTINENT = 11,
       
   589       INTERNAL = 12
       
   590     };
       
   591 
       
   592   public:
       
   593 
       
   594     /// \brief The map for store of embedding
       
   595     typedef typename Graph::template ArcMap<Arc> EmbeddingMap;
       
   596 
       
   597     /// \brief Constructor
       
   598     ///
       
   599     /// \note The graph should be simple, i.e. parallel and loop arc
       
   600     /// free.
       
   601     PlanarEmbedding(const Graph& graph)
       
   602       : _graph(graph), _embedding(_graph), _kuratowski(graph, false) {}
       
   603 
       
   604     /// \brief Runs the algorithm.
       
   605     ///
       
   606     /// Runs the algorithm.
       
   607     /// \param kuratowski If the parameter is false, then the
       
   608     /// algorithm does not compute a Kuratowski subdivision.
       
   609     ///\return %True when the graph is planar.
       
   610     bool run(bool kuratowski = true) {
       
   611       typedef _planarity_bits::PlanarityVisitor<Graph> Visitor;
       
   612 
       
   613       PredMap pred_map(_graph, INVALID);
       
   614       TreeMap tree_map(_graph, false);
       
   615 
       
   616       OrderMap order_map(_graph, -1);
       
   617       OrderList order_list;
       
   618 
       
   619       AncestorMap ancestor_map(_graph, -1);
       
   620       LowMap low_map(_graph, -1);
       
   621 
       
   622       Visitor visitor(_graph, pred_map, tree_map,
       
   623                       order_map, order_list, ancestor_map, low_map);
       
   624       DfsVisit<Graph, Visitor> visit(_graph, visitor);
       
   625       visit.run();
       
   626 
       
   627       ChildLists child_lists(_graph);
       
   628       createChildLists(tree_map, order_map, low_map, child_lists);
       
   629 
       
   630       NodeData node_data(2 * order_list.size());
       
   631 
       
   632       EmbedArc embed_arc(_graph, INVALID);
       
   633 
       
   634       MergeRoots merge_roots(_graph);
       
   635 
       
   636       ArcLists arc_lists(_graph);
       
   637 
       
   638       FlipMap flip_map(_graph, false);
       
   639 
       
   640       for (int i = order_list.size() - 1; i >= 0; --i) {
       
   641 
       
   642         Node node = order_list[i];
       
   643 
       
   644         node_data[i].first = INVALID;
       
   645 
       
   646         Node source = node;
       
   647         for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
   648           Node target = _graph.target(e);
       
   649 
       
   650           if (order_map[source] < order_map[target] && tree_map[e]) {
       
   651             initFace(target, arc_lists, node_data,
       
   652                      pred_map, order_map, order_list);
       
   653           }
       
   654         }
       
   655 
       
   656         for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
   657           Node target = _graph.target(e);
       
   658 
       
   659           if (order_map[source] < order_map[target] && !tree_map[e]) {
       
   660             embed_arc[target] = e;
       
   661             walkUp(target, source, i, pred_map, low_map,
       
   662                    order_map, order_list, node_data, merge_roots);
       
   663           }
       
   664         }
       
   665 
       
   666         for (typename MergeRoots::Value::iterator it =
       
   667                merge_roots[node].begin(); it != merge_roots[node].end(); ++it) {
       
   668           int rn = *it;
       
   669           walkDown(rn, i, node_data, arc_lists, flip_map, order_list,
       
   670                    child_lists, ancestor_map, low_map, embed_arc, merge_roots);
       
   671         }
       
   672         merge_roots[node].clear();
       
   673 
       
   674         for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
   675           Node target = _graph.target(e);
       
   676 
       
   677           if (order_map[source] < order_map[target] && !tree_map[e]) {
       
   678             if (embed_arc[target] != INVALID) {
       
   679               if (kuratowski) {
       
   680                 isolateKuratowski(e, node_data, arc_lists, flip_map,
       
   681                                   order_map, order_list, pred_map, child_lists,
       
   682                                   ancestor_map, low_map,
       
   683                                   embed_arc, merge_roots);
       
   684               }
       
   685               return false;
       
   686             }
       
   687           }
       
   688         }
       
   689       }
       
   690 
       
   691       for (int i = 0; i < int(order_list.size()); ++i) {
       
   692 
       
   693         mergeRemainingFaces(order_list[i], node_data, order_list, order_map,
       
   694                             child_lists, arc_lists);
       
   695         storeEmbedding(order_list[i], node_data, order_map, pred_map,
       
   696                        arc_lists, flip_map);
       
   697       }
       
   698 
       
   699       return true;
       
   700     }
       
   701 
       
   702     /// \brief Gives back the successor of an arc
       
   703     ///
       
   704     /// Gives back the successor of an arc. This function makes
       
   705     /// possible to query the cyclic order of the outgoing arcs from
       
   706     /// a node.
       
   707     Arc next(const Arc& arc) const {
       
   708       return _embedding[arc];
       
   709     }
       
   710 
       
   711     /// \brief Gives back the calculated embedding map
       
   712     ///
       
   713     /// The returned map contains the successor of each arc in the
       
   714     /// graph.
       
   715     const EmbeddingMap& embedding() const {
       
   716       return _embedding;
       
   717     }
       
   718 
       
   719     /// \brief Gives back true if the undirected arc is in the
       
   720     /// kuratowski subdivision
       
   721     ///
       
   722     /// Gives back true if the undirected arc is in the kuratowski
       
   723     /// subdivision
       
   724     /// \note The \c run() had to be called with true value.
       
   725     bool kuratowski(const Edge& edge) {
       
   726       return _kuratowski[edge];
       
   727     }
       
   728 
       
   729   private:
       
   730 
       
   731     void createChildLists(const TreeMap& tree_map, const OrderMap& order_map,
       
   732                           const LowMap& low_map, ChildLists& child_lists) {
       
   733 
       
   734       for (NodeIt n(_graph); n != INVALID; ++n) {
       
   735         Node source = n;
       
   736 
       
   737         std::vector<Node> targets;
       
   738         for (OutArcIt e(_graph, n); e != INVALID; ++e) {
       
   739           Node target = _graph.target(e);
       
   740 
       
   741           if (order_map[source] < order_map[target] && tree_map[e]) {
       
   742             targets.push_back(target);
       
   743           }
       
   744         }
       
   745 
       
   746         if (targets.size() == 0) {
       
   747           child_lists[source].first = INVALID;
       
   748         } else if (targets.size() == 1) {
       
   749           child_lists[source].first = targets[0];
       
   750           child_lists[targets[0]].prev = INVALID;
       
   751           child_lists[targets[0]].next = INVALID;
       
   752         } else {
       
   753           radixSort(targets.begin(), targets.end(), mapToFunctor(low_map));
       
   754           for (int i = 1; i < int(targets.size()); ++i) {
       
   755             child_lists[targets[i]].prev = targets[i - 1];
       
   756             child_lists[targets[i - 1]].next = targets[i];
       
   757           }
       
   758           child_lists[targets.back()].next = INVALID;
       
   759           child_lists[targets.front()].prev = INVALID;
       
   760           child_lists[source].first = targets.front();
       
   761         }
       
   762       }
       
   763     }
       
   764 
       
   765     void walkUp(const Node& node, Node root, int rorder,
       
   766                 const PredMap& pred_map, const LowMap& low_map,
       
   767                 const OrderMap& order_map, const OrderList& order_list,
       
   768                 NodeData& node_data, MergeRoots& merge_roots) {
       
   769 
       
   770       int na, nb;
       
   771       bool da, db;
       
   772 
       
   773       na = nb = order_map[node];
       
   774       da = true; db = false;
       
   775 
       
   776       while (true) {
       
   777 
       
   778         if (node_data[na].visited == rorder) break;
       
   779         if (node_data[nb].visited == rorder) break;
       
   780 
       
   781         node_data[na].visited = rorder;
       
   782         node_data[nb].visited = rorder;
       
   783 
       
   784         int rn = -1;
       
   785 
       
   786         if (na >= int(order_list.size())) {
       
   787           rn = na;
       
   788         } else if (nb >= int(order_list.size())) {
       
   789           rn = nb;
       
   790         }
       
   791 
       
   792         if (rn == -1) {
       
   793           int nn;
       
   794 
       
   795           nn = da ? node_data[na].prev : node_data[na].next;
       
   796           da = node_data[nn].prev != na;
       
   797           na = nn;
       
   798 
       
   799           nn = db ? node_data[nb].prev : node_data[nb].next;
       
   800           db = node_data[nn].prev != nb;
       
   801           nb = nn;
       
   802 
       
   803         } else {
       
   804 
       
   805           Node rep = order_list[rn - order_list.size()];
       
   806           Node parent = _graph.source(pred_map[rep]);
       
   807 
       
   808           if (low_map[rep] < rorder) {
       
   809             merge_roots[parent].push_back(rn);
       
   810           } else {
       
   811             merge_roots[parent].push_front(rn);
       
   812           }
       
   813 
       
   814           if (parent != root) {
       
   815             na = nb = order_map[parent];
       
   816             da = true; db = false;
       
   817           } else {
       
   818             break;
       
   819           }
       
   820         }
       
   821       }
       
   822     }
       
   823 
       
   824     void walkDown(int rn, int rorder, NodeData& node_data,
       
   825                   ArcLists& arc_lists, FlipMap& flip_map,
       
   826                   OrderList& order_list, ChildLists& child_lists,
       
   827                   AncestorMap& ancestor_map, LowMap& low_map,
       
   828                   EmbedArc& embed_arc, MergeRoots& merge_roots) {
       
   829 
       
   830       std::vector<std::pair<int, bool> > merge_stack;
       
   831 
       
   832       for (int di = 0; di < 2; ++di) {
       
   833         bool rd = di == 0;
       
   834         int pn = rn;
       
   835         int n = rd ? node_data[rn].next : node_data[rn].prev;
       
   836 
       
   837         while (n != rn) {
       
   838 
       
   839           Node node = order_list[n];
       
   840 
       
   841           if (embed_arc[node] != INVALID) {
       
   842 
       
   843             // Merging components on the critical path
       
   844             while (!merge_stack.empty()) {
       
   845 
       
   846               // Component root
       
   847               int cn = merge_stack.back().first;
       
   848               bool cd = merge_stack.back().second;
       
   849               merge_stack.pop_back();
       
   850 
       
   851               // Parent of component
       
   852               int dn = merge_stack.back().first;
       
   853               bool dd = merge_stack.back().second;
       
   854               merge_stack.pop_back();
       
   855 
       
   856               Node parent = order_list[dn];
       
   857 
       
   858               // Erasing from merge_roots
       
   859               merge_roots[parent].pop_front();
       
   860 
       
   861               Node child = order_list[cn - order_list.size()];
       
   862 
       
   863               // Erasing from child_lists
       
   864               if (child_lists[child].prev != INVALID) {
       
   865                 child_lists[child_lists[child].prev].next =
       
   866                   child_lists[child].next;
       
   867               } else {
       
   868                 child_lists[parent].first = child_lists[child].next;
       
   869               }
       
   870 
       
   871               if (child_lists[child].next != INVALID) {
       
   872                 child_lists[child_lists[child].next].prev =
       
   873                   child_lists[child].prev;
       
   874               }
       
   875 
       
   876               // Merging arcs + flipping
       
   877               Arc de = node_data[dn].first;
       
   878               Arc ce = node_data[cn].first;
       
   879 
       
   880               flip_map[order_list[cn - order_list.size()]] = cd != dd;
       
   881               if (cd != dd) {
       
   882                 std::swap(arc_lists[ce].prev, arc_lists[ce].next);
       
   883                 ce = arc_lists[ce].prev;
       
   884                 std::swap(arc_lists[ce].prev, arc_lists[ce].next);
       
   885               }
       
   886 
       
   887               {
       
   888                 Arc dne = arc_lists[de].next;
       
   889                 Arc cne = arc_lists[ce].next;
       
   890 
       
   891                 arc_lists[de].next = cne;
       
   892                 arc_lists[ce].next = dne;
       
   893 
       
   894                 arc_lists[dne].prev = ce;
       
   895                 arc_lists[cne].prev = de;
       
   896               }
       
   897 
       
   898               if (dd) {
       
   899                 node_data[dn].first = ce;
       
   900               }
       
   901 
       
   902               // Merging external faces
       
   903               {
       
   904                 int en = cn;
       
   905                 cn = cd ? node_data[cn].prev : node_data[cn].next;
       
   906                 cd = node_data[cn].next == en;
       
   907 
       
   908                  if (node_data[cn].prev == node_data[cn].next &&
       
   909                     node_data[cn].inverted) {
       
   910                    cd = !cd;
       
   911                  }
       
   912               }
       
   913 
       
   914               if (cd) node_data[cn].next = dn; else node_data[cn].prev = dn;
       
   915               if (dd) node_data[dn].prev = cn; else node_data[dn].next = cn;
       
   916 
       
   917             }
       
   918 
       
   919             bool d = pn == node_data[n].prev;
       
   920 
       
   921             if (node_data[n].prev == node_data[n].next &&
       
   922                 node_data[n].inverted) {
       
   923               d = !d;
       
   924             }
       
   925 
       
   926             // Add new arc
       
   927             {
       
   928               Arc arc = embed_arc[node];
       
   929               Arc re = node_data[rn].first;
       
   930 
       
   931               arc_lists[arc_lists[re].next].prev = arc;
       
   932               arc_lists[arc].next = arc_lists[re].next;
       
   933               arc_lists[arc].prev = re;
       
   934               arc_lists[re].next = arc;
       
   935 
       
   936               if (!rd) {
       
   937                 node_data[rn].first = arc;
       
   938               }
       
   939 
       
   940               Arc rev = _graph.oppositeArc(arc);
       
   941               Arc e = node_data[n].first;
       
   942 
       
   943               arc_lists[arc_lists[e].next].prev = rev;
       
   944               arc_lists[rev].next = arc_lists[e].next;
       
   945               arc_lists[rev].prev = e;
       
   946               arc_lists[e].next = rev;
       
   947 
       
   948               if (d) {
       
   949                 node_data[n].first = rev;
       
   950               }
       
   951 
       
   952             }
       
   953 
       
   954             // Embedding arc into external face
       
   955             if (rd) node_data[rn].next = n; else node_data[rn].prev = n;
       
   956             if (d) node_data[n].prev = rn; else node_data[n].next = rn;
       
   957             pn = rn;
       
   958 
       
   959             embed_arc[order_list[n]] = INVALID;
       
   960           }
       
   961 
       
   962           if (!merge_roots[node].empty()) {
       
   963 
       
   964             bool d = pn == node_data[n].prev;
       
   965             if (node_data[n].prev == node_data[n].next &&
       
   966                 node_data[n].inverted) {
       
   967               d = !d;
       
   968             }
       
   969 
       
   970             merge_stack.push_back(std::make_pair(n, d));
       
   971 
       
   972             int rn = merge_roots[node].front();
       
   973 
       
   974             int xn = node_data[rn].next;
       
   975             Node xnode = order_list[xn];
       
   976 
       
   977             int yn = node_data[rn].prev;
       
   978             Node ynode = order_list[yn];
       
   979 
       
   980             bool rd;
       
   981             if (!external(xnode, rorder, child_lists, ancestor_map, low_map)) {
       
   982               rd = true;
       
   983             } else if (!external(ynode, rorder, child_lists,
       
   984                                  ancestor_map, low_map)) {
       
   985               rd = false;
       
   986             } else if (pertinent(xnode, embed_arc, merge_roots)) {
       
   987               rd = true;
       
   988             } else {
       
   989               rd = false;
       
   990             }
       
   991 
       
   992             merge_stack.push_back(std::make_pair(rn, rd));
       
   993 
       
   994             pn = rn;
       
   995             n = rd ? xn : yn;
       
   996 
       
   997           } else if (!external(node, rorder, child_lists,
       
   998                                ancestor_map, low_map)) {
       
   999             int nn = (node_data[n].next != pn ?
       
  1000                       node_data[n].next : node_data[n].prev);
       
  1001 
       
  1002             bool nd = n == node_data[nn].prev;
       
  1003 
       
  1004             if (nd) node_data[nn].prev = pn;
       
  1005             else node_data[nn].next = pn;
       
  1006 
       
  1007             if (n == node_data[pn].prev) node_data[pn].prev = nn;
       
  1008             else node_data[pn].next = nn;
       
  1009 
       
  1010             node_data[nn].inverted =
       
  1011               (node_data[nn].prev == node_data[nn].next && nd != rd);
       
  1012 
       
  1013             n = nn;
       
  1014           }
       
  1015           else break;
       
  1016 
       
  1017         }
       
  1018 
       
  1019         if (!merge_stack.empty() || n == rn) {
       
  1020           break;
       
  1021         }
       
  1022       }
       
  1023     }
       
  1024 
       
  1025     void initFace(const Node& node, ArcLists& arc_lists,
       
  1026                   NodeData& node_data, const PredMap& pred_map,
       
  1027                   const OrderMap& order_map, const OrderList& order_list) {
       
  1028       int n = order_map[node];
       
  1029       int rn = n + order_list.size();
       
  1030 
       
  1031       node_data[n].next = node_data[n].prev = rn;
       
  1032       node_data[rn].next = node_data[rn].prev = n;
       
  1033 
       
  1034       node_data[n].visited = order_list.size();
       
  1035       node_data[rn].visited = order_list.size();
       
  1036 
       
  1037       node_data[n].inverted = false;
       
  1038       node_data[rn].inverted = false;
       
  1039 
       
  1040       Arc arc = pred_map[node];
       
  1041       Arc rev = _graph.oppositeArc(arc);
       
  1042 
       
  1043       node_data[rn].first = arc;
       
  1044       node_data[n].first = rev;
       
  1045 
       
  1046       arc_lists[arc].prev = arc;
       
  1047       arc_lists[arc].next = arc;
       
  1048 
       
  1049       arc_lists[rev].prev = rev;
       
  1050       arc_lists[rev].next = rev;
       
  1051 
       
  1052     }
       
  1053 
       
  1054     void mergeRemainingFaces(const Node& node, NodeData& node_data,
       
  1055                              OrderList& order_list, OrderMap& order_map,
       
  1056                              ChildLists& child_lists, ArcLists& arc_lists) {
       
  1057       while (child_lists[node].first != INVALID) {
       
  1058         int dd = order_map[node];
       
  1059         Node child = child_lists[node].first;
       
  1060         int cd = order_map[child] + order_list.size();
       
  1061         child_lists[node].first = child_lists[child].next;
       
  1062 
       
  1063         Arc de = node_data[dd].first;
       
  1064         Arc ce = node_data[cd].first;
       
  1065 
       
  1066         if (de != INVALID) {
       
  1067           Arc dne = arc_lists[de].next;
       
  1068           Arc cne = arc_lists[ce].next;
       
  1069 
       
  1070           arc_lists[de].next = cne;
       
  1071           arc_lists[ce].next = dne;
       
  1072 
       
  1073           arc_lists[dne].prev = ce;
       
  1074           arc_lists[cne].prev = de;
       
  1075         }
       
  1076 
       
  1077         node_data[dd].first = ce;
       
  1078 
       
  1079       }
       
  1080     }
       
  1081 
       
  1082     void storeEmbedding(const Node& node, NodeData& node_data,
       
  1083                         OrderMap& order_map, PredMap& pred_map,
       
  1084                         ArcLists& arc_lists, FlipMap& flip_map) {
       
  1085 
       
  1086       if (node_data[order_map[node]].first == INVALID) return;
       
  1087 
       
  1088       if (pred_map[node] != INVALID) {
       
  1089         Node source = _graph.source(pred_map[node]);
       
  1090         flip_map[node] = flip_map[node] != flip_map[source];
       
  1091       }
       
  1092 
       
  1093       Arc first = node_data[order_map[node]].first;
       
  1094       Arc prev = first;
       
  1095 
       
  1096       Arc arc = flip_map[node] ?
       
  1097         arc_lists[prev].prev : arc_lists[prev].next;
       
  1098 
       
  1099       _embedding[prev] = arc;
       
  1100 
       
  1101       while (arc != first) {
       
  1102         Arc next = arc_lists[arc].prev == prev ?
       
  1103           arc_lists[arc].next : arc_lists[arc].prev;
       
  1104         prev = arc; arc = next;
       
  1105         _embedding[prev] = arc;
       
  1106       }
       
  1107     }
       
  1108 
       
  1109 
       
  1110     bool external(const Node& node, int rorder,
       
  1111                   ChildLists& child_lists, AncestorMap& ancestor_map,
       
  1112                   LowMap& low_map) {
       
  1113       Node child = child_lists[node].first;
       
  1114 
       
  1115       if (child != INVALID) {
       
  1116         if (low_map[child] < rorder) return true;
       
  1117       }
       
  1118 
       
  1119       if (ancestor_map[node] < rorder) return true;
       
  1120 
       
  1121       return false;
       
  1122     }
       
  1123 
       
  1124     bool pertinent(const Node& node, const EmbedArc& embed_arc,
       
  1125                    const MergeRoots& merge_roots) {
       
  1126       return !merge_roots[node].empty() || embed_arc[node] != INVALID;
       
  1127     }
       
  1128 
       
  1129     int lowPoint(const Node& node, OrderMap& order_map, ChildLists& child_lists,
       
  1130                  AncestorMap& ancestor_map, LowMap& low_map) {
       
  1131       int low_point;
       
  1132 
       
  1133       Node child = child_lists[node].first;
       
  1134 
       
  1135       if (child != INVALID) {
       
  1136         low_point = low_map[child];
       
  1137       } else {
       
  1138         low_point = order_map[node];
       
  1139       }
       
  1140 
       
  1141       if (low_point > ancestor_map[node]) {
       
  1142         low_point = ancestor_map[node];
       
  1143       }
       
  1144 
       
  1145       return low_point;
       
  1146     }
       
  1147 
       
  1148     int findComponentRoot(Node root, Node node, ChildLists& child_lists,
       
  1149                           OrderMap& order_map, OrderList& order_list) {
       
  1150 
       
  1151       int order = order_map[root];
       
  1152       int norder = order_map[node];
       
  1153 
       
  1154       Node child = child_lists[root].first;
       
  1155       while (child != INVALID) {
       
  1156         int corder = order_map[child];
       
  1157         if (corder > order && corder < norder) {
       
  1158           order = corder;
       
  1159         }
       
  1160         child = child_lists[child].next;
       
  1161       }
       
  1162       return order + order_list.size();
       
  1163     }
       
  1164 
       
  1165     Node findPertinent(Node node, OrderMap& order_map, NodeData& node_data,
       
  1166                        EmbedArc& embed_arc, MergeRoots& merge_roots) {
       
  1167       Node wnode =_graph.target(node_data[order_map[node]].first);
       
  1168       while (!pertinent(wnode, embed_arc, merge_roots)) {
       
  1169         wnode = _graph.target(node_data[order_map[wnode]].first);
       
  1170       }
       
  1171       return wnode;
       
  1172     }
       
  1173 
       
  1174 
       
  1175     Node findExternal(Node node, int rorder, OrderMap& order_map,
       
  1176                       ChildLists& child_lists, AncestorMap& ancestor_map,
       
  1177                       LowMap& low_map, NodeData& node_data) {
       
  1178       Node wnode =_graph.target(node_data[order_map[node]].first);
       
  1179       while (!external(wnode, rorder, child_lists, ancestor_map, low_map)) {
       
  1180         wnode = _graph.target(node_data[order_map[wnode]].first);
       
  1181       }
       
  1182       return wnode;
       
  1183     }
       
  1184 
       
  1185     void markCommonPath(Node node, int rorder, Node& wnode, Node& znode,
       
  1186                         OrderList& order_list, OrderMap& order_map,
       
  1187                         NodeData& node_data, ArcLists& arc_lists,
       
  1188                         EmbedArc& embed_arc, MergeRoots& merge_roots,
       
  1189                         ChildLists& child_lists, AncestorMap& ancestor_map,
       
  1190                         LowMap& low_map) {
       
  1191 
       
  1192       Node cnode = node;
       
  1193       Node pred = INVALID;
       
  1194 
       
  1195       while (true) {
       
  1196 
       
  1197         bool pert = pertinent(cnode, embed_arc, merge_roots);
       
  1198         bool ext = external(cnode, rorder, child_lists, ancestor_map, low_map);
       
  1199 
       
  1200         if (pert && ext) {
       
  1201           if (!merge_roots[cnode].empty()) {
       
  1202             int cn = merge_roots[cnode].back();
       
  1203 
       
  1204             if (low_map[order_list[cn - order_list.size()]] < rorder) {
       
  1205               Arc arc = node_data[cn].first;
       
  1206               _kuratowski.set(arc, true);
       
  1207 
       
  1208               pred = cnode;
       
  1209               cnode = _graph.target(arc);
       
  1210 
       
  1211               continue;
       
  1212             }
       
  1213           }
       
  1214           wnode = znode = cnode;
       
  1215           return;
       
  1216 
       
  1217         } else if (pert) {
       
  1218           wnode = cnode;
       
  1219 
       
  1220           while (!external(cnode, rorder, child_lists, ancestor_map, low_map)) {
       
  1221             Arc arc = node_data[order_map[cnode]].first;
       
  1222 
       
  1223             if (_graph.target(arc) == pred) {
       
  1224               arc = arc_lists[arc].next;
       
  1225             }
       
  1226             _kuratowski.set(arc, true);
       
  1227 
       
  1228             Node next = _graph.target(arc);
       
  1229             pred = cnode; cnode = next;
       
  1230           }
       
  1231 
       
  1232           znode = cnode;
       
  1233           return;
       
  1234 
       
  1235         } else if (ext) {
       
  1236           znode = cnode;
       
  1237 
       
  1238           while (!pertinent(cnode, embed_arc, merge_roots)) {
       
  1239             Arc arc = node_data[order_map[cnode]].first;
       
  1240 
       
  1241             if (_graph.target(arc) == pred) {
       
  1242               arc = arc_lists[arc].next;
       
  1243             }
       
  1244             _kuratowski.set(arc, true);
       
  1245 
       
  1246             Node next = _graph.target(arc);
       
  1247             pred = cnode; cnode = next;
       
  1248           }
       
  1249 
       
  1250           wnode = cnode;
       
  1251           return;
       
  1252 
       
  1253         } else {
       
  1254           Arc arc = node_data[order_map[cnode]].first;
       
  1255 
       
  1256           if (_graph.target(arc) == pred) {
       
  1257             arc = arc_lists[arc].next;
       
  1258           }
       
  1259           _kuratowski.set(arc, true);
       
  1260 
       
  1261           Node next = _graph.target(arc);
       
  1262           pred = cnode; cnode = next;
       
  1263         }
       
  1264 
       
  1265       }
       
  1266 
       
  1267     }
       
  1268 
       
  1269     void orientComponent(Node root, int rn, OrderMap& order_map,
       
  1270                          PredMap& pred_map, NodeData& node_data,
       
  1271                          ArcLists& arc_lists, FlipMap& flip_map,
       
  1272                          TypeMap& type_map) {
       
  1273       node_data[order_map[root]].first = node_data[rn].first;
       
  1274       type_map[root] = 1;
       
  1275 
       
  1276       std::vector<Node> st, qu;
       
  1277 
       
  1278       st.push_back(root);
       
  1279       while (!st.empty()) {
       
  1280         Node node = st.back();
       
  1281         st.pop_back();
       
  1282         qu.push_back(node);
       
  1283 
       
  1284         Arc arc = node_data[order_map[node]].first;
       
  1285 
       
  1286         if (type_map[_graph.target(arc)] == 0) {
       
  1287           st.push_back(_graph.target(arc));
       
  1288           type_map[_graph.target(arc)] = 1;
       
  1289         }
       
  1290 
       
  1291         Arc last = arc, pred = arc;
       
  1292         arc = arc_lists[arc].next;
       
  1293         while (arc != last) {
       
  1294 
       
  1295           if (type_map[_graph.target(arc)] == 0) {
       
  1296             st.push_back(_graph.target(arc));
       
  1297             type_map[_graph.target(arc)] = 1;
       
  1298           }
       
  1299 
       
  1300           Arc next = arc_lists[arc].next != pred ?
       
  1301             arc_lists[arc].next : arc_lists[arc].prev;
       
  1302           pred = arc; arc = next;
       
  1303         }
       
  1304 
       
  1305       }
       
  1306 
       
  1307       type_map[root] = 2;
       
  1308       flip_map[root] = false;
       
  1309 
       
  1310       for (int i = 1; i < int(qu.size()); ++i) {
       
  1311 
       
  1312         Node node = qu[i];
       
  1313 
       
  1314         while (type_map[node] != 2) {
       
  1315           st.push_back(node);
       
  1316           type_map[node] = 2;
       
  1317           node = _graph.source(pred_map[node]);
       
  1318         }
       
  1319 
       
  1320         bool flip = flip_map[node];
       
  1321 
       
  1322         while (!st.empty()) {
       
  1323           node = st.back();
       
  1324           st.pop_back();
       
  1325 
       
  1326           flip_map[node] = flip != flip_map[node];
       
  1327           flip = flip_map[node];
       
  1328 
       
  1329           if (flip) {
       
  1330             Arc arc = node_data[order_map[node]].first;
       
  1331             std::swap(arc_lists[arc].prev, arc_lists[arc].next);
       
  1332             arc = arc_lists[arc].prev;
       
  1333             std::swap(arc_lists[arc].prev, arc_lists[arc].next);
       
  1334             node_data[order_map[node]].first = arc;
       
  1335           }
       
  1336         }
       
  1337       }
       
  1338 
       
  1339       for (int i = 0; i < int(qu.size()); ++i) {
       
  1340 
       
  1341         Arc arc = node_data[order_map[qu[i]]].first;
       
  1342         Arc last = arc, pred = arc;
       
  1343 
       
  1344         arc = arc_lists[arc].next;
       
  1345         while (arc != last) {
       
  1346 
       
  1347           if (arc_lists[arc].next == pred) {
       
  1348             std::swap(arc_lists[arc].next, arc_lists[arc].prev);
       
  1349           }
       
  1350           pred = arc; arc = arc_lists[arc].next;
       
  1351         }
       
  1352 
       
  1353       }
       
  1354     }
       
  1355 
       
  1356     void setFaceFlags(Node root, Node wnode, Node ynode, Node xnode,
       
  1357                       OrderMap& order_map, NodeData& node_data,
       
  1358                       TypeMap& type_map) {
       
  1359       Node node = _graph.target(node_data[order_map[root]].first);
       
  1360 
       
  1361       while (node != ynode) {
       
  1362         type_map[node] = HIGHY;
       
  1363         node = _graph.target(node_data[order_map[node]].first);
       
  1364       }
       
  1365 
       
  1366       while (node != wnode) {
       
  1367         type_map[node] = LOWY;
       
  1368         node = _graph.target(node_data[order_map[node]].first);
       
  1369       }
       
  1370 
       
  1371       node = _graph.target(node_data[order_map[wnode]].first);
       
  1372 
       
  1373       while (node != xnode) {
       
  1374         type_map[node] = LOWX;
       
  1375         node = _graph.target(node_data[order_map[node]].first);
       
  1376       }
       
  1377       type_map[node] = LOWX;
       
  1378 
       
  1379       node = _graph.target(node_data[order_map[xnode]].first);
       
  1380       while (node != root) {
       
  1381         type_map[node] = HIGHX;
       
  1382         node = _graph.target(node_data[order_map[node]].first);
       
  1383       }
       
  1384 
       
  1385       type_map[wnode] = PERTINENT;
       
  1386       type_map[root] = ROOT;
       
  1387     }
       
  1388 
       
  1389     void findInternalPath(std::vector<Arc>& ipath,
       
  1390                           Node wnode, Node root, TypeMap& type_map,
       
  1391                           OrderMap& order_map, NodeData& node_data,
       
  1392                           ArcLists& arc_lists) {
       
  1393       std::vector<Arc> st;
       
  1394 
       
  1395       Node node = wnode;
       
  1396 
       
  1397       while (node != root) {
       
  1398         Arc arc = arc_lists[node_data[order_map[node]].first].next;
       
  1399         st.push_back(arc);
       
  1400         node = _graph.target(arc);
       
  1401       }
       
  1402 
       
  1403       while (true) {
       
  1404         Arc arc = st.back();
       
  1405         if (type_map[_graph.target(arc)] == LOWX ||
       
  1406             type_map[_graph.target(arc)] == HIGHX) {
       
  1407           break;
       
  1408         }
       
  1409         if (type_map[_graph.target(arc)] == 2) {
       
  1410           type_map[_graph.target(arc)] = 3;
       
  1411 
       
  1412           arc = arc_lists[_graph.oppositeArc(arc)].next;
       
  1413           st.push_back(arc);
       
  1414         } else {
       
  1415           st.pop_back();
       
  1416           arc = arc_lists[arc].next;
       
  1417 
       
  1418           while (_graph.oppositeArc(arc) == st.back()) {
       
  1419             arc = st.back();
       
  1420             st.pop_back();
       
  1421             arc = arc_lists[arc].next;
       
  1422           }
       
  1423           st.push_back(arc);
       
  1424         }
       
  1425       }
       
  1426 
       
  1427       for (int i = 0; i < int(st.size()); ++i) {
       
  1428         if (type_map[_graph.target(st[i])] != LOWY &&
       
  1429             type_map[_graph.target(st[i])] != HIGHY) {
       
  1430           for (; i < int(st.size()); ++i) {
       
  1431             ipath.push_back(st[i]);
       
  1432           }
       
  1433         }
       
  1434       }
       
  1435     }
       
  1436 
       
  1437     void setInternalFlags(std::vector<Arc>& ipath, TypeMap& type_map) {
       
  1438       for (int i = 1; i < int(ipath.size()); ++i) {
       
  1439         type_map[_graph.source(ipath[i])] = INTERNAL;
       
  1440       }
       
  1441     }
       
  1442 
       
  1443     void findPilePath(std::vector<Arc>& ppath,
       
  1444                       Node root, TypeMap& type_map, OrderMap& order_map,
       
  1445                       NodeData& node_data, ArcLists& arc_lists) {
       
  1446       std::vector<Arc> st;
       
  1447 
       
  1448       st.push_back(_graph.oppositeArc(node_data[order_map[root]].first));
       
  1449       st.push_back(node_data[order_map[root]].first);
       
  1450 
       
  1451       while (st.size() > 1) {
       
  1452         Arc arc = st.back();
       
  1453         if (type_map[_graph.target(arc)] == INTERNAL) {
       
  1454           break;
       
  1455         }
       
  1456         if (type_map[_graph.target(arc)] == 3) {
       
  1457           type_map[_graph.target(arc)] = 4;
       
  1458 
       
  1459           arc = arc_lists[_graph.oppositeArc(arc)].next;
       
  1460           st.push_back(arc);
       
  1461         } else {
       
  1462           st.pop_back();
       
  1463           arc = arc_lists[arc].next;
       
  1464 
       
  1465           while (!st.empty() && _graph.oppositeArc(arc) == st.back()) {
       
  1466             arc = st.back();
       
  1467             st.pop_back();
       
  1468             arc = arc_lists[arc].next;
       
  1469           }
       
  1470           st.push_back(arc);
       
  1471         }
       
  1472       }
       
  1473 
       
  1474       for (int i = 1; i < int(st.size()); ++i) {
       
  1475         ppath.push_back(st[i]);
       
  1476       }
       
  1477     }
       
  1478 
       
  1479 
       
  1480     int markExternalPath(Node node, OrderMap& order_map,
       
  1481                          ChildLists& child_lists, PredMap& pred_map,
       
  1482                          AncestorMap& ancestor_map, LowMap& low_map) {
       
  1483       int lp = lowPoint(node, order_map, child_lists,
       
  1484                         ancestor_map, low_map);
       
  1485 
       
  1486       if (ancestor_map[node] != lp) {
       
  1487         node = child_lists[node].first;
       
  1488         _kuratowski[pred_map[node]] = true;
       
  1489 
       
  1490         while (ancestor_map[node] != lp) {
       
  1491           for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
  1492             Node tnode = _graph.target(e);
       
  1493             if (order_map[tnode] > order_map[node] && low_map[tnode] == lp) {
       
  1494               node = tnode;
       
  1495               _kuratowski[e] = true;
       
  1496               break;
       
  1497             }
       
  1498           }
       
  1499         }
       
  1500       }
       
  1501 
       
  1502       for (OutArcIt e(_graph, node); e != INVALID; ++e) {
       
  1503         if (order_map[_graph.target(e)] == lp) {
       
  1504           _kuratowski[e] = true;
       
  1505           break;
       
  1506         }
       
  1507       }
       
  1508 
       
  1509       return lp;
       
  1510     }
       
  1511 
       
  1512     void markPertinentPath(Node node, OrderMap& order_map,
       
  1513                            NodeData& node_data, ArcLists& arc_lists,
       
  1514                            EmbedArc& embed_arc, MergeRoots& merge_roots) {
       
  1515       while (embed_arc[node] == INVALID) {
       
  1516         int n = merge_roots[node].front();
       
  1517         Arc arc = node_data[n].first;
       
  1518 
       
  1519         _kuratowski.set(arc, true);
       
  1520 
       
  1521         Node pred = node;
       
  1522         node = _graph.target(arc);
       
  1523         while (!pertinent(node, embed_arc, merge_roots)) {
       
  1524           arc = node_data[order_map[node]].first;
       
  1525           if (_graph.target(arc) == pred) {
       
  1526             arc = arc_lists[arc].next;
       
  1527           }
       
  1528           _kuratowski.set(arc, true);
       
  1529           pred = node;
       
  1530           node = _graph.target(arc);
       
  1531         }
       
  1532       }
       
  1533       _kuratowski.set(embed_arc[node], true);
       
  1534     }
       
  1535 
       
  1536     void markPredPath(Node node, Node snode, PredMap& pred_map) {
       
  1537       while (node != snode) {
       
  1538         _kuratowski.set(pred_map[node], true);
       
  1539         node = _graph.source(pred_map[node]);
       
  1540       }
       
  1541     }
       
  1542 
       
  1543     void markFacePath(Node ynode, Node xnode,
       
  1544                       OrderMap& order_map, NodeData& node_data) {
       
  1545       Arc arc = node_data[order_map[ynode]].first;
       
  1546       Node node = _graph.target(arc);
       
  1547       _kuratowski.set(arc, true);
       
  1548 
       
  1549       while (node != xnode) {
       
  1550         arc = node_data[order_map[node]].first;
       
  1551         _kuratowski.set(arc, true);
       
  1552         node = _graph.target(arc);
       
  1553       }
       
  1554     }
       
  1555 
       
  1556     void markInternalPath(std::vector<Arc>& path) {
       
  1557       for (int i = 0; i < int(path.size()); ++i) {
       
  1558         _kuratowski.set(path[i], true);
       
  1559       }
       
  1560     }
       
  1561 
       
  1562     void markPilePath(std::vector<Arc>& path) {
       
  1563       for (int i = 0; i < int(path.size()); ++i) {
       
  1564         _kuratowski.set(path[i], true);
       
  1565       }
       
  1566     }
       
  1567 
       
  1568     void isolateKuratowski(Arc arc, NodeData& node_data,
       
  1569                            ArcLists& arc_lists, FlipMap& flip_map,
       
  1570                            OrderMap& order_map, OrderList& order_list,
       
  1571                            PredMap& pred_map, ChildLists& child_lists,
       
  1572                            AncestorMap& ancestor_map, LowMap& low_map,
       
  1573                            EmbedArc& embed_arc, MergeRoots& merge_roots) {
       
  1574 
       
  1575       Node root = _graph.source(arc);
       
  1576       Node enode = _graph.target(arc);
       
  1577 
       
  1578       int rorder = order_map[root];
       
  1579 
       
  1580       TypeMap type_map(_graph, 0);
       
  1581 
       
  1582       int rn = findComponentRoot(root, enode, child_lists,
       
  1583                                  order_map, order_list);
       
  1584 
       
  1585       Node xnode = order_list[node_data[rn].next];
       
  1586       Node ynode = order_list[node_data[rn].prev];
       
  1587 
       
  1588       // Minor-A
       
  1589       {
       
  1590         while (!merge_roots[xnode].empty() || !merge_roots[ynode].empty()) {
       
  1591 
       
  1592           if (!merge_roots[xnode].empty()) {
       
  1593             root = xnode;
       
  1594             rn = merge_roots[xnode].front();
       
  1595           } else {
       
  1596             root = ynode;
       
  1597             rn = merge_roots[ynode].front();
       
  1598           }
       
  1599 
       
  1600           xnode = order_list[node_data[rn].next];
       
  1601           ynode = order_list[node_data[rn].prev];
       
  1602         }
       
  1603 
       
  1604         if (root != _graph.source(arc)) {
       
  1605           orientComponent(root, rn, order_map, pred_map,
       
  1606                           node_data, arc_lists, flip_map, type_map);
       
  1607           markFacePath(root, root, order_map, node_data);
       
  1608           int xlp = markExternalPath(xnode, order_map, child_lists,
       
  1609                                      pred_map, ancestor_map, low_map);
       
  1610           int ylp = markExternalPath(ynode, order_map, child_lists,
       
  1611                                      pred_map, ancestor_map, low_map);
       
  1612           markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
       
  1613           Node lwnode = findPertinent(ynode, order_map, node_data,
       
  1614                                       embed_arc, merge_roots);
       
  1615 
       
  1616           markPertinentPath(lwnode, order_map, node_data, arc_lists,
       
  1617                             embed_arc, merge_roots);
       
  1618 
       
  1619           return;
       
  1620         }
       
  1621       }
       
  1622 
       
  1623       orientComponent(root, rn, order_map, pred_map,
       
  1624                       node_data, arc_lists, flip_map, type_map);
       
  1625 
       
  1626       Node wnode = findPertinent(ynode, order_map, node_data,
       
  1627                                  embed_arc, merge_roots);
       
  1628       setFaceFlags(root, wnode, ynode, xnode, order_map, node_data, type_map);
       
  1629 
       
  1630 
       
  1631       //Minor-B
       
  1632       if (!merge_roots[wnode].empty()) {
       
  1633         int cn = merge_roots[wnode].back();
       
  1634         Node rep = order_list[cn - order_list.size()];
       
  1635         if (low_map[rep] < rorder) {
       
  1636           markFacePath(root, root, order_map, node_data);
       
  1637           int xlp = markExternalPath(xnode, order_map, child_lists,
       
  1638                                      pred_map, ancestor_map, low_map);
       
  1639           int ylp = markExternalPath(ynode, order_map, child_lists,
       
  1640                                      pred_map, ancestor_map, low_map);
       
  1641 
       
  1642           Node lwnode, lznode;
       
  1643           markCommonPath(wnode, rorder, lwnode, lznode, order_list,
       
  1644                          order_map, node_data, arc_lists, embed_arc,
       
  1645                          merge_roots, child_lists, ancestor_map, low_map);
       
  1646 
       
  1647           markPertinentPath(lwnode, order_map, node_data, arc_lists,
       
  1648                             embed_arc, merge_roots);
       
  1649           int zlp = markExternalPath(lznode, order_map, child_lists,
       
  1650                                      pred_map, ancestor_map, low_map);
       
  1651 
       
  1652           int minlp = xlp < ylp ? xlp : ylp;
       
  1653           if (zlp < minlp) minlp = zlp;
       
  1654 
       
  1655           int maxlp = xlp > ylp ? xlp : ylp;
       
  1656           if (zlp > maxlp) maxlp = zlp;
       
  1657 
       
  1658           markPredPath(order_list[maxlp], order_list[minlp], pred_map);
       
  1659 
       
  1660           return;
       
  1661         }
       
  1662       }
       
  1663 
       
  1664       Node pxnode, pynode;
       
  1665       std::vector<Arc> ipath;
       
  1666       findInternalPath(ipath, wnode, root, type_map, order_map,
       
  1667                        node_data, arc_lists);
       
  1668       setInternalFlags(ipath, type_map);
       
  1669       pynode = _graph.source(ipath.front());
       
  1670       pxnode = _graph.target(ipath.back());
       
  1671 
       
  1672       wnode = findPertinent(pynode, order_map, node_data,
       
  1673                             embed_arc, merge_roots);
       
  1674 
       
  1675       // Minor-C
       
  1676       {
       
  1677         if (type_map[_graph.source(ipath.front())] == HIGHY) {
       
  1678           if (type_map[_graph.target(ipath.back())] == HIGHX) {
       
  1679             markFacePath(xnode, pxnode, order_map, node_data);
       
  1680           }
       
  1681           markFacePath(root, xnode, order_map, node_data);
       
  1682           markPertinentPath(wnode, order_map, node_data, arc_lists,
       
  1683                             embed_arc, merge_roots);
       
  1684           markInternalPath(ipath);
       
  1685           int xlp = markExternalPath(xnode, order_map, child_lists,
       
  1686                                      pred_map, ancestor_map, low_map);
       
  1687           int ylp = markExternalPath(ynode, order_map, child_lists,
       
  1688                                      pred_map, ancestor_map, low_map);
       
  1689           markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
       
  1690           return;
       
  1691         }
       
  1692 
       
  1693         if (type_map[_graph.target(ipath.back())] == HIGHX) {
       
  1694           markFacePath(ynode, root, order_map, node_data);
       
  1695           markPertinentPath(wnode, order_map, node_data, arc_lists,
       
  1696                             embed_arc, merge_roots);
       
  1697           markInternalPath(ipath);
       
  1698           int xlp = markExternalPath(xnode, order_map, child_lists,
       
  1699                                      pred_map, ancestor_map, low_map);
       
  1700           int ylp = markExternalPath(ynode, order_map, child_lists,
       
  1701                                      pred_map, ancestor_map, low_map);
       
  1702           markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
       
  1703           return;
       
  1704         }
       
  1705       }
       
  1706 
       
  1707       std::vector<Arc> ppath;
       
  1708       findPilePath(ppath, root, type_map, order_map, node_data, arc_lists);
       
  1709 
       
  1710       // Minor-D
       
  1711       if (!ppath.empty()) {
       
  1712         markFacePath(ynode, xnode, order_map, node_data);
       
  1713         markPertinentPath(wnode, order_map, node_data, arc_lists,
       
  1714                           embed_arc, merge_roots);
       
  1715         markPilePath(ppath);
       
  1716         markInternalPath(ipath);
       
  1717         int xlp = markExternalPath(xnode, order_map, child_lists,
       
  1718                                    pred_map, ancestor_map, low_map);
       
  1719         int ylp = markExternalPath(ynode, order_map, child_lists,
       
  1720                                    pred_map, ancestor_map, low_map);
       
  1721         markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
       
  1722         return;
       
  1723       }
       
  1724 
       
  1725       // Minor-E*
       
  1726       {
       
  1727 
       
  1728         if (!external(wnode, rorder, child_lists, ancestor_map, low_map)) {
       
  1729           Node znode = findExternal(pynode, rorder, order_map,
       
  1730                                     child_lists, ancestor_map,
       
  1731                                     low_map, node_data);
       
  1732 
       
  1733           if (type_map[znode] == LOWY) {
       
  1734             markFacePath(root, xnode, order_map, node_data);
       
  1735             markPertinentPath(wnode, order_map, node_data, arc_lists,
       
  1736                               embed_arc, merge_roots);
       
  1737             markInternalPath(ipath);
       
  1738             int xlp = markExternalPath(xnode, order_map, child_lists,
       
  1739                                        pred_map, ancestor_map, low_map);
       
  1740             int zlp = markExternalPath(znode, order_map, child_lists,
       
  1741                                        pred_map, ancestor_map, low_map);
       
  1742             markPredPath(root, order_list[xlp < zlp ? xlp : zlp], pred_map);
       
  1743           } else {
       
  1744             markFacePath(ynode, root, order_map, node_data);
       
  1745             markPertinentPath(wnode, order_map, node_data, arc_lists,
       
  1746                               embed_arc, merge_roots);
       
  1747             markInternalPath(ipath);
       
  1748             int ylp = markExternalPath(ynode, order_map, child_lists,
       
  1749                                        pred_map, ancestor_map, low_map);
       
  1750             int zlp = markExternalPath(znode, order_map, child_lists,
       
  1751                                        pred_map, ancestor_map, low_map);
       
  1752             markPredPath(root, order_list[ylp < zlp ? ylp : zlp], pred_map);
       
  1753           }
       
  1754           return;
       
  1755         }
       
  1756 
       
  1757         int xlp = markExternalPath(xnode, order_map, child_lists,
       
  1758                                    pred_map, ancestor_map, low_map);
       
  1759         int ylp = markExternalPath(ynode, order_map, child_lists,
       
  1760                                    pred_map, ancestor_map, low_map);
       
  1761         int wlp = markExternalPath(wnode, order_map, child_lists,
       
  1762                                    pred_map, ancestor_map, low_map);
       
  1763 
       
  1764         if (wlp > xlp && wlp > ylp) {
       
  1765           markFacePath(root, root, order_map, node_data);
       
  1766           markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
       
  1767           return;
       
  1768         }
       
  1769 
       
  1770         markInternalPath(ipath);
       
  1771         markPertinentPath(wnode, order_map, node_data, arc_lists,
       
  1772                           embed_arc, merge_roots);
       
  1773 
       
  1774         if (xlp > ylp && xlp > wlp) {
       
  1775           markFacePath(root, pynode, order_map, node_data);
       
  1776           markFacePath(wnode, xnode, order_map, node_data);
       
  1777           markPredPath(root, order_list[ylp < wlp ? ylp : wlp], pred_map);
       
  1778           return;
       
  1779         }
       
  1780 
       
  1781         if (ylp > xlp && ylp > wlp) {
       
  1782           markFacePath(pxnode, root, order_map, node_data);
       
  1783           markFacePath(ynode, wnode, order_map, node_data);
       
  1784           markPredPath(root, order_list[xlp < wlp ? xlp : wlp], pred_map);
       
  1785           return;
       
  1786         }
       
  1787 
       
  1788         if (pynode != ynode) {
       
  1789           markFacePath(pxnode, wnode, order_map, node_data);
       
  1790 
       
  1791           int minlp = xlp < ylp ? xlp : ylp;
       
  1792           if (wlp < minlp) minlp = wlp;
       
  1793 
       
  1794           int maxlp = xlp > ylp ? xlp : ylp;
       
  1795           if (wlp > maxlp) maxlp = wlp;
       
  1796 
       
  1797           markPredPath(order_list[maxlp], order_list[minlp], pred_map);
       
  1798           return;
       
  1799         }
       
  1800 
       
  1801         if (pxnode != xnode) {
       
  1802           markFacePath(wnode, pynode, order_map, node_data);
       
  1803 
       
  1804           int minlp = xlp < ylp ? xlp : ylp;
       
  1805           if (wlp < minlp) minlp = wlp;
       
  1806 
       
  1807           int maxlp = xlp > ylp ? xlp : ylp;
       
  1808           if (wlp > maxlp) maxlp = wlp;
       
  1809 
       
  1810           markPredPath(order_list[maxlp], order_list[minlp], pred_map);
       
  1811           return;
       
  1812         }
       
  1813 
       
  1814         markFacePath(root, root, order_map, node_data);
       
  1815         int minlp = xlp < ylp ? xlp : ylp;
       
  1816         if (wlp < minlp) minlp = wlp;
       
  1817         markPredPath(root, order_list[minlp], pred_map);
       
  1818         return;
       
  1819       }
       
  1820 
       
  1821     }
       
  1822 
       
  1823   };
       
  1824 
       
  1825   namespace _planarity_bits {
       
  1826 
       
  1827     template <typename Graph, typename EmbeddingMap>
       
  1828     void makeConnected(Graph& graph, EmbeddingMap& embedding) {
       
  1829       DfsVisitor<Graph> null_visitor;
       
  1830       DfsVisit<Graph, DfsVisitor<Graph> > dfs(graph, null_visitor);
       
  1831       dfs.init();
       
  1832 
       
  1833       typename Graph::Node u = INVALID;
       
  1834       for (typename Graph::NodeIt n(graph); n != INVALID; ++n) {
       
  1835         if (!dfs.reached(n)) {
       
  1836           dfs.addSource(n);
       
  1837           dfs.start();
       
  1838           if (u == INVALID) {
       
  1839             u = n;
       
  1840           } else {
       
  1841             typename Graph::Node v = n;
       
  1842 
       
  1843             typename Graph::Arc ue = typename Graph::OutArcIt(graph, u);
       
  1844             typename Graph::Arc ve = typename Graph::OutArcIt(graph, v);
       
  1845 
       
  1846             typename Graph::Arc e = graph.direct(graph.addEdge(u, v), true);
       
  1847 
       
  1848             if (ue != INVALID) {
       
  1849               embedding[e] = embedding[ue];
       
  1850               embedding[ue] = e;
       
  1851             } else {
       
  1852               embedding[e] = e;
       
  1853             }
       
  1854 
       
  1855             if (ve != INVALID) {
       
  1856               embedding[graph.oppositeArc(e)] = embedding[ve];
       
  1857               embedding[ve] = graph.oppositeArc(e);
       
  1858             } else {
       
  1859               embedding[graph.oppositeArc(e)] = graph.oppositeArc(e);
       
  1860             }
       
  1861           }
       
  1862         }
       
  1863       }
       
  1864     }
       
  1865 
       
  1866     template <typename Graph, typename EmbeddingMap>
       
  1867     void makeBiNodeConnected(Graph& graph, EmbeddingMap& embedding) {
       
  1868       typename Graph::template ArcMap<bool> processed(graph);
       
  1869 
       
  1870       std::vector<typename Graph::Arc> arcs;
       
  1871       for (typename Graph::ArcIt e(graph); e != INVALID; ++e) {
       
  1872         arcs.push_back(e);
       
  1873       }
       
  1874 
       
  1875       IterableBoolMap<Graph, typename Graph::Node> visited(graph, false);
       
  1876 
       
  1877       for (int i = 0; i < int(arcs.size()); ++i) {
       
  1878         typename Graph::Arc pp = arcs[i];
       
  1879         if (processed[pp]) continue;
       
  1880 
       
  1881         typename Graph::Arc e = embedding[graph.oppositeArc(pp)];
       
  1882         processed[e] = true;
       
  1883         visited.set(graph.source(e), true);
       
  1884 
       
  1885         typename Graph::Arc p = e, l = e;
       
  1886         e = embedding[graph.oppositeArc(e)];
       
  1887 
       
  1888         while (e != l) {
       
  1889           processed[e] = true;
       
  1890 
       
  1891           if (visited[graph.source(e)]) {
       
  1892 
       
  1893             typename Graph::Arc n =
       
  1894               graph.direct(graph.addEdge(graph.source(p),
       
  1895                                            graph.target(e)), true);
       
  1896             embedding[n] = p;
       
  1897             embedding[graph.oppositeArc(pp)] = n;
       
  1898 
       
  1899             embedding[graph.oppositeArc(n)] =
       
  1900               embedding[graph.oppositeArc(e)];
       
  1901             embedding[graph.oppositeArc(e)] =
       
  1902               graph.oppositeArc(n);
       
  1903 
       
  1904             p = n;
       
  1905             e = embedding[graph.oppositeArc(n)];
       
  1906           } else {
       
  1907             visited.set(graph.source(e), true);
       
  1908             pp = p;
       
  1909             p = e;
       
  1910             e = embedding[graph.oppositeArc(e)];
       
  1911           }
       
  1912         }
       
  1913         visited.setAll(false);
       
  1914       }
       
  1915     }
       
  1916 
       
  1917 
       
  1918     template <typename Graph, typename EmbeddingMap>
       
  1919     void makeMaxPlanar(Graph& graph, EmbeddingMap& embedding) {
       
  1920 
       
  1921       typename Graph::template NodeMap<int> degree(graph);
       
  1922 
       
  1923       for (typename Graph::NodeIt n(graph); n != INVALID; ++n) {
       
  1924         degree[n] = countIncEdges(graph, n);
       
  1925       }
       
  1926 
       
  1927       typename Graph::template ArcMap<bool> processed(graph);
       
  1928       IterableBoolMap<Graph, typename Graph::Node> visited(graph, false);
       
  1929 
       
  1930       std::vector<typename Graph::Arc> arcs;
       
  1931       for (typename Graph::ArcIt e(graph); e != INVALID; ++e) {
       
  1932         arcs.push_back(e);
       
  1933       }
       
  1934 
       
  1935       for (int i = 0; i < int(arcs.size()); ++i) {
       
  1936         typename Graph::Arc e = arcs[i];
       
  1937 
       
  1938         if (processed[e]) continue;
       
  1939         processed[e] = true;
       
  1940 
       
  1941         typename Graph::Arc mine = e;
       
  1942         int mind = degree[graph.source(e)];
       
  1943 
       
  1944         int face_size = 1;
       
  1945 
       
  1946         typename Graph::Arc l = e;
       
  1947         e = embedding[graph.oppositeArc(e)];
       
  1948         while (l != e) {
       
  1949           processed[e] = true;
       
  1950 
       
  1951           ++face_size;
       
  1952 
       
  1953           if (degree[graph.source(e)] < mind) {
       
  1954             mine = e;
       
  1955             mind = degree[graph.source(e)];
       
  1956           }
       
  1957 
       
  1958           e = embedding[graph.oppositeArc(e)];
       
  1959         }
       
  1960 
       
  1961         if (face_size < 4) {
       
  1962           continue;
       
  1963         }
       
  1964 
       
  1965         typename Graph::Node s = graph.source(mine);
       
  1966         for (typename Graph::OutArcIt e(graph, s); e != INVALID; ++e) {
       
  1967           visited.set(graph.target(e), true);
       
  1968         }
       
  1969 
       
  1970         typename Graph::Arc oppe = INVALID;
       
  1971 
       
  1972         e = embedding[graph.oppositeArc(mine)];
       
  1973         e = embedding[graph.oppositeArc(e)];
       
  1974         while (graph.target(e) != s) {
       
  1975           if (visited[graph.source(e)]) {
       
  1976             oppe = e;
       
  1977             break;
       
  1978           }
       
  1979           e = embedding[graph.oppositeArc(e)];
       
  1980         }
       
  1981         visited.setAll(false);
       
  1982 
       
  1983         if (oppe == INVALID) {
       
  1984 
       
  1985           e = embedding[graph.oppositeArc(mine)];
       
  1986           typename Graph::Arc pn = mine, p = e;
       
  1987 
       
  1988           e = embedding[graph.oppositeArc(e)];
       
  1989           while (graph.target(e) != s) {
       
  1990             typename Graph::Arc n =
       
  1991               graph.direct(graph.addEdge(s, graph.source(e)), true);
       
  1992 
       
  1993             embedding[n] = pn;
       
  1994             embedding[graph.oppositeArc(n)] = e;
       
  1995             embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);
       
  1996 
       
  1997             pn = n;
       
  1998 
       
  1999             p = e;
       
  2000             e = embedding[graph.oppositeArc(e)];
       
  2001           }
       
  2002 
       
  2003           embedding[graph.oppositeArc(e)] = pn;
       
  2004 
       
  2005         } else {
       
  2006 
       
  2007           mine = embedding[graph.oppositeArc(mine)];
       
  2008           s = graph.source(mine);
       
  2009           oppe = embedding[graph.oppositeArc(oppe)];
       
  2010           typename Graph::Node t = graph.source(oppe);
       
  2011 
       
  2012           typename Graph::Arc ce = graph.direct(graph.addEdge(s, t), true);
       
  2013           embedding[ce] = mine;
       
  2014           embedding[graph.oppositeArc(ce)] = oppe;
       
  2015 
       
  2016           typename Graph::Arc pn = ce, p = oppe;
       
  2017           e = embedding[graph.oppositeArc(oppe)];
       
  2018           while (graph.target(e) != s) {
       
  2019             typename Graph::Arc n =
       
  2020               graph.direct(graph.addEdge(s, graph.source(e)), true);
       
  2021 
       
  2022             embedding[n] = pn;
       
  2023             embedding[graph.oppositeArc(n)] = e;
       
  2024             embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);
       
  2025 
       
  2026             pn = n;
       
  2027 
       
  2028             p = e;
       
  2029             e = embedding[graph.oppositeArc(e)];
       
  2030 
       
  2031           }
       
  2032           embedding[graph.oppositeArc(e)] = pn;
       
  2033 
       
  2034           pn = graph.oppositeArc(ce), p = mine;
       
  2035           e = embedding[graph.oppositeArc(mine)];
       
  2036           while (graph.target(e) != t) {
       
  2037             typename Graph::Arc n =
       
  2038               graph.direct(graph.addEdge(t, graph.source(e)), true);
       
  2039 
       
  2040             embedding[n] = pn;
       
  2041             embedding[graph.oppositeArc(n)] = e;
       
  2042             embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);
       
  2043 
       
  2044             pn = n;
       
  2045 
       
  2046             p = e;
       
  2047             e = embedding[graph.oppositeArc(e)];
       
  2048 
       
  2049           }
       
  2050           embedding[graph.oppositeArc(e)] = pn;
       
  2051         }
       
  2052       }
       
  2053     }
       
  2054 
       
  2055   }
       
  2056 
       
  2057   /// \ingroup planar
       
  2058   ///
       
  2059   /// \brief Schnyder's planar drawing algorithm
       
  2060   ///
       
  2061   /// The planar drawing algorithm calculates positions for the nodes
       
  2062   /// in the plane which coordinates satisfy that if the arcs are
       
  2063   /// represented with straight lines then they will not intersect
       
  2064   /// each other.
       
  2065   ///
       
  2066   /// Scnyder's algorithm embeds the graph on \c (n-2,n-2) size grid,
       
  2067   /// i.e. each node will be located in the \c [0,n-2]x[0,n-2] square.
       
  2068   /// The time complexity of the algorithm is O(n).
       
  2069   template <typename Graph>
       
  2070   class PlanarDrawing {
       
  2071   public:
       
  2072 
       
  2073     TEMPLATE_GRAPH_TYPEDEFS(Graph);
       
  2074 
       
  2075     /// \brief The point type for store coordinates
       
  2076     typedef dim2::Point<int> Point;
       
  2077     /// \brief The map type for store coordinates
       
  2078     typedef typename Graph::template NodeMap<Point> PointMap;
       
  2079 
       
  2080 
       
  2081     /// \brief Constructor
       
  2082     ///
       
  2083     /// Constructor
       
  2084     /// \pre The graph should be simple, i.e. loop and parallel arc free.
       
  2085     PlanarDrawing(const Graph& graph)
       
  2086       : _graph(graph), _point_map(graph) {}
       
  2087 
       
  2088   private:
       
  2089 
       
  2090     template <typename AuxGraph, typename AuxEmbeddingMap>
       
  2091     void drawing(const AuxGraph& graph,
       
  2092                  const AuxEmbeddingMap& next,
       
  2093                  PointMap& point_map) {
       
  2094       TEMPLATE_GRAPH_TYPEDEFS(AuxGraph);
       
  2095 
       
  2096       typename AuxGraph::template ArcMap<Arc> prev(graph);
       
  2097 
       
  2098       for (NodeIt n(graph); n != INVALID; ++n) {
       
  2099         Arc e = OutArcIt(graph, n);
       
  2100 
       
  2101         Arc p = e, l = e;
       
  2102 
       
  2103         e = next[e];
       
  2104         while (e != l) {
       
  2105           prev[e] = p;
       
  2106           p = e;
       
  2107           e = next[e];
       
  2108         }
       
  2109         prev[e] = p;
       
  2110       }
       
  2111 
       
  2112       Node anode, bnode, cnode;
       
  2113 
       
  2114       {
       
  2115         Arc e = ArcIt(graph);
       
  2116         anode = graph.source(e);
       
  2117         bnode = graph.target(e);
       
  2118         cnode = graph.target(next[graph.oppositeArc(e)]);
       
  2119       }
       
  2120 
       
  2121       IterableBoolMap<AuxGraph, Node> proper(graph, false);
       
  2122       typename AuxGraph::template NodeMap<int> conn(graph, -1);
       
  2123 
       
  2124       conn[anode] = conn[bnode] = -2;
       
  2125       {
       
  2126         for (OutArcIt e(graph, anode); e != INVALID; ++e) {
       
  2127           Node m = graph.target(e);
       
  2128           if (conn[m] == -1) {
       
  2129             conn[m] = 1;
       
  2130           }
       
  2131         }
       
  2132         conn[cnode] = 2;
       
  2133 
       
  2134         for (OutArcIt e(graph, bnode); e != INVALID; ++e) {
       
  2135           Node m = graph.target(e);
       
  2136           if (conn[m] == -1) {
       
  2137             conn[m] = 1;
       
  2138           } else if (conn[m] != -2) {
       
  2139             conn[m] += 1;
       
  2140             Arc pe = graph.oppositeArc(e);
       
  2141             if (conn[graph.target(next[pe])] == -2) {
       
  2142               conn[m] -= 1;
       
  2143             }
       
  2144             if (conn[graph.target(prev[pe])] == -2) {
       
  2145               conn[m] -= 1;
       
  2146             }
       
  2147 
       
  2148             proper.set(m, conn[m] == 1);
       
  2149           }
       
  2150         }
       
  2151       }
       
  2152 
       
  2153 
       
  2154       typename AuxGraph::template ArcMap<int> angle(graph, -1);
       
  2155 
       
  2156       while (proper.trueNum() != 0) {
       
  2157         Node n = typename IterableBoolMap<AuxGraph, Node>::TrueIt(proper);
       
  2158         proper.set(n, false);
       
  2159         conn[n] = -2;
       
  2160 
       
  2161         for (OutArcIt e(graph, n); e != INVALID; ++e) {
       
  2162           Node m = graph.target(e);
       
  2163           if (conn[m] == -1) {
       
  2164             conn[m] = 1;
       
  2165           } else if (conn[m] != -2) {
       
  2166             conn[m] += 1;
       
  2167             Arc pe = graph.oppositeArc(e);
       
  2168             if (conn[graph.target(next[pe])] == -2) {
       
  2169               conn[m] -= 1;
       
  2170             }
       
  2171             if (conn[graph.target(prev[pe])] == -2) {
       
  2172               conn[m] -= 1;
       
  2173             }
       
  2174 
       
  2175             proper.set(m, conn[m] == 1);
       
  2176           }
       
  2177         }
       
  2178 
       
  2179         {
       
  2180           Arc e = OutArcIt(graph, n);
       
  2181           Arc p = e, l = e;
       
  2182 
       
  2183           e = next[e];
       
  2184           while (e != l) {
       
  2185 
       
  2186             if (conn[graph.target(e)] == -2 && conn[graph.target(p)] == -2) {
       
  2187               Arc f = e;
       
  2188               angle[f] = 0;
       
  2189               f = next[graph.oppositeArc(f)];
       
  2190               angle[f] = 1;
       
  2191               f = next[graph.oppositeArc(f)];
       
  2192               angle[f] = 2;
       
  2193             }
       
  2194 
       
  2195             p = e;
       
  2196             e = next[e];
       
  2197           }
       
  2198 
       
  2199           if (conn[graph.target(e)] == -2 && conn[graph.target(p)] == -2) {
       
  2200             Arc f = e;
       
  2201             angle[f] = 0;
       
  2202             f = next[graph.oppositeArc(f)];
       
  2203             angle[f] = 1;
       
  2204             f = next[graph.oppositeArc(f)];
       
  2205             angle[f] = 2;
       
  2206           }
       
  2207         }
       
  2208       }
       
  2209 
       
  2210       typename AuxGraph::template NodeMap<Node> apred(graph, INVALID);
       
  2211       typename AuxGraph::template NodeMap<Node> bpred(graph, INVALID);
       
  2212       typename AuxGraph::template NodeMap<Node> cpred(graph, INVALID);
       
  2213 
       
  2214       typename AuxGraph::template NodeMap<int> apredid(graph, -1);
       
  2215       typename AuxGraph::template NodeMap<int> bpredid(graph, -1);
       
  2216       typename AuxGraph::template NodeMap<int> cpredid(graph, -1);
       
  2217 
       
  2218       for (ArcIt e(graph); e != INVALID; ++e) {
       
  2219         if (angle[e] == angle[next[e]]) {
       
  2220           switch (angle[e]) {
       
  2221           case 2:
       
  2222             apred[graph.target(e)] = graph.source(e);
       
  2223             apredid[graph.target(e)] = graph.id(graph.source(e));
       
  2224             break;
       
  2225           case 1:
       
  2226             bpred[graph.target(e)] = graph.source(e);
       
  2227             bpredid[graph.target(e)] = graph.id(graph.source(e));
       
  2228             break;
       
  2229           case 0:
       
  2230             cpred[graph.target(e)] = graph.source(e);
       
  2231             cpredid[graph.target(e)] = graph.id(graph.source(e));
       
  2232             break;
       
  2233           }
       
  2234         }
       
  2235       }
       
  2236 
       
  2237       cpred[anode] = INVALID;
       
  2238       cpred[bnode] = INVALID;
       
  2239 
       
  2240       std::vector<Node> aorder, border, corder;
       
  2241 
       
  2242       {
       
  2243         typename AuxGraph::template NodeMap<bool> processed(graph, false);
       
  2244         std::vector<Node> st;
       
  2245         for (NodeIt n(graph); n != INVALID; ++n) {
       
  2246           if (!processed[n] && n != bnode && n != cnode) {
       
  2247             st.push_back(n);
       
  2248             processed[n] = true;
       
  2249             Node m = apred[n];
       
  2250             while (m != INVALID && !processed[m]) {
       
  2251               st.push_back(m);
       
  2252               processed[m] = true;
       
  2253               m = apred[m];
       
  2254             }
       
  2255             while (!st.empty()) {
       
  2256               aorder.push_back(st.back());
       
  2257               st.pop_back();
       
  2258             }
       
  2259           }
       
  2260         }
       
  2261       }
       
  2262 
       
  2263       {
       
  2264         typename AuxGraph::template NodeMap<bool> processed(graph, false);
       
  2265         std::vector<Node> st;
       
  2266         for (NodeIt n(graph); n != INVALID; ++n) {
       
  2267           if (!processed[n] && n != cnode && n != anode) {
       
  2268             st.push_back(n);
       
  2269             processed[n] = true;
       
  2270             Node m = bpred[n];
       
  2271             while (m != INVALID && !processed[m]) {
       
  2272               st.push_back(m);
       
  2273               processed[m] = true;
       
  2274               m = bpred[m];
       
  2275             }
       
  2276             while (!st.empty()) {
       
  2277               border.push_back(st.back());
       
  2278               st.pop_back();
       
  2279             }
       
  2280           }
       
  2281         }
       
  2282       }
       
  2283 
       
  2284       {
       
  2285         typename AuxGraph::template NodeMap<bool> processed(graph, false);
       
  2286         std::vector<Node> st;
       
  2287         for (NodeIt n(graph); n != INVALID; ++n) {
       
  2288           if (!processed[n] && n != anode && n != bnode) {
       
  2289             st.push_back(n);
       
  2290             processed[n] = true;
       
  2291             Node m = cpred[n];
       
  2292             while (m != INVALID && !processed[m]) {
       
  2293               st.push_back(m);
       
  2294               processed[m] = true;
       
  2295               m = cpred[m];
       
  2296             }
       
  2297             while (!st.empty()) {
       
  2298               corder.push_back(st.back());
       
  2299               st.pop_back();
       
  2300             }
       
  2301           }
       
  2302         }
       
  2303       }
       
  2304 
       
  2305       typename AuxGraph::template NodeMap<int> atree(graph, 0);
       
  2306       for (int i = aorder.size() - 1; i >= 0; --i) {
       
  2307         Node n = aorder[i];
       
  2308         atree[n] = 1;
       
  2309         for (OutArcIt e(graph, n); e != INVALID; ++e) {
       
  2310           if (apred[graph.target(e)] == n) {
       
  2311             atree[n] += atree[graph.target(e)];
       
  2312           }
       
  2313         }
       
  2314       }
       
  2315 
       
  2316       typename AuxGraph::template NodeMap<int> btree(graph, 0);
       
  2317       for (int i = border.size() - 1; i >= 0; --i) {
       
  2318         Node n = border[i];
       
  2319         btree[n] = 1;
       
  2320         for (OutArcIt e(graph, n); e != INVALID; ++e) {
       
  2321           if (bpred[graph.target(e)] == n) {
       
  2322             btree[n] += btree[graph.target(e)];
       
  2323           }
       
  2324         }
       
  2325       }
       
  2326 
       
  2327       typename AuxGraph::template NodeMap<int> apath(graph, 0);
       
  2328       apath[bnode] = apath[cnode] = 1;
       
  2329       typename AuxGraph::template NodeMap<int> apath_btree(graph, 0);
       
  2330       apath_btree[bnode] = btree[bnode];
       
  2331       for (int i = 1; i < int(aorder.size()); ++i) {
       
  2332         Node n = aorder[i];
       
  2333         apath[n] = apath[apred[n]] + 1;
       
  2334         apath_btree[n] = btree[n] + apath_btree[apred[n]];
       
  2335       }
       
  2336 
       
  2337       typename AuxGraph::template NodeMap<int> bpath_atree(graph, 0);
       
  2338       bpath_atree[anode] = atree[anode];
       
  2339       for (int i = 1; i < int(border.size()); ++i) {
       
  2340         Node n = border[i];
       
  2341         bpath_atree[n] = atree[n] + bpath_atree[bpred[n]];
       
  2342       }
       
  2343 
       
  2344       typename AuxGraph::template NodeMap<int> cpath(graph, 0);
       
  2345       cpath[anode] = cpath[bnode] = 1;
       
  2346       typename AuxGraph::template NodeMap<int> cpath_atree(graph, 0);
       
  2347       cpath_atree[anode] = atree[anode];
       
  2348       typename AuxGraph::template NodeMap<int> cpath_btree(graph, 0);
       
  2349       cpath_btree[bnode] = btree[bnode];
       
  2350       for (int i = 1; i < int(corder.size()); ++i) {
       
  2351         Node n = corder[i];
       
  2352         cpath[n] = cpath[cpred[n]] + 1;
       
  2353         cpath_atree[n] = atree[n] + cpath_atree[cpred[n]];
       
  2354         cpath_btree[n] = btree[n] + cpath_btree[cpred[n]];
       
  2355       }
       
  2356 
       
  2357       typename AuxGraph::template NodeMap<int> third(graph);
       
  2358       for (NodeIt n(graph); n != INVALID; ++n) {
       
  2359         point_map[n].x =
       
  2360           bpath_atree[n] + cpath_atree[n] - atree[n] - cpath[n] + 1;
       
  2361         point_map[n].y =
       
  2362           cpath_btree[n] + apath_btree[n] - btree[n] - apath[n] + 1;
       
  2363       }
       
  2364 
       
  2365     }
       
  2366 
       
  2367   public:
       
  2368 
       
  2369     /// \brief Calculates the node positions
       
  2370     ///
       
  2371     /// This function calculates the node positions.
       
  2372     /// \return %True if the graph is planar.
       
  2373     bool run() {
       
  2374       PlanarEmbedding<Graph> pe(_graph);
       
  2375       if (!pe.run()) return false;
       
  2376 
       
  2377       run(pe);
       
  2378       return true;
       
  2379     }
       
  2380 
       
  2381     /// \brief Calculates the node positions according to a
       
  2382     /// combinatorical embedding
       
  2383     ///
       
  2384     /// This function calculates the node locations. The \c embedding
       
  2385     /// parameter should contain a valid combinatorical embedding, i.e.
       
  2386     /// a valid cyclic order of the arcs.
       
  2387     template <typename EmbeddingMap>
       
  2388     void run(const EmbeddingMap& embedding) {
       
  2389       typedef SmartEdgeSet<Graph> AuxGraph;
       
  2390 
       
  2391       if (3 * countNodes(_graph) - 6 == countEdges(_graph)) {
       
  2392         drawing(_graph, embedding, _point_map);
       
  2393         return;
       
  2394       }
       
  2395 
       
  2396       AuxGraph aux_graph(_graph);
       
  2397       typename AuxGraph::template ArcMap<typename AuxGraph::Arc>
       
  2398         aux_embedding(aux_graph);
       
  2399 
       
  2400       {
       
  2401 
       
  2402         typename Graph::template EdgeMap<typename AuxGraph::Edge>
       
  2403           ref(_graph);
       
  2404 
       
  2405         for (EdgeIt e(_graph); e != INVALID; ++e) {
       
  2406           ref[e] = aux_graph.addEdge(_graph.u(e), _graph.v(e));
       
  2407         }
       
  2408 
       
  2409         for (EdgeIt e(_graph); e != INVALID; ++e) {
       
  2410           Arc ee = embedding[_graph.direct(e, true)];
       
  2411           aux_embedding[aux_graph.direct(ref[e], true)] =
       
  2412             aux_graph.direct(ref[ee], _graph.direction(ee));
       
  2413           ee = embedding[_graph.direct(e, false)];
       
  2414           aux_embedding[aux_graph.direct(ref[e], false)] =
       
  2415             aux_graph.direct(ref[ee], _graph.direction(ee));
       
  2416         }
       
  2417       }
       
  2418       _planarity_bits::makeConnected(aux_graph, aux_embedding);
       
  2419       _planarity_bits::makeBiNodeConnected(aux_graph, aux_embedding);
       
  2420       _planarity_bits::makeMaxPlanar(aux_graph, aux_embedding);
       
  2421       drawing(aux_graph, aux_embedding, _point_map);
       
  2422     }
       
  2423 
       
  2424     /// \brief The coordinate of the given node
       
  2425     ///
       
  2426     /// The coordinate of the given node.
       
  2427     Point operator[](const Node& node) const {
       
  2428       return _point_map[node];
       
  2429     }
       
  2430 
       
  2431     /// \brief Returns the grid embedding in a \e NodeMap.
       
  2432     ///
       
  2433     /// Returns the grid embedding in a \e NodeMap of \c dim2::Point<int> .
       
  2434     const PointMap& coords() const {
       
  2435       return _point_map;
       
  2436     }
       
  2437 
       
  2438   private:
       
  2439 
       
  2440     const Graph& _graph;
       
  2441     PointMap _point_map;
       
  2442 
       
  2443   };
       
  2444 
       
  2445   namespace _planarity_bits {
       
  2446 
       
  2447     template <typename ColorMap>
       
  2448     class KempeFilter {
       
  2449     public:
       
  2450       typedef typename ColorMap::Key Key;
       
  2451       typedef bool Value;
       
  2452 
       
  2453       KempeFilter(const ColorMap& color_map,
       
  2454                   const typename ColorMap::Value& first,
       
  2455                   const typename ColorMap::Value& second)
       
  2456         : _color_map(color_map), _first(first), _second(second) {}
       
  2457 
       
  2458       Value operator[](const Key& key) const {
       
  2459         return _color_map[key] == _first || _color_map[key] == _second;
       
  2460       }
       
  2461 
       
  2462     private:
       
  2463       const ColorMap& _color_map;
       
  2464       typename ColorMap::Value _first, _second;
       
  2465     };
       
  2466   }
       
  2467 
       
  2468   /// \ingroup planar
       
  2469   ///
       
  2470   /// \brief Coloring planar graphs
       
  2471   ///
       
  2472   /// The graph coloring problem is the coloring of the graph nodes
       
  2473   /// that there are not adjacent nodes with the same color. The
       
  2474   /// planar graphs can be always colored with four colors, it is
       
  2475   /// proved by Appel and Haken and their proofs provide a quadratic
       
  2476   /// time algorithm for four coloring, but it could not be used to
       
  2477   /// implement efficient algorithm. The five and six coloring can be
       
  2478   /// made in linear time, but in this class the five coloring has
       
  2479   /// quadratic worst case time complexity. The two coloring (if
       
  2480   /// possible) is solvable with a graph search algorithm and it is
       
  2481   /// implemented in \ref bipartitePartitions() function in LEMON. To
       
  2482   /// decide whether the planar graph is three colorable is
       
  2483   /// NP-complete.
       
  2484   ///
       
  2485   /// This class contains member functions for calculate colorings
       
  2486   /// with five and six colors. The six coloring algorithm is a simple
       
  2487   /// greedy coloring on the backward minimum outgoing order of nodes.
       
  2488   /// This order can be computed as in each phase the node with least
       
  2489   /// outgoing arcs to unprocessed nodes is chosen. This order
       
  2490   /// guarantees that when a node is chosen for coloring it has at
       
  2491   /// most five already colored adjacents. The five coloring algorithm
       
  2492   /// use the same method, but if the greedy approach fails to color
       
  2493   /// with five colors, i.e. the node has five already different
       
  2494   /// colored neighbours, it swaps the colors in one of the connected
       
  2495   /// two colored sets with the Kempe recoloring method.
       
  2496   template <typename Graph>
       
  2497   class PlanarColoring {
       
  2498   public:
       
  2499 
       
  2500     TEMPLATE_GRAPH_TYPEDEFS(Graph);
       
  2501 
       
  2502     /// \brief The map type for store color indexes
       
  2503     typedef typename Graph::template NodeMap<int> IndexMap;
       
  2504     /// \brief The map type for store colors
       
  2505     typedef ComposeMap<Palette, IndexMap> ColorMap;
       
  2506 
       
  2507     /// \brief Constructor
       
  2508     ///
       
  2509     /// Constructor
       
  2510     /// \pre The graph should be simple, i.e. loop and parallel arc free.
       
  2511     PlanarColoring(const Graph& graph)
       
  2512       : _graph(graph), _color_map(graph), _palette(0) {
       
  2513       _palette.add(Color(1,0,0));
       
  2514       _palette.add(Color(0,1,0));
       
  2515       _palette.add(Color(0,0,1));
       
  2516       _palette.add(Color(1,1,0));
       
  2517       _palette.add(Color(1,0,1));
       
  2518       _palette.add(Color(0,1,1));
       
  2519     }
       
  2520 
       
  2521     /// \brief Returns the \e NodeMap of color indexes
       
  2522     ///
       
  2523     /// Returns the \e NodeMap of color indexes. The values are in the
       
  2524     /// range \c [0..4] or \c [0..5] according to the coloring method.
       
  2525     IndexMap colorIndexMap() const {
       
  2526       return _color_map;
       
  2527     }
       
  2528 
       
  2529     /// \brief Returns the \e NodeMap of colors
       
  2530     ///
       
  2531     /// Returns the \e NodeMap of colors. The values are five or six
       
  2532     /// distinct \ref lemon::Color "colors".
       
  2533     ColorMap colorMap() const {
       
  2534       return composeMap(_palette, _color_map);
       
  2535     }
       
  2536 
       
  2537     /// \brief Returns the color index of the node
       
  2538     ///
       
  2539     /// Returns the color index of the node. The values are in the
       
  2540     /// range \c [0..4] or \c [0..5] according to the coloring method.
       
  2541     int colorIndex(const Node& node) const {
       
  2542       return _color_map[node];
       
  2543     }
       
  2544 
       
  2545     /// \brief Returns the color of the node
       
  2546     ///
       
  2547     /// Returns the color of the node. The values are five or six
       
  2548     /// distinct \ref lemon::Color "colors".
       
  2549     Color color(const Node& node) const {
       
  2550       return _palette[_color_map[node]];
       
  2551     }
       
  2552 
       
  2553 
       
  2554     /// \brief Calculates a coloring with at most six colors
       
  2555     ///
       
  2556     /// This function calculates a coloring with at most six colors. The time
       
  2557     /// complexity of this variant is linear in the size of the graph.
       
  2558     /// \return %True when the algorithm could color the graph with six color.
       
  2559     /// If the algorithm fails, then the graph could not be planar.
       
  2560     /// \note This function can return true if the graph is not
       
  2561     /// planar but it can be colored with 6 colors.
       
  2562     bool runSixColoring() {
       
  2563 
       
  2564       typename Graph::template NodeMap<int> heap_index(_graph, -1);
       
  2565       BucketHeap<typename Graph::template NodeMap<int> > heap(heap_index);
       
  2566 
       
  2567       for (NodeIt n(_graph); n != INVALID; ++n) {
       
  2568         _color_map[n] = -2;
       
  2569         heap.push(n, countOutArcs(_graph, n));
       
  2570       }
       
  2571 
       
  2572       std::vector<Node> order;
       
  2573 
       
  2574       while (!heap.empty()) {
       
  2575         Node n = heap.top();
       
  2576         heap.pop();
       
  2577         _color_map[n] = -1;
       
  2578         order.push_back(n);
       
  2579         for (OutArcIt e(_graph, n); e != INVALID; ++e) {
       
  2580           Node t = _graph.runningNode(e);
       
  2581           if (_color_map[t] == -2) {
       
  2582             heap.decrease(t, heap[t] - 1);
       
  2583           }
       
  2584         }
       
  2585       }
       
  2586 
       
  2587       for (int i = order.size() - 1; i >= 0; --i) {
       
  2588         std::vector<bool> forbidden(6, false);
       
  2589         for (OutArcIt e(_graph, order[i]); e != INVALID; ++e) {
       
  2590           Node t = _graph.runningNode(e);
       
  2591           if (_color_map[t] != -1) {
       
  2592             forbidden[_color_map[t]] = true;
       
  2593           }
       
  2594         }
       
  2595                for (int k = 0; k < 6; ++k) {
       
  2596           if (!forbidden[k]) {
       
  2597             _color_map[order[i]] = k;
       
  2598             break;
       
  2599           }
       
  2600         }
       
  2601         if (_color_map[order[i]] == -1) {
       
  2602           return false;
       
  2603         }
       
  2604       }
       
  2605       return true;
       
  2606     }
       
  2607 
       
  2608   private:
       
  2609 
       
  2610     bool recolor(const Node& u, const Node& v) {
       
  2611       int ucolor = _color_map[u];
       
  2612       int vcolor = _color_map[v];
       
  2613       typedef _planarity_bits::KempeFilter<IndexMap> KempeFilter;
       
  2614       KempeFilter filter(_color_map, ucolor, vcolor);
       
  2615 
       
  2616       typedef FilterNodes<const Graph, const KempeFilter> KempeGraph;
       
  2617       KempeGraph kempe_graph(_graph, filter);
       
  2618 
       
  2619       std::vector<Node> comp;
       
  2620       Bfs<KempeGraph> bfs(kempe_graph);
       
  2621       bfs.init();
       
  2622       bfs.addSource(u);
       
  2623       while (!bfs.emptyQueue()) {
       
  2624         Node n = bfs.nextNode();
       
  2625         if (n == v) return false;
       
  2626         comp.push_back(n);
       
  2627         bfs.processNextNode();
       
  2628       }
       
  2629 
       
  2630       int scolor = ucolor + vcolor;
       
  2631       for (int i = 0; i < static_cast<int>(comp.size()); ++i) {
       
  2632         _color_map[comp[i]] = scolor - _color_map[comp[i]];
       
  2633       }
       
  2634 
       
  2635       return true;
       
  2636     }
       
  2637 
       
  2638     template <typename EmbeddingMap>
       
  2639     void kempeRecoloring(const Node& node, const EmbeddingMap& embedding) {
       
  2640       std::vector<Node> nodes;
       
  2641       nodes.reserve(4);
       
  2642 
       
  2643       for (Arc e = OutArcIt(_graph, node); e != INVALID; e = embedding[e]) {
       
  2644         Node t = _graph.target(e);
       
  2645         if (_color_map[t] != -1) {
       
  2646           nodes.push_back(t);
       
  2647           if (nodes.size() == 4) break;
       
  2648         }
       
  2649       }
       
  2650 
       
  2651       int color = _color_map[nodes[0]];
       
  2652       if (recolor(nodes[0], nodes[2])) {
       
  2653         _color_map[node] = color;
       
  2654       } else {
       
  2655         color = _color_map[nodes[1]];
       
  2656         recolor(nodes[1], nodes[3]);
       
  2657         _color_map[node] = color;
       
  2658       }
       
  2659     }
       
  2660 
       
  2661   public:
       
  2662 
       
  2663     /// \brief Calculates a coloring with at most five colors
       
  2664     ///
       
  2665     /// This function calculates a coloring with at most five
       
  2666     /// colors. The worst case time complexity of this variant is
       
  2667     /// quadratic in the size of the graph.
       
  2668     template <typename EmbeddingMap>
       
  2669     void runFiveColoring(const EmbeddingMap& embedding) {
       
  2670 
       
  2671       typename Graph::template NodeMap<int> heap_index(_graph, -1);
       
  2672       BucketHeap<typename Graph::template NodeMap<int> > heap(heap_index);
       
  2673 
       
  2674       for (NodeIt n(_graph); n != INVALID; ++n) {
       
  2675         _color_map[n] = -2;
       
  2676         heap.push(n, countOutArcs(_graph, n));
       
  2677       }
       
  2678 
       
  2679       std::vector<Node> order;
       
  2680 
       
  2681       while (!heap.empty()) {
       
  2682         Node n = heap.top();
       
  2683         heap.pop();
       
  2684         _color_map[n] = -1;
       
  2685         order.push_back(n);
       
  2686         for (OutArcIt e(_graph, n); e != INVALID; ++e) {
       
  2687           Node t = _graph.runningNode(e);
       
  2688           if (_color_map[t] == -2) {
       
  2689             heap.decrease(t, heap[t] - 1);
       
  2690           }
       
  2691         }
       
  2692       }
       
  2693 
       
  2694       for (int i = order.size() - 1; i >= 0; --i) {
       
  2695         std::vector<bool> forbidden(5, false);
       
  2696         for (OutArcIt e(_graph, order[i]); e != INVALID; ++e) {
       
  2697           Node t = _graph.runningNode(e);
       
  2698           if (_color_map[t] != -1) {
       
  2699             forbidden[_color_map[t]] = true;
       
  2700           }
       
  2701         }
       
  2702         for (int k = 0; k < 5; ++k) {
       
  2703           if (!forbidden[k]) {
       
  2704             _color_map[order[i]] = k;
       
  2705             break;
       
  2706           }
       
  2707         }
       
  2708         if (_color_map[order[i]] == -1) {
       
  2709           kempeRecoloring(order[i], embedding);
       
  2710         }
       
  2711       }
       
  2712     }
       
  2713 
       
  2714     /// \brief Calculates a coloring with at most five colors
       
  2715     ///
       
  2716     /// This function calculates a coloring with at most five
       
  2717     /// colors. The worst case time complexity of this variant is
       
  2718     /// quadratic in the size of the graph.
       
  2719     /// \return %True when the graph is planar.
       
  2720     bool runFiveColoring() {
       
  2721       PlanarEmbedding<Graph> pe(_graph);
       
  2722       if (!pe.run()) return false;
       
  2723 
       
  2724       runFiveColoring(pe.embeddingMap());
       
  2725       return true;
       
  2726     }
       
  2727 
       
  2728   private:
       
  2729 
       
  2730     const Graph& _graph;
       
  2731     IndexMap _color_map;
       
  2732     Palette _palette;
       
  2733   };
       
  2734 
       
  2735 }
       
  2736 
       
  2737 #endif