lemon-main: lemon/planarity.h@a7e93de12cbd

     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  */

    19 #ifndef LEMON_PLANARITY_H

    20 #define LEMON_PLANARITY_H

    22 /// \ingroup planar

    23 /// \file

    24 /// \brief Planarity checking, embedding, drawing and coloring

    26 #include <vector>

    27 #include <list>

    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>

    40 namespace lemon {

    42   namespace _planarity_bits {

    44     template <typename Graph>

    45     struct PlanarityVisitor : DfsVisitor<Graph> {

    47       TEMPLATE_GRAPH_TYPEDEFS(Graph);

    49       typedef typename Graph::template NodeMap<Arc> PredMap;

    51       typedef typename Graph::template EdgeMap<bool> TreeMap;

    53       typedef typename Graph::template NodeMap<int> OrderMap;

    54       typedef std::vector<Node> OrderList;

    56       typedef typename Graph::template NodeMap<int> LowMap;

    57       typedef typename Graph::template NodeMap<int> AncestorMap;

    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) {}

    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       }

    74       void discover(const Arc& arc) {

    75         Node source = _graph.source(arc);

    76         Node target = _graph.target(arc);

    78         _tree_map[arc] = true;

    79         _pred_map[target] = arc;

    80       }

    82       void examine(const Arc& arc) {

    83         Node source = _graph.source(arc);

    84         Node target = _graph.target(arc);

    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       }

    96       void backtrack(const Arc& arc) {

    97         Node source = _graph.source(arc);

    98         Node target = _graph.target(arc);

   100         if (_low_map[source] > _low_map[target]) {

   101           _low_map[source] = _low_map[target];

   102         }

   103       }

   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     };

   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     };

   122     template <typename Graph>

   123     struct NodeDataNode<Graph, false> {

   124       int prev, next;

   125       int visited;

   126     };

   128     template <typename Graph>

   129     struct ChildListNode {

   130       typedef typename Graph::Node Node;

   131       Node first;

   132       Node prev, next;

   133     };

   135     template <typename Graph>

   136     struct ArcListNode {

   137       typename Graph::Arc prev, next;

   138     };

   140     template <typename Graph>

   141     class PlanarityChecking {

   142     private:

   144       TEMPLATE_GRAPH_TYPEDEFS(Graph);

   146       const Graph& _graph;

   148     private:

   150       typedef typename Graph::template NodeMap<Arc> PredMap;

   152       typedef typename Graph::template EdgeMap<bool> TreeMap;

   154       typedef typename Graph::template NodeMap<int> OrderMap;

   155       typedef std::vector<Node> OrderList;

   157       typedef typename Graph::template NodeMap<int> LowMap;

   158       typedef typename Graph::template NodeMap<int> AncestorMap;

   160       typedef _planarity_bits::NodeDataNode<Graph> NodeDataNode;

   161       typedef std::vector<NodeDataNode> NodeData;

   163       typedef _planarity_bits::ChildListNode<Graph> ChildListNode;

   164       typedef typename Graph::template NodeMap<ChildListNode> ChildLists;

   166       typedef typename Graph::template NodeMap<std::list<int> > MergeRoots;

   168       typedef typename Graph::template NodeMap<bool> EmbedArc;

   170     public:

   172       PlanarityChecking(const Graph& graph) : _graph(graph) {}

   174       bool run() {

   175         typedef _planarity_bits::PlanarityVisitor<Graph> Visitor;

   177         PredMap pred_map(_graph, INVALID);

   178         TreeMap tree_map(_graph, false);

   180         OrderMap order_map(_graph, -1);

   181         OrderList order_list;

   183         AncestorMap ancestor_map(_graph, -1);

   184         LowMap low_map(_graph, -1);

   186         Visitor visitor(_graph, pred_map, tree_map,

   187                         order_map, order_list, ancestor_map, low_map);

   188         DfsVisit<Graph, Visitor> visit(_graph, visitor);

   189         visit.run();

   191         ChildLists child_lists(_graph);

   192         createChildLists(tree_map, order_map, low_map, child_lists);

   194         NodeData node_data(2 * order_list.size());

   196         EmbedArc embed_arc(_graph, false);

   198         MergeRoots merge_roots(_graph);

   200         for (int i = order_list.size() - 1; i >= 0; --i) {

   202           Node node = order_list[i];

   204           Node source = node;

   205           for (OutArcIt e(_graph, node); e != INVALID; ++e) {

   206             Node target = _graph.target(e);

   208             if (order_map[source] < order_map[target] && tree_map[e]) {

   209               initFace(target, node_data, order_map, order_list);

   210             }

   211           }

   213           for (OutArcIt e(_graph, node); e != INVALID; ++e) {

   214             Node target = _graph.target(e);

   216             if (order_map[source] < order_map[target] && !tree_map[e]) {

   217               embed_arc[target] = true;

   218               walkUp(target, source, i, pred_map, low_map,

   219                      order_map, order_list, node_data, merge_roots);

   220             }

   221           }

   223           for (typename MergeRoots::Value::iterator it =

   224                  merge_roots[node].begin();

   225                it != merge_roots[node].end(); ++it) {

   226             int rn = *it;

   227             walkDown(rn, i, node_data, order_list, child_lists,

   228                      ancestor_map, low_map, embed_arc, merge_roots);

   229           }

   230           merge_roots[node].clear();

   232           for (OutArcIt e(_graph, node); e != INVALID; ++e) {

   233             Node target = _graph.target(e);

   235             if (order_map[source] < order_map[target] && !tree_map[e]) {

   236               if (embed_arc[target]) {

   237                 return false;

   238               }

   239             }

   240           }

   241         }

   243         return true;

   244       }

   246     private:

   248       void createChildLists(const TreeMap& tree_map, const OrderMap& order_map,

   249                             const LowMap& low_map, ChildLists& child_lists) {

   251         for (NodeIt n(_graph); n != INVALID; ++n) {

   252           Node source = n;

   254           std::vector<Node> targets;

   255           for (OutArcIt e(_graph, n); e != INVALID; ++e) {

   256             Node target = _graph.target(e);

   258             if (order_map[source] < order_map[target] && tree_map[e]) {

   259               targets.push_back(target);

   260             }

   261           }

   263           if (targets.size() == 0) {

   264             child_lists[source].first = INVALID;

   265           } else if (targets.size() == 1) {

   266             child_lists[source].first = targets[0];

   267             child_lists[targets[0]].prev = INVALID;

   268             child_lists[targets[0]].next = INVALID;

   269           } else {

   270             radixSort(targets.begin(), targets.end(), mapToFunctor(low_map));

   271             for (int i = 1; i < int(targets.size()); ++i) {

   272               child_lists[targets[i]].prev = targets[i - 1];

   273               child_lists[targets[i - 1]].next = targets[i];

   274             }

   275             child_lists[targets.back()].next = INVALID;

   276             child_lists[targets.front()].prev = INVALID;

   277             child_lists[source].first = targets.front();

   278           }

   279         }

   280       }

   282       void walkUp(const Node& node, Node root, int rorder,

   283                   const PredMap& pred_map, const LowMap& low_map,

   284                   const OrderMap& order_map, const OrderList& order_list,

   285                   NodeData& node_data, MergeRoots& merge_roots) {

   287         int na, nb;

   288         bool da, db;

   290         na = nb = order_map[node];

   291         da = true; db = false;

   293         while (true) {

   295           if (node_data[na].visited == rorder) break;

   296           if (node_data[nb].visited == rorder) break;

   298           node_data[na].visited = rorder;

   299           node_data[nb].visited = rorder;

   301           int rn = -1;

   303           if (na >= int(order_list.size())) {

   304             rn = na;

   305           } else if (nb >= int(order_list.size())) {

   306             rn = nb;

   307           }

   309           if (rn == -1) {

   310             int nn;

   312             nn = da ? node_data[na].prev : node_data[na].next;

   313             da = node_data[nn].prev != na;

   314             na = nn;

   316             nn = db ? node_data[nb].prev : node_data[nb].next;

   317             db = node_data[nn].prev != nb;

   318             nb = nn;

   320           } else {

   322             Node rep = order_list[rn - order_list.size()];

   323             Node parent = _graph.source(pred_map[rep]);

   325             if (low_map[rep] < rorder) {

   326               merge_roots[parent].push_back(rn);

   327             } else {

   328               merge_roots[parent].push_front(rn);

   329             }

   331             if (parent != root) {

   332               na = nb = order_map[parent];

   333               da = true; db = false;

   334             } else {

   335               break;

   336             }

   337           }

   338         }

   339       }

   341       void walkDown(int rn, int rorder, NodeData& node_data,

   342                     OrderList& order_list, ChildLists& child_lists,

   343                     AncestorMap& ancestor_map, LowMap& low_map,

   344                     EmbedArc& embed_arc, MergeRoots& merge_roots) {

   346         std::vector<std::pair<int, bool> > merge_stack;

   348         for (int di = 0; di < 2; ++di) {

   349           bool rd = di == 0;

   350           int pn = rn;

   351           int n = rd ? node_data[rn].next : node_data[rn].prev;

   353           while (n != rn) {

   355             Node node = order_list[n];

   357             if (embed_arc[node]) {

   359               // Merging components on the critical path

   360               while (!merge_stack.empty()) {

   362                 // Component root

   363                 int cn = merge_stack.back().first;

   364                 bool cd = merge_stack.back().second;

   365                 merge_stack.pop_back();

   367                 // Parent of component

   368                 int dn = merge_stack.back().first;

   369                 bool dd = merge_stack.back().second;

   370                 merge_stack.pop_back();

   372                 Node parent = order_list[dn];

   374                 // Erasing from merge_roots

   375                 merge_roots[parent].pop_front();

   377                 Node child = order_list[cn - order_list.size()];

   379                 // Erasing from child_lists

   380                 if (child_lists[child].prev != INVALID) {

   381                   child_lists[child_lists[child].prev].next =

   382                     child_lists[child].next;

   383                 } else {

   384                   child_lists[parent].first = child_lists[child].next;

   385                 }

   387                 if (child_lists[child].next != INVALID) {

   388                   child_lists[child_lists[child].next].prev =

   389                     child_lists[child].prev;

   390                 }

   392                 // Merging external faces

   393                 {

   394                   int en = cn;

   395                   cn = cd ? node_data[cn].prev : node_data[cn].next;

   396                   cd = node_data[cn].next == en;

   398                 }

   400                 if (cd) node_data[cn].next = dn; else node_data[cn].prev = dn;

   401                 if (dd) node_data[dn].prev = cn; else node_data[dn].next = cn;

   403               }

   405               bool d = pn == node_data[n].prev;

   407               if (node_data[n].prev == node_data[n].next &&

   408                   node_data[n].inverted) {

   409                 d = !d;

   410               }

   412               // Embedding arc into external face

   413               if (rd) node_data[rn].next = n; else node_data[rn].prev = n;

   414               if (d) node_data[n].prev = rn; else node_data[n].next = rn;

   415               pn = rn;

   417               embed_arc[order_list[n]] = false;

   418             }

   420             if (!merge_roots[node].empty()) {

   422               bool d = pn == node_data[n].prev;

   424               merge_stack.push_back(std::make_pair(n, d));

   426               int rn = merge_roots[node].front();

   428               int xn = node_data[rn].next;

   429               Node xnode = order_list[xn];

   431               int yn = node_data[rn].prev;

   432               Node ynode = order_list[yn];

   434               bool rd;

   435               if (!external(xnode, rorder, child_lists,

   436                             ancestor_map, low_map)) {

   437                 rd = true;

   438               } else if (!external(ynode, rorder, child_lists,

   439                                    ancestor_map, low_map)) {

   440                 rd = false;

   441               } else if (pertinent(xnode, embed_arc, merge_roots)) {

   442                 rd = true;

   443               } else {

   444                 rd = false;

   445               }

   447               merge_stack.push_back(std::make_pair(rn, rd));

   449               pn = rn;

   450               n = rd ? xn : yn;

   452             } else if (!external(node, rorder, child_lists,

   453                                  ancestor_map, low_map)) {

   454               int nn = (node_data[n].next != pn ?

   455                         node_data[n].next : node_data[n].prev);

   457               bool nd = n == node_data[nn].prev;

   459               if (nd) node_data[nn].prev = pn;

   460               else node_data[nn].next = pn;

   462               if (n == node_data[pn].prev) node_data[pn].prev = nn;

   463               else node_data[pn].next = nn;

   465               node_data[nn].inverted =

   466                 (node_data[nn].prev == node_data[nn].next && nd != rd);

   468               n = nn;

   469             }

   470             else break;

   472           }

   474           if (!merge_stack.empty() || n == rn) {

   475             break;

   476           }

   477         }

   478       }

   480       void initFace(const Node& node, NodeData& node_data,

   481                     const OrderMap& order_map, const OrderList& order_list) {

   482         int n = order_map[node];

   483         int rn = n + order_list.size();

   485         node_data[n].next = node_data[n].prev = rn;

   486         node_data[rn].next = node_data[rn].prev = n;

   488         node_data[n].visited = order_list.size();

   489         node_data[rn].visited = order_list.size();

   491       }

   493       bool external(const Node& node, int rorder,

   494                     ChildLists& child_lists, AncestorMap& ancestor_map,

   495                     LowMap& low_map) {

   496         Node child = child_lists[node].first;

   498         if (child != INVALID) {

   499           if (low_map[child] < rorder) return true;

   500         }

   502         if (ancestor_map[node] < rorder) return true;

   504         return false;

   505       }

   507       bool pertinent(const Node& node, const EmbedArc& embed_arc,

   508                      const MergeRoots& merge_roots) {

   509         return !merge_roots[node].empty() || embed_arc[node];

   510       }

   512     };

   514   }

   516   /// \ingroup planar

   517   ///

   518   /// \brief Planarity checking of an undirected simple graph

   519   ///

   520   /// This function implements the Boyer-Myrvold algorithm for

   521   /// planarity checking of an undirected graph. It is a simplified

   522   /// version of the PlanarEmbedding algorithm class because neither

   523   /// the embedding nor the kuratowski subdivisons are not computed.

   524   template <typename GR>

   525   bool checkPlanarity(const GR& graph) {

   526     _planarity_bits::PlanarityChecking<GR> pc(graph);

   527     return pc.run();

   528   }

   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:

   549     TEMPLATE_GRAPH_TYPEDEFS(Graph);

   551     const Graph& _graph;

   552     typename Graph::template ArcMap<Arc> _embedding;

   554     typename Graph::template EdgeMap<bool> _kuratowski;

   556   private:

   558     typedef typename Graph::template NodeMap<Arc> PredMap;

   560     typedef typename Graph::template EdgeMap<bool> TreeMap;

   562     typedef typename Graph::template NodeMap<int> OrderMap;

   563     typedef std::vector<Node> OrderList;

   565     typedef typename Graph::template NodeMap<int> LowMap;

   566     typedef typename Graph::template NodeMap<int> AncestorMap;

   568     typedef _planarity_bits::NodeDataNode<Graph> NodeDataNode;

   569     typedef std::vector<NodeDataNode> NodeData;

   571     typedef _planarity_bits::ChildListNode<Graph> ChildListNode;

   572     typedef typename Graph::template NodeMap<ChildListNode> ChildLists;

   574     typedef typename Graph::template NodeMap<std::list<int> > MergeRoots;

   576     typedef typename Graph::template NodeMap<Arc> EmbedArc;

   578     typedef _planarity_bits::ArcListNode<Graph> ArcListNode;

   579     typedef typename Graph::template ArcMap<ArcListNode> ArcLists;

   581     typedef typename Graph::template NodeMap<bool> FlipMap;

   583     typedef typename Graph::template NodeMap<int> TypeMap;

   585     enum IsolatorNodeType {

   586       HIGHX = 6, LOWX = 7,

   587       HIGHY = 8, LOWY = 9,

   588       ROOT = 10, PERTINENT = 11,

   589       INTERNAL = 12

   590     };

   592   public:

   594     /// \brief The map for store of embedding

   595     typedef typename Graph::template ArcMap<Arc> EmbeddingMap;

   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) {}

   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;

   613       PredMap pred_map(_graph, INVALID);

   614       TreeMap tree_map(_graph, false);

   616       OrderMap order_map(_graph, -1);

   617       OrderList order_list;

   619       AncestorMap ancestor_map(_graph, -1);

   620       LowMap low_map(_graph, -1);

   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();

   627       ChildLists child_lists(_graph);

   628       createChildLists(tree_map, order_map, low_map, child_lists);

   630       NodeData node_data(2 * order_list.size());

   632       EmbedArc embed_arc(_graph, INVALID);

   634       MergeRoots merge_roots(_graph);

   636       ArcLists arc_lists(_graph);

   638       FlipMap flip_map(_graph, false);

   640       for (int i = order_list.size() - 1; i >= 0; --i) {

   642         Node node = order_list[i];

   644         node_data[i].first = INVALID;

   646         Node source = node;

   647         for (OutArcIt e(_graph, node); e != INVALID; ++e) {

   648           Node target = _graph.target(e);

   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         }

   656         for (OutArcIt e(_graph, node); e != INVALID; ++e) {

   657           Node target = _graph.target(e);

   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         }

   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();

   674         for (OutArcIt e(_graph, node); e != INVALID; ++e) {

   675           Node target = _graph.target(e);

   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       }

   691       for (int i = 0; i < int(order_list.size()); ++i) {

   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       }

   699       return true;

   700     }

   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     }

   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& embeddingMap() const {

   716       return _embedding;

   717     }

   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     }

   729   private:

   731     void createChildLists(const TreeMap& tree_map, const OrderMap& order_map,

   732                           const LowMap& low_map, ChildLists& child_lists) {

   734       for (NodeIt n(_graph); n != INVALID; ++n) {

   735         Node source = n;

   737         std::vector<Node> targets;

   738         for (OutArcIt e(_graph, n); e != INVALID; ++e) {

   739           Node target = _graph.target(e);

   741           if (order_map[source] < order_map[target] && tree_map[e]) {

   742             targets.push_back(target);

   743           }

   744         }

   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     }

   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) {

   770       int na, nb;

   771       bool da, db;

   773       na = nb = order_map[node];

   774       da = true; db = false;

   776       while (true) {

   778         if (node_data[na].visited == rorder) break;

   779         if (node_data[nb].visited == rorder) break;

   781         node_data[na].visited = rorder;

   782         node_data[nb].visited = rorder;

   784         int rn = -1;

   786         if (na >= int(order_list.size())) {

   787           rn = na;

   788         } else if (nb >= int(order_list.size())) {

   789           rn = nb;

   790         }

   792         if (rn == -1) {

   793           int nn;

   795           nn = da ? node_data[na].prev : node_data[na].next;

   796           da = node_data[nn].prev != na;

   797           na = nn;

   799           nn = db ? node_data[nb].prev : node_data[nb].next;

   800           db = node_data[nn].prev != nb;

   801           nb = nn;

   803         } else {

   805           Node rep = order_list[rn - order_list.size()];

   806           Node parent = _graph.source(pred_map[rep]);

   808           if (low_map[rep] < rorder) {

   809             merge_roots[parent].push_back(rn);

   810           } else {

   811             merge_roots[parent].push_front(rn);

   812           }

   814           if (parent != root) {

   815             na = nb = order_map[parent];

   816             da = true; db = false;

   817           } else {

   818             break;

   819           }

   820         }

   821       }

   822     }

   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) {

   830       std::vector<std::pair<int, bool> > merge_stack;

   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;

   837         while (n != rn) {

   839           Node node = order_list[n];

   841           if (embed_arc[node] != INVALID) {

   843             // Merging components on the critical path

   844             while (!merge_stack.empty()) {

   846               // Component root

   847               int cn = merge_stack.back().first;

   848               bool cd = merge_stack.back().second;

   849               merge_stack.pop_back();

   851               // Parent of component

   852               int dn = merge_stack.back().first;

   853               bool dd = merge_stack.back().second;

   854               merge_stack.pop_back();

   856               Node parent = order_list[dn];

   858               // Erasing from merge_roots

   859               merge_roots[parent].pop_front();

   861               Node child = order_list[cn - order_list.size()];

   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               }

   871               if (child_lists[child].next != INVALID) {

   872                 child_lists[child_lists[child].next].prev =

   873                   child_lists[child].prev;

   874               }

   876               // Merging arcs + flipping

   877               Arc de = node_data[dn].first;

   878               Arc ce = node_data[cn].first;

   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               }

   887               {

   888                 Arc dne = arc_lists[de].next;

   889                 Arc cne = arc_lists[ce].next;

   891                 arc_lists[de].next = cne;

   892                 arc_lists[ce].next = dne;

   894                 arc_lists[dne].prev = ce;

   895                 arc_lists[cne].prev = de;

   896               }

   898               if (dd) {

   899                 node_data[dn].first = ce;

   900               }

   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;

   908                  if (node_data[cn].prev == node_data[cn].next &&

   909                     node_data[cn].inverted) {

   910                    cd = !cd;

   911                  }

   912               }

   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;

   917             }

   919             bool d = pn == node_data[n].prev;

   921             if (node_data[n].prev == node_data[n].next &&

   922                 node_data[n].inverted) {

   923               d = !d;

   924             }

   926             // Add new arc

   927             {

   928               Arc arc = embed_arc[node];

   929               Arc re = node_data[rn].first;

   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;

   936               if (!rd) {

   937                 node_data[rn].first = arc;

   938               }

   940               Arc rev = _graph.oppositeArc(arc);

   941               Arc e = node_data[n].first;

   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;

   948               if (d) {

   949                 node_data[n].first = rev;

   950               }

   952             }

   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;

   959             embed_arc[order_list[n]] = INVALID;

   960           }

   962           if (!merge_roots[node].empty()) {

   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             }

   970             merge_stack.push_back(std::make_pair(n, d));

   972             int rn = merge_roots[node].front();

   974             int xn = node_data[rn].next;

   975             Node xnode = order_list[xn];

   977             int yn = node_data[rn].prev;

   978             Node ynode = order_list[yn];

   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             }

   992             merge_stack.push_back(std::make_pair(rn, rd));

   994             pn = rn;

   995             n = rd ? xn : yn;

   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);

  1002             bool nd = n == node_data[nn].prev;

  1004             if (nd) node_data[nn].prev = pn;

  1005             else node_data[nn].next = pn;

  1007             if (n == node_data[pn].prev) node_data[pn].prev = nn;

  1008             else node_data[pn].next = nn;

  1010             node_data[nn].inverted =

  1011               (node_data[nn].prev == node_data[nn].next && nd != rd);

  1013             n = nn;

  1014           }

  1015           else break;

  1017         }

  1019         if (!merge_stack.empty() || n == rn) {

  1020           break;

  1021         }

  1022       }

  1023     }

  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();

  1031       node_data[n].next = node_data[n].prev = rn;

  1032       node_data[rn].next = node_data[rn].prev = n;

  1034       node_data[n].visited = order_list.size();

  1035       node_data[rn].visited = order_list.size();

  1037       node_data[n].inverted = false;

  1038       node_data[rn].inverted = false;

  1040       Arc arc = pred_map[node];

  1041       Arc rev = _graph.oppositeArc(arc);

  1043       node_data[rn].first = arc;

  1044       node_data[n].first = rev;

  1046       arc_lists[arc].prev = arc;

  1047       arc_lists[arc].next = arc;

  1049       arc_lists[rev].prev = rev;

  1050       arc_lists[rev].next = rev;

  1052     }

  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;

  1063         Arc de = node_data[dd].first;

  1064         Arc ce = node_data[cd].first;

  1066         if (de != INVALID) {

  1067           Arc dne = arc_lists[de].next;

  1068           Arc cne = arc_lists[ce].next;

  1070           arc_lists[de].next = cne;

  1071           arc_lists[ce].next = dne;

  1073           arc_lists[dne].prev = ce;

  1074           arc_lists[cne].prev = de;

  1075         }

  1077         node_data[dd].first = ce;

  1079       }

  1080     }

  1082     void storeEmbedding(const Node& node, NodeData& node_data,

  1083                         OrderMap& order_map, PredMap& pred_map,

  1084                         ArcLists& arc_lists, FlipMap& flip_map) {

  1086       if (node_data[order_map[node]].first == INVALID) return;

  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       }

  1093       Arc first = node_data[order_map[node]].first;

  1094       Arc prev = first;

  1096       Arc arc = flip_map[node] ?

  1097         arc_lists[prev].prev : arc_lists[prev].next;

  1099       _embedding[prev] = arc;

  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     }

  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;

  1115       if (child != INVALID) {

  1116         if (low_map[child] < rorder) return true;

  1117       }

  1119       if (ancestor_map[node] < rorder) return true;

  1121       return false;

  1122     }

  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     }

  1129     int lowPoint(const Node& node, OrderMap& order_map, ChildLists& child_lists,

  1130                  AncestorMap& ancestor_map, LowMap& low_map) {

  1131       int low_point;

  1133       Node child = child_lists[node].first;

  1135       if (child != INVALID) {

  1136         low_point = low_map[child];

  1137       } else {

  1138         low_point = order_map[node];

  1139       }

  1141       if (low_point > ancestor_map[node]) {

  1142         low_point = ancestor_map[node];

  1143       }

  1145       return low_point;

  1146     }

  1148     int findComponentRoot(Node root, Node node, ChildLists& child_lists,

  1149                           OrderMap& order_map, OrderList& order_list) {

  1151       int order = order_map[root];

  1152       int norder = order_map[node];

  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     }

  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     }

  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     }

  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) {

  1192       Node cnode = node;

  1193       Node pred = INVALID;

  1195       while (true) {

  1197         bool pert = pertinent(cnode, embed_arc, merge_roots);

  1198         bool ext = external(cnode, rorder, child_lists, ancestor_map, low_map);

  1200         if (pert && ext) {

  1201           if (!merge_roots[cnode].empty()) {

  1202             int cn = merge_roots[cnode].back();

  1204             if (low_map[order_list[cn - order_list.size()]] < rorder) {

  1205               Arc arc = node_data[cn].first;

  1206               _kuratowski.set(arc, true);

  1208               pred = cnode;

  1209               cnode = _graph.target(arc);

  1211               continue;

  1212             }

  1213           }

  1214           wnode = znode = cnode;

  1215           return;

  1217         } else if (pert) {

  1218           wnode = cnode;

  1220           while (!external(cnode, rorder, child_lists, ancestor_map, low_map)) {

  1221             Arc arc = node_data[order_map[cnode]].first;

  1223             if (_graph.target(arc) == pred) {

  1224               arc = arc_lists[arc].next;

  1225             }

  1226             _kuratowski.set(arc, true);

  1228             Node next = _graph.target(arc);

  1229             pred = cnode; cnode = next;

  1230           }

  1232           znode = cnode;

  1233           return;

  1235         } else if (ext) {

  1236           znode = cnode;

  1238           while (!pertinent(cnode, embed_arc, merge_roots)) {

  1239             Arc arc = node_data[order_map[cnode]].first;

  1241             if (_graph.target(arc) == pred) {

  1242               arc = arc_lists[arc].next;

  1243             }

  1244             _kuratowski.set(arc, true);

  1246             Node next = _graph.target(arc);

  1247             pred = cnode; cnode = next;

  1248           }

  1250           wnode = cnode;

  1251           return;

  1253         } else {

  1254           Arc arc = node_data[order_map[cnode]].first;

  1256           if (_graph.target(arc) == pred) {

  1257             arc = arc_lists[arc].next;

  1258           }

  1259           _kuratowski.set(arc, true);

  1261           Node next = _graph.target(arc);

  1262           pred = cnode; cnode = next;

  1263         }

  1265       }

  1267     }

  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;

  1276       std::vector<Node> st, qu;

  1278       st.push_back(root);

  1279       while (!st.empty()) {

  1280         Node node = st.back();

  1281         st.pop_back();

  1282         qu.push_back(node);

  1284         Arc arc = node_data[order_map[node]].first;

  1286         if (type_map[_graph.target(arc)] == 0) {

  1287           st.push_back(_graph.target(arc));

  1288           type_map[_graph.target(arc)] = 1;

  1289         }

  1291         Arc last = arc, pred = arc;

  1292         arc = arc_lists[arc].next;

  1293         while (arc != last) {

  1295           if (type_map[_graph.target(arc)] == 0) {

  1296             st.push_back(_graph.target(arc));

  1297             type_map[_graph.target(arc)] = 1;

  1298           }

  1300           Arc next = arc_lists[arc].next != pred ?

  1301             arc_lists[arc].next : arc_lists[arc].prev;

  1302           pred = arc; arc = next;

  1303         }

  1305       }

  1307       type_map[root] = 2;

  1308       flip_map[root] = false;

  1310       for (int i = 1; i < int(qu.size()); ++i) {

  1312         Node node = qu[i];

  1314         while (type_map[node] != 2) {

  1315           st.push_back(node);

  1316           type_map[node] = 2;

  1317           node = _graph.source(pred_map[node]);

  1318         }

  1320         bool flip = flip_map[node];

  1322         while (!st.empty()) {

  1323           node = st.back();

  1324           st.pop_back();

  1326           flip_map[node] = flip != flip_map[node];

  1327           flip = flip_map[node];

  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       }

  1339       for (int i = 0; i < int(qu.size()); ++i) {

  1341         Arc arc = node_data[order_map[qu[i]]].first;

  1342         Arc last = arc, pred = arc;

  1344         arc = arc_lists[arc].next;

  1345         while (arc != last) {

  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         }

  1353       }

  1354     }

  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);

  1361       while (node != ynode) {

  1362         type_map[node] = HIGHY;

  1363         node = _graph.target(node_data[order_map[node]].first);

  1364       }

  1366       while (node != wnode) {

  1367         type_map[node] = LOWY;

  1368         node = _graph.target(node_data[order_map[node]].first);

  1369       }

  1371       node = _graph.target(node_data[order_map[wnode]].first);

  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;

  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       }

  1385       type_map[wnode] = PERTINENT;

  1386       type_map[root] = ROOT;

  1387     }

  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;

  1395       Node node = wnode;

  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       }

  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;

  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;

  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       }

  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     }

  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     }

  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;

  1448       st.push_back(_graph.oppositeArc(node_data[order_map[root]].first));

  1449       st.push_back(node_data[order_map[root]].first);

  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;

  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;

  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       }

  1474       for (int i = 1; i < int(st.size()); ++i) {

  1475         ppath.push_back(st[i]);

  1476       }

  1477     }

  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);

  1486       if (ancestor_map[node] != lp) {

  1487         node = child_lists[node].first;

  1488         _kuratowski[pred_map[node]] = true;

  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       }

  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       }

  1509       return lp;

  1510     }

  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;

  1519         _kuratowski.set(arc, true);

  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     }

  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     }

  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);

  1549       while (node != xnode) {

  1550         arc = node_data[order_map[node]].first;

  1551         _kuratowski.set(arc, true);

  1552         node = _graph.target(arc);

  1553       }

  1554     }

  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     }

  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     }

  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) {

  1575       Node root = _graph.source(arc);

  1576       Node enode = _graph.target(arc);

  1578       int rorder = order_map[root];

  1580       TypeMap type_map(_graph, 0);

  1582       int rn = findComponentRoot(root, enode, child_lists,

  1583                                  order_map, order_list);

  1585       Node xnode = order_list[node_data[rn].next];

  1586       Node ynode = order_list[node_data[rn].prev];

  1588       // Minor-A

  1589       {

  1590         while (!merge_roots[xnode].empty() || !merge_roots[ynode].empty()) {

  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           }

  1600           xnode = order_list[node_data[rn].next];

  1601           ynode = order_list[node_data[rn].prev];

  1602         }

  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);

  1616           markPertinentPath(lwnode, order_map, node_data, arc_lists,

  1617                             embed_arc, merge_roots);

  1619           return;

  1620         }

  1621       }

  1623       orientComponent(root, rn, order_map, pred_map,

  1624                       node_data, arc_lists, flip_map, type_map);

  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);

  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);

  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);

  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);

  1652           int minlp = xlp < ylp ? xlp : ylp;

  1653           if (zlp < minlp) minlp = zlp;

  1655           int maxlp = xlp > ylp ? xlp : ylp;

  1656           if (zlp > maxlp) maxlp = zlp;

  1658           markPredPath(order_list[maxlp], order_list[minlp], pred_map);

  1660           return;

  1661         }

  1662       }

  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());

  1672       wnode = findPertinent(pynode, order_map, node_data,

  1673                             embed_arc, merge_roots);

  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         }

  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       }

  1707       std::vector<Arc> ppath;

  1708       findPilePath(ppath, root, type_map, order_map, node_data, arc_lists);

  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       }

  1725       // Minor-E*

  1726       {

  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);

  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         }

  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);

  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         }

  1770         markInternalPath(ipath);

  1771         markPertinentPath(wnode, order_map, node_data, arc_lists,

  1772                           embed_arc, merge_roots);

  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         }

  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         }

  1788         if (pynode != ynode) {

  1789           markFacePath(pxnode, wnode, order_map, node_data);

  1791           int minlp = xlp < ylp ? xlp : ylp;

  1792           if (wlp < minlp) minlp = wlp;

  1794           int maxlp = xlp > ylp ? xlp : ylp;

  1795           if (wlp > maxlp) maxlp = wlp;

  1797           markPredPath(order_list[maxlp], order_list[minlp], pred_map);

  1798           return;

  1799         }

  1801         if (pxnode != xnode) {

  1802           markFacePath(wnode, pynode, order_map, node_data);

  1804           int minlp = xlp < ylp ? xlp : ylp;

  1805           if (wlp < minlp) minlp = wlp;

  1807           int maxlp = xlp > ylp ? xlp : ylp;

  1808           if (wlp > maxlp) maxlp = wlp;

  1810           markPredPath(order_list[maxlp], order_list[minlp], pred_map);

  1811           return;

  1812         }

  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       }

  1821     }

  1823   };

  1825   namespace _planarity_bits {

  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();

  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;

  1843             typename Graph::Arc ue = typename Graph::OutArcIt(graph, u);

  1844             typename Graph::Arc ve = typename Graph::OutArcIt(graph, v);

  1846             typename Graph::Arc e = graph.direct(graph.addEdge(u, v), true);

  1848             if (ue != INVALID) {

  1849               embedding[e] = embedding[ue];

  1850               embedding[ue] = e;

  1851             } else {

  1852               embedding[e] = e;

  1853             }

  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     }

  1866     template <typename Graph, typename EmbeddingMap>

  1867     void makeBiNodeConnected(Graph& graph, EmbeddingMap& embedding) {

  1868       typename Graph::template ArcMap<bool> processed(graph);

  1870       std::vector<typename Graph::Arc> arcs;

  1871       for (typename Graph::ArcIt e(graph); e != INVALID; ++e) {

  1872         arcs.push_back(e);

  1873       }

  1875       IterableBoolMap<Graph, typename Graph::Node> visited(graph, false);

  1877       for (int i = 0; i < int(arcs.size()); ++i) {

  1878         typename Graph::Arc pp = arcs[i];

  1879         if (processed[pp]) continue;

  1881         typename Graph::Arc e = embedding[graph.oppositeArc(pp)];

  1882         processed[e] = true;

  1883         visited.set(graph.source(e), true);

  1885         typename Graph::Arc p = e, l = e;

  1886         e = embedding[graph.oppositeArc(e)];

  1888         while (e != l) {

  1889           processed[e] = true;

  1891           if (visited[graph.source(e)]) {

  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;

  1899             embedding[graph.oppositeArc(n)] =

  1900               embedding[graph.oppositeArc(e)];

  1901             embedding[graph.oppositeArc(e)] =

  1902               graph.oppositeArc(n);

  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     }

  1918     template <typename Graph, typename EmbeddingMap>

  1919     void makeMaxPlanar(Graph& graph, EmbeddingMap& embedding) {

  1921       typename Graph::template NodeMap<int> degree(graph);

  1923       for (typename Graph::NodeIt n(graph); n != INVALID; ++n) {

  1924         degree[n] = countIncEdges(graph, n);

  1925       }

  1927       typename Graph::template ArcMap<bool> processed(graph);

  1928       IterableBoolMap<Graph, typename Graph::Node> visited(graph, false);

  1930       std::vector<typename Graph::Arc> arcs;

  1931       for (typename Graph::ArcIt e(graph); e != INVALID; ++e) {

  1932         arcs.push_back(e);

  1933       }

  1935       for (int i = 0; i < int(arcs.size()); ++i) {

  1936         typename Graph::Arc e = arcs[i];

  1938         if (processed[e]) continue;

  1939         processed[e] = true;

  1941         typename Graph::Arc mine = e;

  1942         int mind = degree[graph.source(e)];

  1944         int face_size = 1;

  1946         typename Graph::Arc l = e;

  1947         e = embedding[graph.oppositeArc(e)];

  1948         while (l != e) {

  1949           processed[e] = true;

  1951           ++face_size;

  1953           if (degree[graph.source(e)] < mind) {

  1954             mine = e;

  1955             mind = degree[graph.source(e)];

  1956           }

  1958           e = embedding[graph.oppositeArc(e)];

  1959         }

  1961         if (face_size < 4) {

  1962           continue;

  1963         }

  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         }

  1970         typename Graph::Arc oppe = INVALID;

  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);

  1983         if (oppe == INVALID) {

  1985           e = embedding[graph.oppositeArc(mine)];

  1986           typename Graph::Arc pn = mine, p = e;

  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);

  1993             embedding[n] = pn;

  1994             embedding[graph.oppositeArc(n)] = e;

  1995             embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);

  1997             pn = n;

  1999             p = e;

  2000             e = embedding[graph.oppositeArc(e)];

  2001           }

  2003           embedding[graph.oppositeArc(e)] = pn;

  2005         } else {

  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);

  2012           typename Graph::Arc ce = graph.direct(graph.addEdge(s, t), true);

  2013           embedding[ce] = mine;

  2014           embedding[graph.oppositeArc(ce)] = oppe;

  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);

  2022             embedding[n] = pn;

  2023             embedding[graph.oppositeArc(n)] = e;

  2024             embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);

  2026             pn = n;

  2028             p = e;

  2029             e = embedding[graph.oppositeArc(e)];

  2031           }

  2032           embedding[graph.oppositeArc(e)] = pn;

  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);

  2040             embedding[n] = pn;

  2041             embedding[graph.oppositeArc(n)] = e;

  2042             embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);

  2044             pn = n;

  2046             p = e;

  2047             e = embedding[graph.oppositeArc(e)];

  2049           }

  2050           embedding[graph.oppositeArc(e)] = pn;

  2051         }

  2052       }

  2053     }

  2055   }

  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:

  2073     TEMPLATE_GRAPH_TYPEDEFS(Graph);

  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;

  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) {}

  2088   private:

  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);

  2096       typename AuxGraph::template ArcMap<Arc> prev(graph);

  2098       for (NodeIt n(graph); n != INVALID; ++n) {

  2099         Arc e = OutArcIt(graph, n);

  2101         Arc p = e, l = e;

  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       }

  2112       Node anode, bnode, cnode;

  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       }

  2121       IterableBoolMap<AuxGraph, Node> proper(graph, false);

  2122       typename AuxGraph::template NodeMap<int> conn(graph, -1);

  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;

  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             }

  2148             proper.set(m, conn[m] == 1);

  2149           }

  2150         }

  2151       }

  2154       typename AuxGraph::template ArcMap<int> angle(graph, -1);

  2156       while (proper.trueNum() != 0) {

  2157         Node n = typename IterableBoolMap<AuxGraph, Node>::TrueIt(proper);

  2158         proper.set(n, false);

  2159         conn[n] = -2;

  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             }

  2175             proper.set(m, conn[m] == 1);

  2176           }

  2177         }

  2179         {

  2180           Arc e = OutArcIt(graph, n);

  2181           Arc p = e, l = e;

  2183           e = next[e];

  2184           while (e != l) {

  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             }

  2195             p = e;

  2196             e = next[e];

  2197           }

  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       }

  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);

  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);

  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       }

  2237       cpred[anode] = INVALID;

  2238       cpred[bnode] = INVALID;

  2240       std::vector<Node> aorder, border, corder;

  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       }

  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       }

  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       }

  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       }

  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       }

  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       }

  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       }

  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       }

  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       }

  2365     }

  2367   public:

  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;

  2377       run(pe);

  2378       return true;

  2379     }

  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;

  2391       if (3 * countNodes(_graph) - 6 == countEdges(_graph)) {

  2392         drawing(_graph, embedding, _point_map);

  2393         return;

  2394       }

  2396       AuxGraph aux_graph(_graph);

  2397       typename AuxGraph::template ArcMap<typename AuxGraph::Arc>

  2398         aux_embedding(aux_graph);

  2400       {

  2402         typename Graph::template EdgeMap<typename AuxGraph::Edge>

  2403           ref(_graph);

  2405         for (EdgeIt e(_graph); e != INVALID; ++e) {

  2406           ref[e] = aux_graph.addEdge(_graph.u(e), _graph.v(e));

  2407         }

  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     }

  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     }

  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     }

  2438   private:

  2440     const Graph& _graph;

  2441     PointMap _point_map;

  2443   };

  2445   namespace _planarity_bits {

  2447     template <typename ColorMap>

  2448     class KempeFilter {

  2449     public:

  2450       typedef typename ColorMap::Key Key;

  2451       typedef bool Value;

  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) {}

  2458       Value operator[](const Key& key) const {

  2459         return _color_map[key] == _first || _color_map[key] == _second;

  2460       }

  2462     private:

  2463       const ColorMap& _color_map;

  2464       typename ColorMap::Value _first, _second;

  2465     };

  2466   }

  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:

  2500     TEMPLATE_GRAPH_TYPEDEFS(Graph);

  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;

  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     }

  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     }

  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     }

  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     }

  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     }

  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() {

  2564       typename Graph::template NodeMap<int> heap_index(_graph, -1);

  2565       BucketHeap<typename Graph::template NodeMap<int> > heap(heap_index);

  2567       for (NodeIt n(_graph); n != INVALID; ++n) {

  2568         _color_map[n] = -2;

  2569         heap.push(n, countOutArcs(_graph, n));

  2570       }

  2572       std::vector<Node> order;

  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       }

  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     }

  2608   private:

  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);

  2616       typedef FilterNodes<const Graph, const KempeFilter> KempeGraph;

  2617       KempeGraph kempe_graph(_graph, filter);

  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       }

  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       }

  2635       return true;

  2636     }

  2638     template <typename EmbeddingMap>

  2639     void kempeRecoloring(const Node& node, const EmbeddingMap& embedding) {

  2640       std::vector<Node> nodes;

  2641       nodes.reserve(4);

  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       }

  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     }

  2661   public:

  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) {

  2671       typename Graph::template NodeMap<int> heap_index(_graph, -1);

  2672       BucketHeap<typename Graph::template NodeMap<int> > heap(heap_index);

  2674       for (NodeIt n(_graph); n != INVALID; ++n) {

  2675         _color_map[n] = -2;

  2676         heap.push(n, countOutArcs(_graph, n));

  2677       }

  2679       std::vector<Node> order;

  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       }

  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     }

  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;

  2724       runFiveColoring(pe.embeddingMap());

  2725       return true;

  2726     }

  2728   private:

  2730     const Graph& _graph;

  2731     IndexMap _color_map;

  2732     Palette _palette;

  2733   };

  2735 }

  2737 #endif

author	Peter Kovacs <kpeter@inf.elte.hu>
	Tue, 09 Feb 2010 23:29:51 +0100
changeset 823	a7e93de12cbd
parent 797	30cb42e3e43a
child 828	5fd7fafc4470
permissions	-rw-r--r--