Merge #180 and a bugfix in #51
authorAlpar Juttner <alpar@cs.elte.hu>
Mon, 14 Dec 2009 06:07:52 +0100
changeset 822f903263902f6
parent 821 072ec8120958
parent 803 1b89e29c9fc7
child 823 a7e93de12cbd
child 848 e05b2b48515a
child 1154 43647f48e971
Merge #180 and a bugfix in #51
lemon/Makefile.am
     1.1 --- a/lemon/Makefile.am	Fri Nov 13 00:39:28 2009 +0100
     1.2 +++ b/lemon/Makefile.am	Mon Dec 14 06:07:52 2009 +0100
     1.3 @@ -110,6 +110,7 @@
     1.4  	lemon/network_simplex.h \
     1.5  	lemon/pairing_heap.h \
     1.6  	lemon/path.h \
     1.7 +	lemon/planarity.h \
     1.8  	lemon/preflow.h \
     1.9  	lemon/radix_heap.h \
    1.10  	lemon/radix_sort.h \
     2.1 --- a/lemon/bits/map_extender.h	Fri Nov 13 00:39:28 2009 +0100
     2.2 +++ b/lemon/bits/map_extender.h	Mon Dec 14 06:07:52 2009 +0100
     2.3 @@ -84,36 +84,36 @@
     2.4  
     2.5        typedef typename Map::Value Value;
     2.6  
     2.7 -      MapIt() {}
     2.8 +      MapIt() : map(NULL) {}
     2.9  
    2.10 -      MapIt(Invalid i) : Parent(i) { }
    2.11 +      MapIt(Invalid i) : Parent(i), map(NULL) {}
    2.12  
    2.13 -      explicit MapIt(Map& _map) : map(_map) {
    2.14 -        map.notifier()->first(*this);
    2.15 +      explicit MapIt(Map& _map) : map(&_map) {
    2.16 +        map->notifier()->first(*this);
    2.17        }
    2.18  
    2.19        MapIt(const Map& _map, const Item& item)
    2.20 -        : Parent(item), map(_map) {}
    2.21 +        : Parent(item), map(&_map) {}
    2.22  
    2.23        MapIt& operator++() {
    2.24 -        map.notifier()->next(*this);
    2.25 +        map->notifier()->next(*this);
    2.26          return *this;
    2.27        }
    2.28  
    2.29        typename MapTraits<Map>::ConstReturnValue operator*() const {
    2.30 -        return map[*this];
    2.31 +        return (*map)[*this];
    2.32        }
    2.33  
    2.34        typename MapTraits<Map>::ReturnValue operator*() {
    2.35 -        return map[*this];
    2.36 +        return (*map)[*this];
    2.37        }
    2.38  
    2.39        void set(const Value& value) {
    2.40 -        map.set(*this, value);
    2.41 +        map->set(*this, value);
    2.42        }
    2.43  
    2.44      protected:
    2.45 -      Map& map;
    2.46 +      Map* map;
    2.47  
    2.48      };
    2.49  
    2.50 @@ -124,19 +124,19 @@
    2.51  
    2.52        typedef typename Map::Value Value;
    2.53  
    2.54 -      ConstMapIt() {}
    2.55 +      ConstMapIt() : map(NULL) {}
    2.56  
    2.57 -      ConstMapIt(Invalid i) : Parent(i) { }
    2.58 +      ConstMapIt(Invalid i) : Parent(i), map(NULL) {}
    2.59  
    2.60 -      explicit ConstMapIt(Map& _map) : map(_map) {
    2.61 -        map.notifier()->first(*this);
    2.62 +      explicit ConstMapIt(Map& _map) : map(&_map) {
    2.63 +        map->notifier()->first(*this);
    2.64        }
    2.65  
    2.66        ConstMapIt(const Map& _map, const Item& item)
    2.67          : Parent(item), map(_map) {}
    2.68  
    2.69        ConstMapIt& operator++() {
    2.70 -        map.notifier()->next(*this);
    2.71 +        map->notifier()->next(*this);
    2.72          return *this;
    2.73        }
    2.74  
    2.75 @@ -145,32 +145,32 @@
    2.76        }
    2.77  
    2.78      protected:
    2.79 -      const Map& map;
    2.80 +      const Map* map;
    2.81      };
    2.82  
    2.83      class ItemIt : public Item {
    2.84        typedef Item Parent;
    2.85  
    2.86      public:
    2.87 +      ItemIt() : map(NULL) {}
    2.88  
    2.89 -      ItemIt() {}
    2.90  
    2.91 -      ItemIt(Invalid i) : Parent(i) { }
    2.92 +      ItemIt(Invalid i) : Parent(i), map(NULL) {}
    2.93  
    2.94 -      explicit ItemIt(Map& _map) : map(_map) {
    2.95 -        map.notifier()->first(*this);
    2.96 +      explicit ItemIt(Map& _map) : map(&_map) {
    2.97 +        map->notifier()->first(*this);
    2.98        }
    2.99  
   2.100        ItemIt(const Map& _map, const Item& item)
   2.101 -        : Parent(item), map(_map) {}
   2.102 +        : Parent(item), map(&_map) {}
   2.103  
   2.104        ItemIt& operator++() {
   2.105 -        map.notifier()->next(*this);
   2.106 +        map->notifier()->next(*this);
   2.107          return *this;
   2.108        }
   2.109  
   2.110      protected:
   2.111 -      const Map& map;
   2.112 +      const Map* map;
   2.113  
   2.114      };
   2.115    };
   2.116 @@ -231,36 +231,36 @@
   2.117      public:
   2.118        typedef typename Map::Value Value;
   2.119  
   2.120 -      MapIt() {}
   2.121 +      MapIt() : map(NULL) {}
   2.122  
   2.123 -      MapIt(Invalid i) : Parent(i) { }
   2.124 +      MapIt(Invalid i) : Parent(i), map(NULL) { }
   2.125  
   2.126 -      explicit MapIt(Map& _map) : map(_map) {
   2.127 -        map.graph.first(*this);
   2.128 +      explicit MapIt(Map& _map) : map(&_map) {
   2.129 +        map->graph.first(*this);
   2.130        }
   2.131  
   2.132        MapIt(const Map& _map, const Item& item)
   2.133 -        : Parent(item), map(_map) {}
   2.134 +        : Parent(item), map(&_map) {}
   2.135  
   2.136        MapIt& operator++() {
   2.137 -        map.graph.next(*this);
   2.138 +        map->graph.next(*this);
   2.139          return *this;
   2.140        }
   2.141  
   2.142        typename MapTraits<Map>::ConstReturnValue operator*() const {
   2.143 -        return map[*this];
   2.144 +        return (*map)[*this];
   2.145        }
   2.146  
   2.147        typename MapTraits<Map>::ReturnValue operator*() {
   2.148 -        return map[*this];
   2.149 +        return (*map)[*this];
   2.150        }
   2.151  
   2.152        void set(const Value& value) {
   2.153 -        map.set(*this, value);
   2.154 +        map->set(*this, value);
   2.155        }
   2.156  
   2.157      protected:
   2.158 -      Map& map;
   2.159 +      Map* map;
   2.160  
   2.161      };
   2.162  
   2.163 @@ -271,53 +271,53 @@
   2.164  
   2.165        typedef typename Map::Value Value;
   2.166  
   2.167 -      ConstMapIt() {}
   2.168 +      ConstMapIt() : map(NULL) {}
   2.169  
   2.170 -      ConstMapIt(Invalid i) : Parent(i) { }
   2.171 +      ConstMapIt(Invalid i) : Parent(i), map(NULL) { }
   2.172  
   2.173 -      explicit ConstMapIt(Map& _map) : map(_map) {
   2.174 -        map.graph.first(*this);
   2.175 +      explicit ConstMapIt(Map& _map) : map(&_map) {
   2.176 +        map->graph.first(*this);
   2.177        }
   2.178  
   2.179        ConstMapIt(const Map& _map, const Item& item)
   2.180 -        : Parent(item), map(_map) {}
   2.181 +        : Parent(item), map(&_map) {}
   2.182  
   2.183        ConstMapIt& operator++() {
   2.184 -        map.graph.next(*this);
   2.185 +        map->graph.next(*this);
   2.186          return *this;
   2.187        }
   2.188  
   2.189        typename MapTraits<Map>::ConstReturnValue operator*() const {
   2.190 -        return map[*this];
   2.191 +        return (*map)[*this];
   2.192        }
   2.193  
   2.194      protected:
   2.195 -      const Map& map;
   2.196 +      const Map* map;
   2.197      };
   2.198  
   2.199      class ItemIt : public Item {
   2.200        typedef Item Parent;
   2.201  
   2.202      public:
   2.203 +      ItemIt() : map(NULL) {}
   2.204  
   2.205 -      ItemIt() {}
   2.206  
   2.207 -      ItemIt(Invalid i) : Parent(i) { }
   2.208 +      ItemIt(Invalid i) : Parent(i), map(NULL) { }
   2.209  
   2.210 -      explicit ItemIt(Map& _map) : map(_map) {
   2.211 -        map.graph.first(*this);
   2.212 +      explicit ItemIt(Map& _map) : map(&_map) {
   2.213 +        map->graph.first(*this);
   2.214        }
   2.215  
   2.216        ItemIt(const Map& _map, const Item& item)
   2.217 -        : Parent(item), map(_map) {}
   2.218 +        : Parent(item), map(&_map) {}
   2.219  
   2.220        ItemIt& operator++() {
   2.221 -        map.graph.next(*this);
   2.222 +        map->graph.next(*this);
   2.223          return *this;
   2.224        }
   2.225  
   2.226      protected:
   2.227 -      const Map& map;
   2.228 +      const Map* map;
   2.229  
   2.230      };
   2.231  
     3.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     3.2 +++ b/lemon/planarity.h	Mon Dec 14 06:07:52 2009 +0100
     3.3 @@ -0,0 +1,2737 @@
     3.4 +/* -*- mode: C++; indent-tabs-mode: nil; -*-
     3.5 + *
     3.6 + * This file is a part of LEMON, a generic C++ optimization library.
     3.7 + *
     3.8 + * Copyright (C) 2003-2009
     3.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    3.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    3.11 + *
    3.12 + * Permission to use, modify and distribute this software is granted
    3.13 + * provided that this copyright notice appears in all copies. For
    3.14 + * precise terms see the accompanying LICENSE file.
    3.15 + *
    3.16 + * This software is provided "AS IS" with no warranty of any kind,
    3.17 + * express or implied, and with no claim as to its suitability for any
    3.18 + * purpose.
    3.19 + *
    3.20 + */
    3.21 +
    3.22 +#ifndef LEMON_PLANARITY_H
    3.23 +#define LEMON_PLANARITY_H
    3.24 +
    3.25 +/// \ingroup planar
    3.26 +/// \file
    3.27 +/// \brief Planarity checking, embedding, drawing and coloring
    3.28 +
    3.29 +#include <vector>
    3.30 +#include <list>
    3.31 +
    3.32 +#include <lemon/dfs.h>
    3.33 +#include <lemon/bfs.h>
    3.34 +#include <lemon/radix_sort.h>
    3.35 +#include <lemon/maps.h>
    3.36 +#include <lemon/path.h>
    3.37 +#include <lemon/bucket_heap.h>
    3.38 +#include <lemon/adaptors.h>
    3.39 +#include <lemon/edge_set.h>
    3.40 +#include <lemon/color.h>
    3.41 +#include <lemon/dim2.h>
    3.42 +
    3.43 +namespace lemon {
    3.44 +
    3.45 +  namespace _planarity_bits {
    3.46 +
    3.47 +    template <typename Graph>
    3.48 +    struct PlanarityVisitor : DfsVisitor<Graph> {
    3.49 +
    3.50 +      TEMPLATE_GRAPH_TYPEDEFS(Graph);
    3.51 +
    3.52 +      typedef typename Graph::template NodeMap<Arc> PredMap;
    3.53 +
    3.54 +      typedef typename Graph::template EdgeMap<bool> TreeMap;
    3.55 +
    3.56 +      typedef typename Graph::template NodeMap<int> OrderMap;
    3.57 +      typedef std::vector<Node> OrderList;
    3.58 +
    3.59 +      typedef typename Graph::template NodeMap<int> LowMap;
    3.60 +      typedef typename Graph::template NodeMap<int> AncestorMap;
    3.61 +
    3.62 +      PlanarityVisitor(const Graph& graph,
    3.63 +                       PredMap& pred_map, TreeMap& tree_map,
    3.64 +                       OrderMap& order_map, OrderList& order_list,
    3.65 +                       AncestorMap& ancestor_map, LowMap& low_map)
    3.66 +        : _graph(graph), _pred_map(pred_map), _tree_map(tree_map),
    3.67 +          _order_map(order_map), _order_list(order_list),
    3.68 +          _ancestor_map(ancestor_map), _low_map(low_map) {}
    3.69 +
    3.70 +      void reach(const Node& node) {
    3.71 +        _order_map[node] = _order_list.size();
    3.72 +        _low_map[node] = _order_list.size();
    3.73 +        _ancestor_map[node] = _order_list.size();
    3.74 +        _order_list.push_back(node);
    3.75 +      }
    3.76 +
    3.77 +      void discover(const Arc& arc) {
    3.78 +        Node source = _graph.source(arc);
    3.79 +        Node target = _graph.target(arc);
    3.80 +
    3.81 +        _tree_map[arc] = true;
    3.82 +        _pred_map[target] = arc;
    3.83 +      }
    3.84 +
    3.85 +      void examine(const Arc& arc) {
    3.86 +        Node source = _graph.source(arc);
    3.87 +        Node target = _graph.target(arc);
    3.88 +
    3.89 +        if (_order_map[target] < _order_map[source] && !_tree_map[arc]) {
    3.90 +          if (_low_map[source] > _order_map[target]) {
    3.91 +            _low_map[source] = _order_map[target];
    3.92 +          }
    3.93 +          if (_ancestor_map[source] > _order_map[target]) {
    3.94 +            _ancestor_map[source] = _order_map[target];
    3.95 +          }
    3.96 +        }
    3.97 +      }
    3.98 +
    3.99 +      void backtrack(const Arc& arc) {
   3.100 +        Node source = _graph.source(arc);
   3.101 +        Node target = _graph.target(arc);
   3.102 +
   3.103 +        if (_low_map[source] > _low_map[target]) {
   3.104 +          _low_map[source] = _low_map[target];
   3.105 +        }
   3.106 +      }
   3.107 +
   3.108 +      const Graph& _graph;
   3.109 +      PredMap& _pred_map;
   3.110 +      TreeMap& _tree_map;
   3.111 +      OrderMap& _order_map;
   3.112 +      OrderList& _order_list;
   3.113 +      AncestorMap& _ancestor_map;
   3.114 +      LowMap& _low_map;
   3.115 +    };
   3.116 +
   3.117 +    template <typename Graph, bool embedding = true>
   3.118 +    struct NodeDataNode {
   3.119 +      int prev, next;
   3.120 +      int visited;
   3.121 +      typename Graph::Arc first;
   3.122 +      bool inverted;
   3.123 +    };
   3.124 +
   3.125 +    template <typename Graph>
   3.126 +    struct NodeDataNode<Graph, false> {
   3.127 +      int prev, next;
   3.128 +      int visited;
   3.129 +    };
   3.130 +
   3.131 +    template <typename Graph>
   3.132 +    struct ChildListNode {
   3.133 +      typedef typename Graph::Node Node;
   3.134 +      Node first;
   3.135 +      Node prev, next;
   3.136 +    };
   3.137 +
   3.138 +    template <typename Graph>
   3.139 +    struct ArcListNode {
   3.140 +      typename Graph::Arc prev, next;
   3.141 +    };
   3.142 +
   3.143 +    template <typename Graph>
   3.144 +    class PlanarityChecking {
   3.145 +    private:
   3.146 +      
   3.147 +      TEMPLATE_GRAPH_TYPEDEFS(Graph);
   3.148 +
   3.149 +      const Graph& _graph;
   3.150 +
   3.151 +    private:
   3.152 +      
   3.153 +      typedef typename Graph::template NodeMap<Arc> PredMap;
   3.154 +      
   3.155 +      typedef typename Graph::template EdgeMap<bool> TreeMap;
   3.156 +      
   3.157 +      typedef typename Graph::template NodeMap<int> OrderMap;
   3.158 +      typedef std::vector<Node> OrderList;
   3.159 +
   3.160 +      typedef typename Graph::template NodeMap<int> LowMap;
   3.161 +      typedef typename Graph::template NodeMap<int> AncestorMap;
   3.162 +
   3.163 +      typedef _planarity_bits::NodeDataNode<Graph> NodeDataNode;
   3.164 +      typedef std::vector<NodeDataNode> NodeData;
   3.165 +
   3.166 +      typedef _planarity_bits::ChildListNode<Graph> ChildListNode;
   3.167 +      typedef typename Graph::template NodeMap<ChildListNode> ChildLists;
   3.168 +
   3.169 +      typedef typename Graph::template NodeMap<std::list<int> > MergeRoots;
   3.170 +
   3.171 +      typedef typename Graph::template NodeMap<bool> EmbedArc;
   3.172 +
   3.173 +    public:
   3.174 +
   3.175 +      PlanarityChecking(const Graph& graph) : _graph(graph) {}
   3.176 +
   3.177 +      bool run() {
   3.178 +        typedef _planarity_bits::PlanarityVisitor<Graph> Visitor;
   3.179 +
   3.180 +        PredMap pred_map(_graph, INVALID);
   3.181 +        TreeMap tree_map(_graph, false);
   3.182 +
   3.183 +        OrderMap order_map(_graph, -1);
   3.184 +        OrderList order_list;
   3.185 +
   3.186 +        AncestorMap ancestor_map(_graph, -1);
   3.187 +        LowMap low_map(_graph, -1);
   3.188 +
   3.189 +        Visitor visitor(_graph, pred_map, tree_map,
   3.190 +                        order_map, order_list, ancestor_map, low_map);
   3.191 +        DfsVisit<Graph, Visitor> visit(_graph, visitor);
   3.192 +        visit.run();
   3.193 +
   3.194 +        ChildLists child_lists(_graph);
   3.195 +        createChildLists(tree_map, order_map, low_map, child_lists);
   3.196 +
   3.197 +        NodeData node_data(2 * order_list.size());
   3.198 +
   3.199 +        EmbedArc embed_arc(_graph, false);
   3.200 +
   3.201 +        MergeRoots merge_roots(_graph);
   3.202 +
   3.203 +        for (int i = order_list.size() - 1; i >= 0; --i) {
   3.204 +
   3.205 +          Node node = order_list[i];
   3.206 +
   3.207 +          Node source = node;
   3.208 +          for (OutArcIt e(_graph, node); e != INVALID; ++e) {
   3.209 +            Node target = _graph.target(e);
   3.210 +
   3.211 +            if (order_map[source] < order_map[target] && tree_map[e]) {
   3.212 +              initFace(target, node_data, order_map, order_list);
   3.213 +            }
   3.214 +          }
   3.215 +
   3.216 +          for (OutArcIt e(_graph, node); e != INVALID; ++e) {
   3.217 +            Node target = _graph.target(e);
   3.218 +
   3.219 +            if (order_map[source] < order_map[target] && !tree_map[e]) {
   3.220 +              embed_arc[target] = true;
   3.221 +              walkUp(target, source, i, pred_map, low_map,
   3.222 +                     order_map, order_list, node_data, merge_roots);
   3.223 +            }
   3.224 +          }
   3.225 +
   3.226 +          for (typename MergeRoots::Value::iterator it =
   3.227 +                 merge_roots[node].begin(); 
   3.228 +               it != merge_roots[node].end(); ++it) {
   3.229 +            int rn = *it;
   3.230 +            walkDown(rn, i, node_data, order_list, child_lists,
   3.231 +                     ancestor_map, low_map, embed_arc, merge_roots);
   3.232 +          }
   3.233 +          merge_roots[node].clear();
   3.234 +
   3.235 +          for (OutArcIt e(_graph, node); e != INVALID; ++e) {
   3.236 +            Node target = _graph.target(e);
   3.237 +
   3.238 +            if (order_map[source] < order_map[target] && !tree_map[e]) {
   3.239 +              if (embed_arc[target]) {
   3.240 +                return false;
   3.241 +              }
   3.242 +            }
   3.243 +          }
   3.244 +        }
   3.245 +
   3.246 +        return true;
   3.247 +      }
   3.248 +
   3.249 +    private:
   3.250 +
   3.251 +      void createChildLists(const TreeMap& tree_map, const OrderMap& order_map,
   3.252 +                            const LowMap& low_map, ChildLists& child_lists) {
   3.253 +
   3.254 +        for (NodeIt n(_graph); n != INVALID; ++n) {
   3.255 +          Node source = n;
   3.256 +
   3.257 +          std::vector<Node> targets;
   3.258 +          for (OutArcIt e(_graph, n); e != INVALID; ++e) {
   3.259 +            Node target = _graph.target(e);
   3.260 +
   3.261 +            if (order_map[source] < order_map[target] && tree_map[e]) {
   3.262 +              targets.push_back(target);
   3.263 +            }
   3.264 +          }
   3.265 +
   3.266 +          if (targets.size() == 0) {
   3.267 +            child_lists[source].first = INVALID;
   3.268 +          } else if (targets.size() == 1) {
   3.269 +            child_lists[source].first = targets[0];
   3.270 +            child_lists[targets[0]].prev = INVALID;
   3.271 +            child_lists[targets[0]].next = INVALID;
   3.272 +          } else {
   3.273 +            radixSort(targets.begin(), targets.end(), mapToFunctor(low_map));
   3.274 +            for (int i = 1; i < int(targets.size()); ++i) {
   3.275 +              child_lists[targets[i]].prev = targets[i - 1];
   3.276 +              child_lists[targets[i - 1]].next = targets[i];
   3.277 +            }
   3.278 +            child_lists[targets.back()].next = INVALID;
   3.279 +            child_lists[targets.front()].prev = INVALID;
   3.280 +            child_lists[source].first = targets.front();
   3.281 +          }
   3.282 +        }
   3.283 +      }
   3.284 +
   3.285 +      void walkUp(const Node& node, Node root, int rorder,
   3.286 +                  const PredMap& pred_map, const LowMap& low_map,
   3.287 +                  const OrderMap& order_map, const OrderList& order_list,
   3.288 +                  NodeData& node_data, MergeRoots& merge_roots) {
   3.289 +
   3.290 +        int na, nb;
   3.291 +        bool da, db;
   3.292 +
   3.293 +        na = nb = order_map[node];
   3.294 +        da = true; db = false;
   3.295 +
   3.296 +        while (true) {
   3.297 +
   3.298 +          if (node_data[na].visited == rorder) break;
   3.299 +          if (node_data[nb].visited == rorder) break;
   3.300 +
   3.301 +          node_data[na].visited = rorder;
   3.302 +          node_data[nb].visited = rorder;
   3.303 +
   3.304 +          int rn = -1;
   3.305 +
   3.306 +          if (na >= int(order_list.size())) {
   3.307 +            rn = na;
   3.308 +          } else if (nb >= int(order_list.size())) {
   3.309 +            rn = nb;
   3.310 +          }
   3.311 +
   3.312 +          if (rn == -1) {
   3.313 +            int nn;
   3.314 +
   3.315 +            nn = da ? node_data[na].prev : node_data[na].next;
   3.316 +            da = node_data[nn].prev != na;
   3.317 +            na = nn;
   3.318 +
   3.319 +            nn = db ? node_data[nb].prev : node_data[nb].next;
   3.320 +            db = node_data[nn].prev != nb;
   3.321 +            nb = nn;
   3.322 +
   3.323 +          } else {
   3.324 +
   3.325 +            Node rep = order_list[rn - order_list.size()];
   3.326 +            Node parent = _graph.source(pred_map[rep]);
   3.327 +
   3.328 +            if (low_map[rep] < rorder) {
   3.329 +              merge_roots[parent].push_back(rn);
   3.330 +            } else {
   3.331 +              merge_roots[parent].push_front(rn);
   3.332 +            }
   3.333 +
   3.334 +            if (parent != root) {
   3.335 +              na = nb = order_map[parent];
   3.336 +              da = true; db = false;
   3.337 +            } else {
   3.338 +              break;
   3.339 +            }
   3.340 +          }
   3.341 +        }
   3.342 +      }
   3.343 +
   3.344 +      void walkDown(int rn, int rorder, NodeData& node_data,
   3.345 +                    OrderList& order_list, ChildLists& child_lists,
   3.346 +                    AncestorMap& ancestor_map, LowMap& low_map,
   3.347 +                    EmbedArc& embed_arc, MergeRoots& merge_roots) {
   3.348 +
   3.349 +        std::vector<std::pair<int, bool> > merge_stack;
   3.350 +
   3.351 +        for (int di = 0; di < 2; ++di) {
   3.352 +          bool rd = di == 0;
   3.353 +          int pn = rn;
   3.354 +          int n = rd ? node_data[rn].next : node_data[rn].prev;
   3.355 +
   3.356 +          while (n != rn) {
   3.357 +
   3.358 +            Node node = order_list[n];
   3.359 +
   3.360 +            if (embed_arc[node]) {
   3.361 +
   3.362 +              // Merging components on the critical path
   3.363 +              while (!merge_stack.empty()) {
   3.364 +
   3.365 +                // Component root
   3.366 +                int cn = merge_stack.back().first;
   3.367 +                bool cd = merge_stack.back().second;
   3.368 +                merge_stack.pop_back();
   3.369 +
   3.370 +                // Parent of component
   3.371 +                int dn = merge_stack.back().first;
   3.372 +                bool dd = merge_stack.back().second;
   3.373 +                merge_stack.pop_back();
   3.374 +
   3.375 +                Node parent = order_list[dn];
   3.376 +
   3.377 +                // Erasing from merge_roots
   3.378 +                merge_roots[parent].pop_front();
   3.379 +
   3.380 +                Node child = order_list[cn - order_list.size()];
   3.381 +
   3.382 +                // Erasing from child_lists
   3.383 +                if (child_lists[child].prev != INVALID) {
   3.384 +                  child_lists[child_lists[child].prev].next =
   3.385 +                    child_lists[child].next;
   3.386 +                } else {
   3.387 +                  child_lists[parent].first = child_lists[child].next;
   3.388 +                }
   3.389 +
   3.390 +                if (child_lists[child].next != INVALID) {
   3.391 +                  child_lists[child_lists[child].next].prev =
   3.392 +                    child_lists[child].prev;
   3.393 +                }
   3.394 +
   3.395 +                // Merging external faces
   3.396 +                {
   3.397 +                  int en = cn;
   3.398 +                  cn = cd ? node_data[cn].prev : node_data[cn].next;
   3.399 +                  cd = node_data[cn].next == en;
   3.400 +
   3.401 +                }
   3.402 +
   3.403 +                if (cd) node_data[cn].next = dn; else node_data[cn].prev = dn;
   3.404 +                if (dd) node_data[dn].prev = cn; else node_data[dn].next = cn;
   3.405 +
   3.406 +              }
   3.407 +
   3.408 +              bool d = pn == node_data[n].prev;
   3.409 +
   3.410 +              if (node_data[n].prev == node_data[n].next &&
   3.411 +                  node_data[n].inverted) {
   3.412 +                d = !d;
   3.413 +              }
   3.414 +
   3.415 +              // Embedding arc into external face
   3.416 +              if (rd) node_data[rn].next = n; else node_data[rn].prev = n;
   3.417 +              if (d) node_data[n].prev = rn; else node_data[n].next = rn;
   3.418 +              pn = rn;
   3.419 +
   3.420 +              embed_arc[order_list[n]] = false;
   3.421 +            }
   3.422 +
   3.423 +            if (!merge_roots[node].empty()) {
   3.424 +
   3.425 +              bool d = pn == node_data[n].prev;
   3.426 +
   3.427 +              merge_stack.push_back(std::make_pair(n, d));
   3.428 +
   3.429 +              int rn = merge_roots[node].front();
   3.430 +
   3.431 +              int xn = node_data[rn].next;
   3.432 +              Node xnode = order_list[xn];
   3.433 +
   3.434 +              int yn = node_data[rn].prev;
   3.435 +              Node ynode = order_list[yn];
   3.436 +
   3.437 +              bool rd;
   3.438 +              if (!external(xnode, rorder, child_lists, 
   3.439 +                            ancestor_map, low_map)) {
   3.440 +                rd = true;
   3.441 +              } else if (!external(ynode, rorder, child_lists,
   3.442 +                                   ancestor_map, low_map)) {
   3.443 +                rd = false;
   3.444 +              } else if (pertinent(xnode, embed_arc, merge_roots)) {
   3.445 +                rd = true;
   3.446 +              } else {
   3.447 +                rd = false;
   3.448 +              }
   3.449 +
   3.450 +              merge_stack.push_back(std::make_pair(rn, rd));
   3.451 +
   3.452 +              pn = rn;
   3.453 +              n = rd ? xn : yn;
   3.454 +
   3.455 +            } else if (!external(node, rorder, child_lists,
   3.456 +                                 ancestor_map, low_map)) {
   3.457 +              int nn = (node_data[n].next != pn ?
   3.458 +                        node_data[n].next : node_data[n].prev);
   3.459 +
   3.460 +              bool nd = n == node_data[nn].prev;
   3.461 +
   3.462 +              if (nd) node_data[nn].prev = pn;
   3.463 +              else node_data[nn].next = pn;
   3.464 +
   3.465 +              if (n == node_data[pn].prev) node_data[pn].prev = nn;
   3.466 +              else node_data[pn].next = nn;
   3.467 +
   3.468 +              node_data[nn].inverted =
   3.469 +                (node_data[nn].prev == node_data[nn].next && nd != rd);
   3.470 +
   3.471 +              n = nn;
   3.472 +            }
   3.473 +            else break;
   3.474 +
   3.475 +          }
   3.476 +
   3.477 +          if (!merge_stack.empty() || n == rn) {
   3.478 +            break;
   3.479 +          }
   3.480 +        }
   3.481 +      }
   3.482 +
   3.483 +      void initFace(const Node& node, NodeData& node_data,
   3.484 +                    const OrderMap& order_map, const OrderList& order_list) {
   3.485 +        int n = order_map[node];
   3.486 +        int rn = n + order_list.size();
   3.487 +
   3.488 +        node_data[n].next = node_data[n].prev = rn;
   3.489 +        node_data[rn].next = node_data[rn].prev = n;
   3.490 +
   3.491 +        node_data[n].visited = order_list.size();
   3.492 +        node_data[rn].visited = order_list.size();
   3.493 +
   3.494 +      }
   3.495 +
   3.496 +      bool external(const Node& node, int rorder,
   3.497 +                    ChildLists& child_lists, AncestorMap& ancestor_map,
   3.498 +                    LowMap& low_map) {
   3.499 +        Node child = child_lists[node].first;
   3.500 +
   3.501 +        if (child != INVALID) {
   3.502 +          if (low_map[child] < rorder) return true;
   3.503 +        }
   3.504 +
   3.505 +        if (ancestor_map[node] < rorder) return true;
   3.506 +
   3.507 +        return false;
   3.508 +      }
   3.509 +
   3.510 +      bool pertinent(const Node& node, const EmbedArc& embed_arc,
   3.511 +                     const MergeRoots& merge_roots) {
   3.512 +        return !merge_roots[node].empty() || embed_arc[node];
   3.513 +      }
   3.514 +
   3.515 +    };
   3.516 +
   3.517 +  }
   3.518 +
   3.519 +  /// \ingroup planar
   3.520 +  ///
   3.521 +  /// \brief Planarity checking of an undirected simple graph
   3.522 +  ///
   3.523 +  /// This function implements the Boyer-Myrvold algorithm for
   3.524 +  /// planarity checking of an undirected graph. It is a simplified
   3.525 +  /// version of the PlanarEmbedding algorithm class because neither
   3.526 +  /// the embedding nor the kuratowski subdivisons are not computed.
   3.527 +  template <typename GR>
   3.528 +  bool checkPlanarity(const GR& graph) {
   3.529 +    _planarity_bits::PlanarityChecking<GR> pc(graph);
   3.530 +    return pc.run();
   3.531 +  }
   3.532 +
   3.533 +  /// \ingroup planar
   3.534 +  ///
   3.535 +  /// \brief Planar embedding of an undirected simple graph
   3.536 +  ///
   3.537 +  /// This class implements the Boyer-Myrvold algorithm for planar
   3.538 +  /// embedding of an undirected graph. The planar embedding is an
   3.539 +  /// ordering of the outgoing edges of the nodes, which is a possible
   3.540 +  /// configuration to draw the graph in the plane. If there is not
   3.541 +  /// such ordering then the graph contains a \f$ K_5 \f$ (full graph
   3.542 +  /// with 5 nodes) or a \f$ K_{3,3} \f$ (complete bipartite graph on
   3.543 +  /// 3 ANode and 3 BNode) subdivision.
   3.544 +  ///
   3.545 +  /// The current implementation calculates either an embedding or a
   3.546 +  /// Kuratowski subdivision. The running time of the algorithm is 
   3.547 +  /// \f$ O(n) \f$.
   3.548 +  template <typename Graph>
   3.549 +  class PlanarEmbedding {
   3.550 +  private:
   3.551 +
   3.552 +    TEMPLATE_GRAPH_TYPEDEFS(Graph);
   3.553 +
   3.554 +    const Graph& _graph;
   3.555 +    typename Graph::template ArcMap<Arc> _embedding;
   3.556 +
   3.557 +    typename Graph::template EdgeMap<bool> _kuratowski;
   3.558 +
   3.559 +  private:
   3.560 +
   3.561 +    typedef typename Graph::template NodeMap<Arc> PredMap;
   3.562 +
   3.563 +    typedef typename Graph::template EdgeMap<bool> TreeMap;
   3.564 +
   3.565 +    typedef typename Graph::template NodeMap<int> OrderMap;
   3.566 +    typedef std::vector<Node> OrderList;
   3.567 +
   3.568 +    typedef typename Graph::template NodeMap<int> LowMap;
   3.569 +    typedef typename Graph::template NodeMap<int> AncestorMap;
   3.570 +
   3.571 +    typedef _planarity_bits::NodeDataNode<Graph> NodeDataNode;
   3.572 +    typedef std::vector<NodeDataNode> NodeData;
   3.573 +
   3.574 +    typedef _planarity_bits::ChildListNode<Graph> ChildListNode;
   3.575 +    typedef typename Graph::template NodeMap<ChildListNode> ChildLists;
   3.576 +
   3.577 +    typedef typename Graph::template NodeMap<std::list<int> > MergeRoots;
   3.578 +
   3.579 +    typedef typename Graph::template NodeMap<Arc> EmbedArc;
   3.580 +
   3.581 +    typedef _planarity_bits::ArcListNode<Graph> ArcListNode;
   3.582 +    typedef typename Graph::template ArcMap<ArcListNode> ArcLists;
   3.583 +
   3.584 +    typedef typename Graph::template NodeMap<bool> FlipMap;
   3.585 +
   3.586 +    typedef typename Graph::template NodeMap<int> TypeMap;
   3.587 +
   3.588 +    enum IsolatorNodeType {
   3.589 +      HIGHX = 6, LOWX = 7,
   3.590 +      HIGHY = 8, LOWY = 9,
   3.591 +      ROOT = 10, PERTINENT = 11,
   3.592 +      INTERNAL = 12
   3.593 +    };
   3.594 +
   3.595 +  public:
   3.596 +
   3.597 +    /// \brief The map for store of embedding
   3.598 +    typedef typename Graph::template ArcMap<Arc> EmbeddingMap;
   3.599 +
   3.600 +    /// \brief Constructor
   3.601 +    ///
   3.602 +    /// \note The graph should be simple, i.e. parallel and loop arc
   3.603 +    /// free.
   3.604 +    PlanarEmbedding(const Graph& graph)
   3.605 +      : _graph(graph), _embedding(_graph), _kuratowski(graph, false) {}
   3.606 +
   3.607 +    /// \brief Runs the algorithm.
   3.608 +    ///
   3.609 +    /// Runs the algorithm.
   3.610 +    /// \param kuratowski If the parameter is false, then the
   3.611 +    /// algorithm does not compute a Kuratowski subdivision.
   3.612 +    ///\return %True when the graph is planar.
   3.613 +    bool run(bool kuratowski = true) {
   3.614 +      typedef _planarity_bits::PlanarityVisitor<Graph> Visitor;
   3.615 +
   3.616 +      PredMap pred_map(_graph, INVALID);
   3.617 +      TreeMap tree_map(_graph, false);
   3.618 +
   3.619 +      OrderMap order_map(_graph, -1);
   3.620 +      OrderList order_list;
   3.621 +
   3.622 +      AncestorMap ancestor_map(_graph, -1);
   3.623 +      LowMap low_map(_graph, -1);
   3.624 +
   3.625 +      Visitor visitor(_graph, pred_map, tree_map,
   3.626 +                      order_map, order_list, ancestor_map, low_map);
   3.627 +      DfsVisit<Graph, Visitor> visit(_graph, visitor);
   3.628 +      visit.run();
   3.629 +
   3.630 +      ChildLists child_lists(_graph);
   3.631 +      createChildLists(tree_map, order_map, low_map, child_lists);
   3.632 +
   3.633 +      NodeData node_data(2 * order_list.size());
   3.634 +
   3.635 +      EmbedArc embed_arc(_graph, INVALID);
   3.636 +
   3.637 +      MergeRoots merge_roots(_graph);
   3.638 +
   3.639 +      ArcLists arc_lists(_graph);
   3.640 +
   3.641 +      FlipMap flip_map(_graph, false);
   3.642 +
   3.643 +      for (int i = order_list.size() - 1; i >= 0; --i) {
   3.644 +
   3.645 +        Node node = order_list[i];
   3.646 +
   3.647 +        node_data[i].first = INVALID;
   3.648 +
   3.649 +        Node source = node;
   3.650 +        for (OutArcIt e(_graph, node); e != INVALID; ++e) {
   3.651 +          Node target = _graph.target(e);
   3.652 +
   3.653 +          if (order_map[source] < order_map[target] && tree_map[e]) {
   3.654 +            initFace(target, arc_lists, node_data,
   3.655 +                     pred_map, order_map, order_list);
   3.656 +          }
   3.657 +        }
   3.658 +
   3.659 +        for (OutArcIt e(_graph, node); e != INVALID; ++e) {
   3.660 +          Node target = _graph.target(e);
   3.661 +
   3.662 +          if (order_map[source] < order_map[target] && !tree_map[e]) {
   3.663 +            embed_arc[target] = e;
   3.664 +            walkUp(target, source, i, pred_map, low_map,
   3.665 +                   order_map, order_list, node_data, merge_roots);
   3.666 +          }
   3.667 +        }
   3.668 +
   3.669 +        for (typename MergeRoots::Value::iterator it =
   3.670 +               merge_roots[node].begin(); it != merge_roots[node].end(); ++it) {
   3.671 +          int rn = *it;
   3.672 +          walkDown(rn, i, node_data, arc_lists, flip_map, order_list,
   3.673 +                   child_lists, ancestor_map, low_map, embed_arc, merge_roots);
   3.674 +        }
   3.675 +        merge_roots[node].clear();
   3.676 +
   3.677 +        for (OutArcIt e(_graph, node); e != INVALID; ++e) {
   3.678 +          Node target = _graph.target(e);
   3.679 +
   3.680 +          if (order_map[source] < order_map[target] && !tree_map[e]) {
   3.681 +            if (embed_arc[target] != INVALID) {
   3.682 +              if (kuratowski) {
   3.683 +                isolateKuratowski(e, node_data, arc_lists, flip_map,
   3.684 +                                  order_map, order_list, pred_map, child_lists,
   3.685 +                                  ancestor_map, low_map,
   3.686 +                                  embed_arc, merge_roots);
   3.687 +              }
   3.688 +              return false;
   3.689 +            }
   3.690 +          }
   3.691 +        }
   3.692 +      }
   3.693 +
   3.694 +      for (int i = 0; i < int(order_list.size()); ++i) {
   3.695 +
   3.696 +        mergeRemainingFaces(order_list[i], node_data, order_list, order_map,
   3.697 +                            child_lists, arc_lists);
   3.698 +        storeEmbedding(order_list[i], node_data, order_map, pred_map,
   3.699 +                       arc_lists, flip_map);
   3.700 +      }
   3.701 +
   3.702 +      return true;
   3.703 +    }
   3.704 +
   3.705 +    /// \brief Gives back the successor of an arc
   3.706 +    ///
   3.707 +    /// Gives back the successor of an arc. This function makes
   3.708 +    /// possible to query the cyclic order of the outgoing arcs from
   3.709 +    /// a node.
   3.710 +    Arc next(const Arc& arc) const {
   3.711 +      return _embedding[arc];
   3.712 +    }
   3.713 +
   3.714 +    /// \brief Gives back the calculated embedding map
   3.715 +    ///
   3.716 +    /// The returned map contains the successor of each arc in the
   3.717 +    /// graph.
   3.718 +    const EmbeddingMap& embeddingMap() const {
   3.719 +      return _embedding;
   3.720 +    }
   3.721 +
   3.722 +    /// \brief Gives back true if the undirected arc is in the
   3.723 +    /// kuratowski subdivision
   3.724 +    ///
   3.725 +    /// Gives back true if the undirected arc is in the kuratowski
   3.726 +    /// subdivision
   3.727 +    /// \note The \c run() had to be called with true value.
   3.728 +    bool kuratowski(const Edge& edge) {
   3.729 +      return _kuratowski[edge];
   3.730 +    }
   3.731 +
   3.732 +  private:
   3.733 +
   3.734 +    void createChildLists(const TreeMap& tree_map, const OrderMap& order_map,
   3.735 +                          const LowMap& low_map, ChildLists& child_lists) {
   3.736 +
   3.737 +      for (NodeIt n(_graph); n != INVALID; ++n) {
   3.738 +        Node source = n;
   3.739 +
   3.740 +        std::vector<Node> targets;
   3.741 +        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
   3.742 +          Node target = _graph.target(e);
   3.743 +
   3.744 +          if (order_map[source] < order_map[target] && tree_map[e]) {
   3.745 +            targets.push_back(target);
   3.746 +          }
   3.747 +        }
   3.748 +
   3.749 +        if (targets.size() == 0) {
   3.750 +          child_lists[source].first = INVALID;
   3.751 +        } else if (targets.size() == 1) {
   3.752 +          child_lists[source].first = targets[0];
   3.753 +          child_lists[targets[0]].prev = INVALID;
   3.754 +          child_lists[targets[0]].next = INVALID;
   3.755 +        } else {
   3.756 +          radixSort(targets.begin(), targets.end(), mapToFunctor(low_map));
   3.757 +          for (int i = 1; i < int(targets.size()); ++i) {
   3.758 +            child_lists[targets[i]].prev = targets[i - 1];
   3.759 +            child_lists[targets[i - 1]].next = targets[i];
   3.760 +          }
   3.761 +          child_lists[targets.back()].next = INVALID;
   3.762 +          child_lists[targets.front()].prev = INVALID;
   3.763 +          child_lists[source].first = targets.front();
   3.764 +        }
   3.765 +      }
   3.766 +    }
   3.767 +
   3.768 +    void walkUp(const Node& node, Node root, int rorder,
   3.769 +                const PredMap& pred_map, const LowMap& low_map,
   3.770 +                const OrderMap& order_map, const OrderList& order_list,
   3.771 +                NodeData& node_data, MergeRoots& merge_roots) {
   3.772 +
   3.773 +      int na, nb;
   3.774 +      bool da, db;
   3.775 +
   3.776 +      na = nb = order_map[node];
   3.777 +      da = true; db = false;
   3.778 +
   3.779 +      while (true) {
   3.780 +
   3.781 +        if (node_data[na].visited == rorder) break;
   3.782 +        if (node_data[nb].visited == rorder) break;
   3.783 +
   3.784 +        node_data[na].visited = rorder;
   3.785 +        node_data[nb].visited = rorder;
   3.786 +
   3.787 +        int rn = -1;
   3.788 +
   3.789 +        if (na >= int(order_list.size())) {
   3.790 +          rn = na;
   3.791 +        } else if (nb >= int(order_list.size())) {
   3.792 +          rn = nb;
   3.793 +        }
   3.794 +
   3.795 +        if (rn == -1) {
   3.796 +          int nn;
   3.797 +
   3.798 +          nn = da ? node_data[na].prev : node_data[na].next;
   3.799 +          da = node_data[nn].prev != na;
   3.800 +          na = nn;
   3.801 +
   3.802 +          nn = db ? node_data[nb].prev : node_data[nb].next;
   3.803 +          db = node_data[nn].prev != nb;
   3.804 +          nb = nn;
   3.805 +
   3.806 +        } else {
   3.807 +
   3.808 +          Node rep = order_list[rn - order_list.size()];
   3.809 +          Node parent = _graph.source(pred_map[rep]);
   3.810 +
   3.811 +          if (low_map[rep] < rorder) {
   3.812 +            merge_roots[parent].push_back(rn);
   3.813 +          } else {
   3.814 +            merge_roots[parent].push_front(rn);
   3.815 +          }
   3.816 +
   3.817 +          if (parent != root) {
   3.818 +            na = nb = order_map[parent];
   3.819 +            da = true; db = false;
   3.820 +          } else {
   3.821 +            break;
   3.822 +          }
   3.823 +        }
   3.824 +      }
   3.825 +    }
   3.826 +
   3.827 +    void walkDown(int rn, int rorder, NodeData& node_data,
   3.828 +                  ArcLists& arc_lists, FlipMap& flip_map,
   3.829 +                  OrderList& order_list, ChildLists& child_lists,
   3.830 +                  AncestorMap& ancestor_map, LowMap& low_map,
   3.831 +                  EmbedArc& embed_arc, MergeRoots& merge_roots) {
   3.832 +
   3.833 +      std::vector<std::pair<int, bool> > merge_stack;
   3.834 +
   3.835 +      for (int di = 0; di < 2; ++di) {
   3.836 +        bool rd = di == 0;
   3.837 +        int pn = rn;
   3.838 +        int n = rd ? node_data[rn].next : node_data[rn].prev;
   3.839 +
   3.840 +        while (n != rn) {
   3.841 +
   3.842 +          Node node = order_list[n];
   3.843 +
   3.844 +          if (embed_arc[node] != INVALID) {
   3.845 +
   3.846 +            // Merging components on the critical path
   3.847 +            while (!merge_stack.empty()) {
   3.848 +
   3.849 +              // Component root
   3.850 +              int cn = merge_stack.back().first;
   3.851 +              bool cd = merge_stack.back().second;
   3.852 +              merge_stack.pop_back();
   3.853 +
   3.854 +              // Parent of component
   3.855 +              int dn = merge_stack.back().first;
   3.856 +              bool dd = merge_stack.back().second;
   3.857 +              merge_stack.pop_back();
   3.858 +
   3.859 +              Node parent = order_list[dn];
   3.860 +
   3.861 +              // Erasing from merge_roots
   3.862 +              merge_roots[parent].pop_front();
   3.863 +
   3.864 +              Node child = order_list[cn - order_list.size()];
   3.865 +
   3.866 +              // Erasing from child_lists
   3.867 +              if (child_lists[child].prev != INVALID) {
   3.868 +                child_lists[child_lists[child].prev].next =
   3.869 +                  child_lists[child].next;
   3.870 +              } else {
   3.871 +                child_lists[parent].first = child_lists[child].next;
   3.872 +              }
   3.873 +
   3.874 +              if (child_lists[child].next != INVALID) {
   3.875 +                child_lists[child_lists[child].next].prev =
   3.876 +                  child_lists[child].prev;
   3.877 +              }
   3.878 +
   3.879 +              // Merging arcs + flipping
   3.880 +              Arc de = node_data[dn].first;
   3.881 +              Arc ce = node_data[cn].first;
   3.882 +
   3.883 +              flip_map[order_list[cn - order_list.size()]] = cd != dd;
   3.884 +              if (cd != dd) {
   3.885 +                std::swap(arc_lists[ce].prev, arc_lists[ce].next);
   3.886 +                ce = arc_lists[ce].prev;
   3.887 +                std::swap(arc_lists[ce].prev, arc_lists[ce].next);
   3.888 +              }
   3.889 +
   3.890 +              {
   3.891 +                Arc dne = arc_lists[de].next;
   3.892 +                Arc cne = arc_lists[ce].next;
   3.893 +
   3.894 +                arc_lists[de].next = cne;
   3.895 +                arc_lists[ce].next = dne;
   3.896 +
   3.897 +                arc_lists[dne].prev = ce;
   3.898 +                arc_lists[cne].prev = de;
   3.899 +              }
   3.900 +
   3.901 +              if (dd) {
   3.902 +                node_data[dn].first = ce;
   3.903 +              }
   3.904 +
   3.905 +              // Merging external faces
   3.906 +              {
   3.907 +                int en = cn;
   3.908 +                cn = cd ? node_data[cn].prev : node_data[cn].next;
   3.909 +                cd = node_data[cn].next == en;
   3.910 +
   3.911 +                 if (node_data[cn].prev == node_data[cn].next &&
   3.912 +                    node_data[cn].inverted) {
   3.913 +                   cd = !cd;
   3.914 +                 }
   3.915 +              }
   3.916 +
   3.917 +              if (cd) node_data[cn].next = dn; else node_data[cn].prev = dn;
   3.918 +              if (dd) node_data[dn].prev = cn; else node_data[dn].next = cn;
   3.919 +
   3.920 +            }
   3.921 +
   3.922 +            bool d = pn == node_data[n].prev;
   3.923 +
   3.924 +            if (node_data[n].prev == node_data[n].next &&
   3.925 +                node_data[n].inverted) {
   3.926 +              d = !d;
   3.927 +            }
   3.928 +
   3.929 +            // Add new arc
   3.930 +            {
   3.931 +              Arc arc = embed_arc[node];
   3.932 +              Arc re = node_data[rn].first;
   3.933 +
   3.934 +              arc_lists[arc_lists[re].next].prev = arc;
   3.935 +              arc_lists[arc].next = arc_lists[re].next;
   3.936 +              arc_lists[arc].prev = re;
   3.937 +              arc_lists[re].next = arc;
   3.938 +
   3.939 +              if (!rd) {
   3.940 +                node_data[rn].first = arc;
   3.941 +              }
   3.942 +
   3.943 +              Arc rev = _graph.oppositeArc(arc);
   3.944 +              Arc e = node_data[n].first;
   3.945 +
   3.946 +              arc_lists[arc_lists[e].next].prev = rev;
   3.947 +              arc_lists[rev].next = arc_lists[e].next;
   3.948 +              arc_lists[rev].prev = e;
   3.949 +              arc_lists[e].next = rev;
   3.950 +
   3.951 +              if (d) {
   3.952 +                node_data[n].first = rev;
   3.953 +              }
   3.954 +
   3.955 +            }
   3.956 +
   3.957 +            // Embedding arc into external face
   3.958 +            if (rd) node_data[rn].next = n; else node_data[rn].prev = n;
   3.959 +            if (d) node_data[n].prev = rn; else node_data[n].next = rn;
   3.960 +            pn = rn;
   3.961 +
   3.962 +            embed_arc[order_list[n]] = INVALID;
   3.963 +          }
   3.964 +
   3.965 +          if (!merge_roots[node].empty()) {
   3.966 +
   3.967 +            bool d = pn == node_data[n].prev;
   3.968 +            if (node_data[n].prev == node_data[n].next &&
   3.969 +                node_data[n].inverted) {
   3.970 +              d = !d;
   3.971 +            }
   3.972 +
   3.973 +            merge_stack.push_back(std::make_pair(n, d));
   3.974 +
   3.975 +            int rn = merge_roots[node].front();
   3.976 +
   3.977 +            int xn = node_data[rn].next;
   3.978 +            Node xnode = order_list[xn];
   3.979 +
   3.980 +            int yn = node_data[rn].prev;
   3.981 +            Node ynode = order_list[yn];
   3.982 +
   3.983 +            bool rd;
   3.984 +            if (!external(xnode, rorder, child_lists, ancestor_map, low_map)) {
   3.985 +              rd = true;
   3.986 +            } else if (!external(ynode, rorder, child_lists,
   3.987 +                                 ancestor_map, low_map)) {
   3.988 +              rd = false;
   3.989 +            } else if (pertinent(xnode, embed_arc, merge_roots)) {
   3.990 +              rd = true;
   3.991 +            } else {
   3.992 +              rd = false;
   3.993 +            }
   3.994 +
   3.995 +            merge_stack.push_back(std::make_pair(rn, rd));
   3.996 +
   3.997 +            pn = rn;
   3.998 +            n = rd ? xn : yn;
   3.999 +
  3.1000 +          } else if (!external(node, rorder, child_lists,
  3.1001 +                               ancestor_map, low_map)) {
  3.1002 +            int nn = (node_data[n].next != pn ?
  3.1003 +                      node_data[n].next : node_data[n].prev);
  3.1004 +
  3.1005 +            bool nd = n == node_data[nn].prev;
  3.1006 +
  3.1007 +            if (nd) node_data[nn].prev = pn;
  3.1008 +            else node_data[nn].next = pn;
  3.1009 +
  3.1010 +            if (n == node_data[pn].prev) node_data[pn].prev = nn;
  3.1011 +            else node_data[pn].next = nn;
  3.1012 +
  3.1013 +            node_data[nn].inverted =
  3.1014 +              (node_data[nn].prev == node_data[nn].next && nd != rd);
  3.1015 +
  3.1016 +            n = nn;
  3.1017 +          }
  3.1018 +          else break;
  3.1019 +
  3.1020 +        }
  3.1021 +
  3.1022 +        if (!merge_stack.empty() || n == rn) {
  3.1023 +          break;
  3.1024 +        }
  3.1025 +      }
  3.1026 +    }
  3.1027 +
  3.1028 +    void initFace(const Node& node, ArcLists& arc_lists,
  3.1029 +                  NodeData& node_data, const PredMap& pred_map,
  3.1030 +                  const OrderMap& order_map, const OrderList& order_list) {
  3.1031 +      int n = order_map[node];
  3.1032 +      int rn = n + order_list.size();
  3.1033 +
  3.1034 +      node_data[n].next = node_data[n].prev = rn;
  3.1035 +      node_data[rn].next = node_data[rn].prev = n;
  3.1036 +
  3.1037 +      node_data[n].visited = order_list.size();
  3.1038 +      node_data[rn].visited = order_list.size();
  3.1039 +
  3.1040 +      node_data[n].inverted = false;
  3.1041 +      node_data[rn].inverted = false;
  3.1042 +
  3.1043 +      Arc arc = pred_map[node];
  3.1044 +      Arc rev = _graph.oppositeArc(arc);
  3.1045 +
  3.1046 +      node_data[rn].first = arc;
  3.1047 +      node_data[n].first = rev;
  3.1048 +
  3.1049 +      arc_lists[arc].prev = arc;
  3.1050 +      arc_lists[arc].next = arc;
  3.1051 +
  3.1052 +      arc_lists[rev].prev = rev;
  3.1053 +      arc_lists[rev].next = rev;
  3.1054 +
  3.1055 +    }
  3.1056 +
  3.1057 +    void mergeRemainingFaces(const Node& node, NodeData& node_data,
  3.1058 +                             OrderList& order_list, OrderMap& order_map,
  3.1059 +                             ChildLists& child_lists, ArcLists& arc_lists) {
  3.1060 +      while (child_lists[node].first != INVALID) {
  3.1061 +        int dd = order_map[node];
  3.1062 +        Node child = child_lists[node].first;
  3.1063 +        int cd = order_map[child] + order_list.size();
  3.1064 +        child_lists[node].first = child_lists[child].next;
  3.1065 +
  3.1066 +        Arc de = node_data[dd].first;
  3.1067 +        Arc ce = node_data[cd].first;
  3.1068 +
  3.1069 +        if (de != INVALID) {
  3.1070 +          Arc dne = arc_lists[de].next;
  3.1071 +          Arc cne = arc_lists[ce].next;
  3.1072 +
  3.1073 +          arc_lists[de].next = cne;
  3.1074 +          arc_lists[ce].next = dne;
  3.1075 +
  3.1076 +          arc_lists[dne].prev = ce;
  3.1077 +          arc_lists[cne].prev = de;
  3.1078 +        }
  3.1079 +
  3.1080 +        node_data[dd].first = ce;
  3.1081 +
  3.1082 +      }
  3.1083 +    }
  3.1084 +
  3.1085 +    void storeEmbedding(const Node& node, NodeData& node_data,
  3.1086 +                        OrderMap& order_map, PredMap& pred_map,
  3.1087 +                        ArcLists& arc_lists, FlipMap& flip_map) {
  3.1088 +
  3.1089 +      if (node_data[order_map[node]].first == INVALID) return;
  3.1090 +
  3.1091 +      if (pred_map[node] != INVALID) {
  3.1092 +        Node source = _graph.source(pred_map[node]);
  3.1093 +        flip_map[node] = flip_map[node] != flip_map[source];
  3.1094 +      }
  3.1095 +
  3.1096 +      Arc first = node_data[order_map[node]].first;
  3.1097 +      Arc prev = first;
  3.1098 +
  3.1099 +      Arc arc = flip_map[node] ?
  3.1100 +        arc_lists[prev].prev : arc_lists[prev].next;
  3.1101 +
  3.1102 +      _embedding[prev] = arc;
  3.1103 +
  3.1104 +      while (arc != first) {
  3.1105 +        Arc next = arc_lists[arc].prev == prev ?
  3.1106 +          arc_lists[arc].next : arc_lists[arc].prev;
  3.1107 +        prev = arc; arc = next;
  3.1108 +        _embedding[prev] = arc;
  3.1109 +      }
  3.1110 +    }
  3.1111 +
  3.1112 +
  3.1113 +    bool external(const Node& node, int rorder,
  3.1114 +                  ChildLists& child_lists, AncestorMap& ancestor_map,
  3.1115 +                  LowMap& low_map) {
  3.1116 +      Node child = child_lists[node].first;
  3.1117 +
  3.1118 +      if (child != INVALID) {
  3.1119 +        if (low_map[child] < rorder) return true;
  3.1120 +      }
  3.1121 +
  3.1122 +      if (ancestor_map[node] < rorder) return true;
  3.1123 +
  3.1124 +      return false;
  3.1125 +    }
  3.1126 +
  3.1127 +    bool pertinent(const Node& node, const EmbedArc& embed_arc,
  3.1128 +                   const MergeRoots& merge_roots) {
  3.1129 +      return !merge_roots[node].empty() || embed_arc[node] != INVALID;
  3.1130 +    }
  3.1131 +
  3.1132 +    int lowPoint(const Node& node, OrderMap& order_map, ChildLists& child_lists,
  3.1133 +                 AncestorMap& ancestor_map, LowMap& low_map) {
  3.1134 +      int low_point;
  3.1135 +
  3.1136 +      Node child = child_lists[node].first;
  3.1137 +
  3.1138 +      if (child != INVALID) {
  3.1139 +        low_point = low_map[child];
  3.1140 +      } else {
  3.1141 +        low_point = order_map[node];
  3.1142 +      }
  3.1143 +
  3.1144 +      if (low_point > ancestor_map[node]) {
  3.1145 +        low_point = ancestor_map[node];
  3.1146 +      }
  3.1147 +
  3.1148 +      return low_point;
  3.1149 +    }
  3.1150 +
  3.1151 +    int findComponentRoot(Node root, Node node, ChildLists& child_lists,
  3.1152 +                          OrderMap& order_map, OrderList& order_list) {
  3.1153 +
  3.1154 +      int order = order_map[root];
  3.1155 +      int norder = order_map[node];
  3.1156 +
  3.1157 +      Node child = child_lists[root].first;
  3.1158 +      while (child != INVALID) {
  3.1159 +        int corder = order_map[child];
  3.1160 +        if (corder > order && corder < norder) {
  3.1161 +          order = corder;
  3.1162 +        }
  3.1163 +        child = child_lists[child].next;
  3.1164 +      }
  3.1165 +      return order + order_list.size();
  3.1166 +    }
  3.1167 +
  3.1168 +    Node findPertinent(Node node, OrderMap& order_map, NodeData& node_data,
  3.1169 +                       EmbedArc& embed_arc, MergeRoots& merge_roots) {
  3.1170 +      Node wnode =_graph.target(node_data[order_map[node]].first);
  3.1171 +      while (!pertinent(wnode, embed_arc, merge_roots)) {
  3.1172 +        wnode = _graph.target(node_data[order_map[wnode]].first);
  3.1173 +      }
  3.1174 +      return wnode;
  3.1175 +    }
  3.1176 +
  3.1177 +
  3.1178 +    Node findExternal(Node node, int rorder, OrderMap& order_map,
  3.1179 +                      ChildLists& child_lists, AncestorMap& ancestor_map,
  3.1180 +                      LowMap& low_map, NodeData& node_data) {
  3.1181 +      Node wnode =_graph.target(node_data[order_map[node]].first);
  3.1182 +      while (!external(wnode, rorder, child_lists, ancestor_map, low_map)) {
  3.1183 +        wnode = _graph.target(node_data[order_map[wnode]].first);
  3.1184 +      }
  3.1185 +      return wnode;
  3.1186 +    }
  3.1187 +
  3.1188 +    void markCommonPath(Node node, int rorder, Node& wnode, Node& znode,
  3.1189 +                        OrderList& order_list, OrderMap& order_map,
  3.1190 +                        NodeData& node_data, ArcLists& arc_lists,
  3.1191 +                        EmbedArc& embed_arc, MergeRoots& merge_roots,
  3.1192 +                        ChildLists& child_lists, AncestorMap& ancestor_map,
  3.1193 +                        LowMap& low_map) {
  3.1194 +
  3.1195 +      Node cnode = node;
  3.1196 +      Node pred = INVALID;
  3.1197 +
  3.1198 +      while (true) {
  3.1199 +
  3.1200 +        bool pert = pertinent(cnode, embed_arc, merge_roots);
  3.1201 +        bool ext = external(cnode, rorder, child_lists, ancestor_map, low_map);
  3.1202 +
  3.1203 +        if (pert && ext) {
  3.1204 +          if (!merge_roots[cnode].empty()) {
  3.1205 +            int cn = merge_roots[cnode].back();
  3.1206 +
  3.1207 +            if (low_map[order_list[cn - order_list.size()]] < rorder) {
  3.1208 +              Arc arc = node_data[cn].first;
  3.1209 +              _kuratowski.set(arc, true);
  3.1210 +
  3.1211 +              pred = cnode;
  3.1212 +              cnode = _graph.target(arc);
  3.1213 +
  3.1214 +              continue;
  3.1215 +            }
  3.1216 +          }
  3.1217 +          wnode = znode = cnode;
  3.1218 +          return;
  3.1219 +
  3.1220 +        } else if (pert) {
  3.1221 +          wnode = cnode;
  3.1222 +
  3.1223 +          while (!external(cnode, rorder, child_lists, ancestor_map, low_map)) {
  3.1224 +            Arc arc = node_data[order_map[cnode]].first;
  3.1225 +
  3.1226 +            if (_graph.target(arc) == pred) {
  3.1227 +              arc = arc_lists[arc].next;
  3.1228 +            }
  3.1229 +            _kuratowski.set(arc, true);
  3.1230 +
  3.1231 +            Node next = _graph.target(arc);
  3.1232 +            pred = cnode; cnode = next;
  3.1233 +          }
  3.1234 +
  3.1235 +          znode = cnode;
  3.1236 +          return;
  3.1237 +
  3.1238 +        } else if (ext) {
  3.1239 +          znode = cnode;
  3.1240 +
  3.1241 +          while (!pertinent(cnode, embed_arc, merge_roots)) {
  3.1242 +            Arc arc = node_data[order_map[cnode]].first;
  3.1243 +
  3.1244 +            if (_graph.target(arc) == pred) {
  3.1245 +              arc = arc_lists[arc].next;
  3.1246 +            }
  3.1247 +            _kuratowski.set(arc, true);
  3.1248 +
  3.1249 +            Node next = _graph.target(arc);
  3.1250 +            pred = cnode; cnode = next;
  3.1251 +          }
  3.1252 +
  3.1253 +          wnode = cnode;
  3.1254 +          return;
  3.1255 +
  3.1256 +        } else {
  3.1257 +          Arc arc = node_data[order_map[cnode]].first;
  3.1258 +
  3.1259 +          if (_graph.target(arc) == pred) {
  3.1260 +            arc = arc_lists[arc].next;
  3.1261 +          }
  3.1262 +          _kuratowski.set(arc, true);
  3.1263 +
  3.1264 +          Node next = _graph.target(arc);
  3.1265 +          pred = cnode; cnode = next;
  3.1266 +        }
  3.1267 +
  3.1268 +      }
  3.1269 +
  3.1270 +    }
  3.1271 +
  3.1272 +    void orientComponent(Node root, int rn, OrderMap& order_map,
  3.1273 +                         PredMap& pred_map, NodeData& node_data,
  3.1274 +                         ArcLists& arc_lists, FlipMap& flip_map,
  3.1275 +                         TypeMap& type_map) {
  3.1276 +      node_data[order_map[root]].first = node_data[rn].first;
  3.1277 +      type_map[root] = 1;
  3.1278 +
  3.1279 +      std::vector<Node> st, qu;
  3.1280 +
  3.1281 +      st.push_back(root);
  3.1282 +      while (!st.empty()) {
  3.1283 +        Node node = st.back();
  3.1284 +        st.pop_back();
  3.1285 +        qu.push_back(node);
  3.1286 +
  3.1287 +        Arc arc = node_data[order_map[node]].first;
  3.1288 +
  3.1289 +        if (type_map[_graph.target(arc)] == 0) {
  3.1290 +          st.push_back(_graph.target(arc));
  3.1291 +          type_map[_graph.target(arc)] = 1;
  3.1292 +        }
  3.1293 +
  3.1294 +        Arc last = arc, pred = arc;
  3.1295 +        arc = arc_lists[arc].next;
  3.1296 +        while (arc != last) {
  3.1297 +
  3.1298 +          if (type_map[_graph.target(arc)] == 0) {
  3.1299 +            st.push_back(_graph.target(arc));
  3.1300 +            type_map[_graph.target(arc)] = 1;
  3.1301 +          }
  3.1302 +
  3.1303 +          Arc next = arc_lists[arc].next != pred ?
  3.1304 +            arc_lists[arc].next : arc_lists[arc].prev;
  3.1305 +          pred = arc; arc = next;
  3.1306 +        }
  3.1307 +
  3.1308 +      }
  3.1309 +
  3.1310 +      type_map[root] = 2;
  3.1311 +      flip_map[root] = false;
  3.1312 +
  3.1313 +      for (int i = 1; i < int(qu.size()); ++i) {
  3.1314 +
  3.1315 +        Node node = qu[i];
  3.1316 +
  3.1317 +        while (type_map[node] != 2) {
  3.1318 +          st.push_back(node);
  3.1319 +          type_map[node] = 2;
  3.1320 +          node = _graph.source(pred_map[node]);
  3.1321 +        }
  3.1322 +
  3.1323 +        bool flip = flip_map[node];
  3.1324 +
  3.1325 +        while (!st.empty()) {
  3.1326 +          node = st.back();
  3.1327 +          st.pop_back();
  3.1328 +
  3.1329 +          flip_map[node] = flip != flip_map[node];
  3.1330 +          flip = flip_map[node];
  3.1331 +
  3.1332 +          if (flip) {
  3.1333 +            Arc arc = node_data[order_map[node]].first;
  3.1334 +            std::swap(arc_lists[arc].prev, arc_lists[arc].next);
  3.1335 +            arc = arc_lists[arc].prev;
  3.1336 +            std::swap(arc_lists[arc].prev, arc_lists[arc].next);
  3.1337 +            node_data[order_map[node]].first = arc;
  3.1338 +          }
  3.1339 +        }
  3.1340 +      }
  3.1341 +
  3.1342 +      for (int i = 0; i < int(qu.size()); ++i) {
  3.1343 +
  3.1344 +        Arc arc = node_data[order_map[qu[i]]].first;
  3.1345 +        Arc last = arc, pred = arc;
  3.1346 +
  3.1347 +        arc = arc_lists[arc].next;
  3.1348 +        while (arc != last) {
  3.1349 +
  3.1350 +          if (arc_lists[arc].next == pred) {
  3.1351 +            std::swap(arc_lists[arc].next, arc_lists[arc].prev);
  3.1352 +          }
  3.1353 +          pred = arc; arc = arc_lists[arc].next;
  3.1354 +        }
  3.1355 +
  3.1356 +      }
  3.1357 +    }
  3.1358 +
  3.1359 +    void setFaceFlags(Node root, Node wnode, Node ynode, Node xnode,
  3.1360 +                      OrderMap& order_map, NodeData& node_data,
  3.1361 +                      TypeMap& type_map) {
  3.1362 +      Node node = _graph.target(node_data[order_map[root]].first);
  3.1363 +
  3.1364 +      while (node != ynode) {
  3.1365 +        type_map[node] = HIGHY;
  3.1366 +        node = _graph.target(node_data[order_map[node]].first);
  3.1367 +      }
  3.1368 +
  3.1369 +      while (node != wnode) {
  3.1370 +        type_map[node] = LOWY;
  3.1371 +        node = _graph.target(node_data[order_map[node]].first);
  3.1372 +      }
  3.1373 +
  3.1374 +      node = _graph.target(node_data[order_map[wnode]].first);
  3.1375 +
  3.1376 +      while (node != xnode) {
  3.1377 +        type_map[node] = LOWX;
  3.1378 +        node = _graph.target(node_data[order_map[node]].first);
  3.1379 +      }
  3.1380 +      type_map[node] = LOWX;
  3.1381 +
  3.1382 +      node = _graph.target(node_data[order_map[xnode]].first);
  3.1383 +      while (node != root) {
  3.1384 +        type_map[node] = HIGHX;
  3.1385 +        node = _graph.target(node_data[order_map[node]].first);
  3.1386 +      }
  3.1387 +
  3.1388 +      type_map[wnode] = PERTINENT;
  3.1389 +      type_map[root] = ROOT;
  3.1390 +    }
  3.1391 +
  3.1392 +    void findInternalPath(std::vector<Arc>& ipath,
  3.1393 +                          Node wnode, Node root, TypeMap& type_map,
  3.1394 +                          OrderMap& order_map, NodeData& node_data,
  3.1395 +                          ArcLists& arc_lists) {
  3.1396 +      std::vector<Arc> st;
  3.1397 +
  3.1398 +      Node node = wnode;
  3.1399 +
  3.1400 +      while (node != root) {
  3.1401 +        Arc arc = arc_lists[node_data[order_map[node]].first].next;
  3.1402 +        st.push_back(arc);
  3.1403 +        node = _graph.target(arc);
  3.1404 +      }
  3.1405 +
  3.1406 +      while (true) {
  3.1407 +        Arc arc = st.back();
  3.1408 +        if (type_map[_graph.target(arc)] == LOWX ||
  3.1409 +            type_map[_graph.target(arc)] == HIGHX) {
  3.1410 +          break;
  3.1411 +        }
  3.1412 +        if (type_map[_graph.target(arc)] == 2) {
  3.1413 +          type_map[_graph.target(arc)] = 3;
  3.1414 +
  3.1415 +          arc = arc_lists[_graph.oppositeArc(arc)].next;
  3.1416 +          st.push_back(arc);
  3.1417 +        } else {
  3.1418 +          st.pop_back();
  3.1419 +          arc = arc_lists[arc].next;
  3.1420 +
  3.1421 +          while (_graph.oppositeArc(arc) == st.back()) {
  3.1422 +            arc = st.back();
  3.1423 +            st.pop_back();
  3.1424 +            arc = arc_lists[arc].next;
  3.1425 +          }
  3.1426 +          st.push_back(arc);
  3.1427 +        }
  3.1428 +      }
  3.1429 +
  3.1430 +      for (int i = 0; i < int(st.size()); ++i) {
  3.1431 +        if (type_map[_graph.target(st[i])] != LOWY &&
  3.1432 +            type_map[_graph.target(st[i])] != HIGHY) {
  3.1433 +          for (; i < int(st.size()); ++i) {
  3.1434 +            ipath.push_back(st[i]);
  3.1435 +          }
  3.1436 +        }
  3.1437 +      }
  3.1438 +    }
  3.1439 +
  3.1440 +    void setInternalFlags(std::vector<Arc>& ipath, TypeMap& type_map) {
  3.1441 +      for (int i = 1; i < int(ipath.size()); ++i) {
  3.1442 +        type_map[_graph.source(ipath[i])] = INTERNAL;
  3.1443 +      }
  3.1444 +    }
  3.1445 +
  3.1446 +    void findPilePath(std::vector<Arc>& ppath,
  3.1447 +                      Node root, TypeMap& type_map, OrderMap& order_map,
  3.1448 +                      NodeData& node_data, ArcLists& arc_lists) {
  3.1449 +      std::vector<Arc> st;
  3.1450 +
  3.1451 +      st.push_back(_graph.oppositeArc(node_data[order_map[root]].first));
  3.1452 +      st.push_back(node_data[order_map[root]].first);
  3.1453 +
  3.1454 +      while (st.size() > 1) {
  3.1455 +        Arc arc = st.back();
  3.1456 +        if (type_map[_graph.target(arc)] == INTERNAL) {
  3.1457 +          break;
  3.1458 +        }
  3.1459 +        if (type_map[_graph.target(arc)] == 3) {
  3.1460 +          type_map[_graph.target(arc)] = 4;
  3.1461 +
  3.1462 +          arc = arc_lists[_graph.oppositeArc(arc)].next;
  3.1463 +          st.push_back(arc);
  3.1464 +        } else {
  3.1465 +          st.pop_back();
  3.1466 +          arc = arc_lists[arc].next;
  3.1467 +
  3.1468 +          while (!st.empty() && _graph.oppositeArc(arc) == st.back()) {
  3.1469 +            arc = st.back();
  3.1470 +            st.pop_back();
  3.1471 +            arc = arc_lists[arc].next;
  3.1472 +          }
  3.1473 +          st.push_back(arc);
  3.1474 +        }
  3.1475 +      }
  3.1476 +
  3.1477 +      for (int i = 1; i < int(st.size()); ++i) {
  3.1478 +        ppath.push_back(st[i]);
  3.1479 +      }
  3.1480 +    }
  3.1481 +
  3.1482 +
  3.1483 +    int markExternalPath(Node node, OrderMap& order_map,
  3.1484 +                         ChildLists& child_lists, PredMap& pred_map,
  3.1485 +                         AncestorMap& ancestor_map, LowMap& low_map) {
  3.1486 +      int lp = lowPoint(node, order_map, child_lists,
  3.1487 +                        ancestor_map, low_map);
  3.1488 +
  3.1489 +      if (ancestor_map[node] != lp) {
  3.1490 +        node = child_lists[node].first;
  3.1491 +        _kuratowski[pred_map[node]] = true;
  3.1492 +
  3.1493 +        while (ancestor_map[node] != lp) {
  3.1494 +          for (OutArcIt e(_graph, node); e != INVALID; ++e) {
  3.1495 +            Node tnode = _graph.target(e);
  3.1496 +            if (order_map[tnode] > order_map[node] && low_map[tnode] == lp) {
  3.1497 +              node = tnode;
  3.1498 +              _kuratowski[e] = true;
  3.1499 +              break;
  3.1500 +            }
  3.1501 +          }
  3.1502 +        }
  3.1503 +      }
  3.1504 +
  3.1505 +      for (OutArcIt e(_graph, node); e != INVALID; ++e) {
  3.1506 +        if (order_map[_graph.target(e)] == lp) {
  3.1507 +          _kuratowski[e] = true;
  3.1508 +          break;
  3.1509 +        }
  3.1510 +      }
  3.1511 +
  3.1512 +      return lp;
  3.1513 +    }
  3.1514 +
  3.1515 +    void markPertinentPath(Node node, OrderMap& order_map,
  3.1516 +                           NodeData& node_data, ArcLists& arc_lists,
  3.1517 +                           EmbedArc& embed_arc, MergeRoots& merge_roots) {
  3.1518 +      while (embed_arc[node] == INVALID) {
  3.1519 +        int n = merge_roots[node].front();
  3.1520 +        Arc arc = node_data[n].first;
  3.1521 +
  3.1522 +        _kuratowski.set(arc, true);
  3.1523 +
  3.1524 +        Node pred = node;
  3.1525 +        node = _graph.target(arc);
  3.1526 +        while (!pertinent(node, embed_arc, merge_roots)) {
  3.1527 +          arc = node_data[order_map[node]].first;
  3.1528 +          if (_graph.target(arc) == pred) {
  3.1529 +            arc = arc_lists[arc].next;
  3.1530 +          }
  3.1531 +          _kuratowski.set(arc, true);
  3.1532 +          pred = node;
  3.1533 +          node = _graph.target(arc);
  3.1534 +        }
  3.1535 +      }
  3.1536 +      _kuratowski.set(embed_arc[node], true);
  3.1537 +    }
  3.1538 +
  3.1539 +    void markPredPath(Node node, Node snode, PredMap& pred_map) {
  3.1540 +      while (node != snode) {
  3.1541 +        _kuratowski.set(pred_map[node], true);
  3.1542 +        node = _graph.source(pred_map[node]);
  3.1543 +      }
  3.1544 +    }
  3.1545 +
  3.1546 +    void markFacePath(Node ynode, Node xnode,
  3.1547 +                      OrderMap& order_map, NodeData& node_data) {
  3.1548 +      Arc arc = node_data[order_map[ynode]].first;
  3.1549 +      Node node = _graph.target(arc);
  3.1550 +      _kuratowski.set(arc, true);
  3.1551 +
  3.1552 +      while (node != xnode) {
  3.1553 +        arc = node_data[order_map[node]].first;
  3.1554 +        _kuratowski.set(arc, true);
  3.1555 +        node = _graph.target(arc);
  3.1556 +      }
  3.1557 +    }
  3.1558 +
  3.1559 +    void markInternalPath(std::vector<Arc>& path) {
  3.1560 +      for (int i = 0; i < int(path.size()); ++i) {
  3.1561 +        _kuratowski.set(path[i], true);
  3.1562 +      }
  3.1563 +    }
  3.1564 +
  3.1565 +    void markPilePath(std::vector<Arc>& path) {
  3.1566 +      for (int i = 0; i < int(path.size()); ++i) {
  3.1567 +        _kuratowski.set(path[i], true);
  3.1568 +      }
  3.1569 +    }
  3.1570 +
  3.1571 +    void isolateKuratowski(Arc arc, NodeData& node_data,
  3.1572 +                           ArcLists& arc_lists, FlipMap& flip_map,
  3.1573 +                           OrderMap& order_map, OrderList& order_list,
  3.1574 +                           PredMap& pred_map, ChildLists& child_lists,
  3.1575 +                           AncestorMap& ancestor_map, LowMap& low_map,
  3.1576 +                           EmbedArc& embed_arc, MergeRoots& merge_roots) {
  3.1577 +
  3.1578 +      Node root = _graph.source(arc);
  3.1579 +      Node enode = _graph.target(arc);
  3.1580 +
  3.1581 +      int rorder = order_map[root];
  3.1582 +
  3.1583 +      TypeMap type_map(_graph, 0);
  3.1584 +
  3.1585 +      int rn = findComponentRoot(root, enode, child_lists,
  3.1586 +                                 order_map, order_list);
  3.1587 +
  3.1588 +      Node xnode = order_list[node_data[rn].next];
  3.1589 +      Node ynode = order_list[node_data[rn].prev];
  3.1590 +
  3.1591 +      // Minor-A
  3.1592 +      {
  3.1593 +        while (!merge_roots[xnode].empty() || !merge_roots[ynode].empty()) {
  3.1594 +
  3.1595 +          if (!merge_roots[xnode].empty()) {
  3.1596 +            root = xnode;
  3.1597 +            rn = merge_roots[xnode].front();
  3.1598 +          } else {
  3.1599 +            root = ynode;
  3.1600 +            rn = merge_roots[ynode].front();
  3.1601 +          }
  3.1602 +
  3.1603 +          xnode = order_list[node_data[rn].next];
  3.1604 +          ynode = order_list[node_data[rn].prev];
  3.1605 +        }
  3.1606 +
  3.1607 +        if (root != _graph.source(arc)) {
  3.1608 +          orientComponent(root, rn, order_map, pred_map,
  3.1609 +                          node_data, arc_lists, flip_map, type_map);
  3.1610 +          markFacePath(root, root, order_map, node_data);
  3.1611 +          int xlp = markExternalPath(xnode, order_map, child_lists,
  3.1612 +                                     pred_map, ancestor_map, low_map);
  3.1613 +          int ylp = markExternalPath(ynode, order_map, child_lists,
  3.1614 +                                     pred_map, ancestor_map, low_map);
  3.1615 +          markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
  3.1616 +          Node lwnode = findPertinent(ynode, order_map, node_data,
  3.1617 +                                      embed_arc, merge_roots);
  3.1618 +
  3.1619 +          markPertinentPath(lwnode, order_map, node_data, arc_lists,
  3.1620 +                            embed_arc, merge_roots);
  3.1621 +
  3.1622 +          return;
  3.1623 +        }
  3.1624 +      }
  3.1625 +
  3.1626 +      orientComponent(root, rn, order_map, pred_map,
  3.1627 +                      node_data, arc_lists, flip_map, type_map);
  3.1628 +
  3.1629 +      Node wnode = findPertinent(ynode, order_map, node_data,
  3.1630 +                                 embed_arc, merge_roots);
  3.1631 +      setFaceFlags(root, wnode, ynode, xnode, order_map, node_data, type_map);
  3.1632 +
  3.1633 +
  3.1634 +      //Minor-B
  3.1635 +      if (!merge_roots[wnode].empty()) {
  3.1636 +        int cn = merge_roots[wnode].back();
  3.1637 +        Node rep = order_list[cn - order_list.size()];
  3.1638 +        if (low_map[rep] < rorder) {
  3.1639 +          markFacePath(root, root, order_map, node_data);
  3.1640 +          int xlp = markExternalPath(xnode, order_map, child_lists,
  3.1641 +                                     pred_map, ancestor_map, low_map);
  3.1642 +          int ylp = markExternalPath(ynode, order_map, child_lists,
  3.1643 +                                     pred_map, ancestor_map, low_map);
  3.1644 +
  3.1645 +          Node lwnode, lznode;
  3.1646 +          markCommonPath(wnode, rorder, lwnode, lznode, order_list,
  3.1647 +                         order_map, node_data, arc_lists, embed_arc,
  3.1648 +                         merge_roots, child_lists, ancestor_map, low_map);
  3.1649 +
  3.1650 +          markPertinentPath(lwnode, order_map, node_data, arc_lists,
  3.1651 +                            embed_arc, merge_roots);
  3.1652 +          int zlp = markExternalPath(lznode, order_map, child_lists,
  3.1653 +                                     pred_map, ancestor_map, low_map);
  3.1654 +
  3.1655 +          int minlp = xlp < ylp ? xlp : ylp;
  3.1656 +          if (zlp < minlp) minlp = zlp;
  3.1657 +
  3.1658 +          int maxlp = xlp > ylp ? xlp : ylp;
  3.1659 +          if (zlp > maxlp) maxlp = zlp;
  3.1660 +
  3.1661 +          markPredPath(order_list[maxlp], order_list[minlp], pred_map);
  3.1662 +
  3.1663 +          return;
  3.1664 +        }
  3.1665 +      }
  3.1666 +
  3.1667 +      Node pxnode, pynode;
  3.1668 +      std::vector<Arc> ipath;
  3.1669 +      findInternalPath(ipath, wnode, root, type_map, order_map,
  3.1670 +                       node_data, arc_lists);
  3.1671 +      setInternalFlags(ipath, type_map);
  3.1672 +      pynode = _graph.source(ipath.front());
  3.1673 +      pxnode = _graph.target(ipath.back());
  3.1674 +
  3.1675 +      wnode = findPertinent(pynode, order_map, node_data,
  3.1676 +                            embed_arc, merge_roots);
  3.1677 +
  3.1678 +      // Minor-C
  3.1679 +      {
  3.1680 +        if (type_map[_graph.source(ipath.front())] == HIGHY) {
  3.1681 +          if (type_map[_graph.target(ipath.back())] == HIGHX) {
  3.1682 +            markFacePath(xnode, pxnode, order_map, node_data);
  3.1683 +          }
  3.1684 +          markFacePath(root, xnode, order_map, node_data);
  3.1685 +          markPertinentPath(wnode, order_map, node_data, arc_lists,
  3.1686 +                            embed_arc, merge_roots);
  3.1687 +          markInternalPath(ipath);
  3.1688 +          int xlp = markExternalPath(xnode, order_map, child_lists,
  3.1689 +                                     pred_map, ancestor_map, low_map);
  3.1690 +          int ylp = markExternalPath(ynode, order_map, child_lists,
  3.1691 +                                     pred_map, ancestor_map, low_map);
  3.1692 +          markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
  3.1693 +          return;
  3.1694 +        }
  3.1695 +
  3.1696 +        if (type_map[_graph.target(ipath.back())] == HIGHX) {
  3.1697 +          markFacePath(ynode, root, order_map, node_data);
  3.1698 +          markPertinentPath(wnode, order_map, node_data, arc_lists,
  3.1699 +                            embed_arc, merge_roots);
  3.1700 +          markInternalPath(ipath);
  3.1701 +          int xlp = markExternalPath(xnode, order_map, child_lists,
  3.1702 +                                     pred_map, ancestor_map, low_map);
  3.1703 +          int ylp = markExternalPath(ynode, order_map, child_lists,
  3.1704 +                                     pred_map, ancestor_map, low_map);
  3.1705 +          markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
  3.1706 +          return;
  3.1707 +        }
  3.1708 +      }
  3.1709 +
  3.1710 +      std::vector<Arc> ppath;
  3.1711 +      findPilePath(ppath, root, type_map, order_map, node_data, arc_lists);
  3.1712 +
  3.1713 +      // Minor-D
  3.1714 +      if (!ppath.empty()) {
  3.1715 +        markFacePath(ynode, xnode, order_map, node_data);
  3.1716 +        markPertinentPath(wnode, order_map, node_data, arc_lists,
  3.1717 +                          embed_arc, merge_roots);
  3.1718 +        markPilePath(ppath);
  3.1719 +        markInternalPath(ipath);
  3.1720 +        int xlp = markExternalPath(xnode, order_map, child_lists,
  3.1721 +                                   pred_map, ancestor_map, low_map);
  3.1722 +        int ylp = markExternalPath(ynode, order_map, child_lists,
  3.1723 +                                   pred_map, ancestor_map, low_map);
  3.1724 +        markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
  3.1725 +        return;
  3.1726 +      }
  3.1727 +
  3.1728 +      // Minor-E*
  3.1729 +      {
  3.1730 +
  3.1731 +        if (!external(wnode, rorder, child_lists, ancestor_map, low_map)) {
  3.1732 +          Node znode = findExternal(pynode, rorder, order_map,
  3.1733 +                                    child_lists, ancestor_map,
  3.1734 +                                    low_map, node_data);
  3.1735 +
  3.1736 +          if (type_map[znode] == LOWY) {
  3.1737 +            markFacePath(root, xnode, order_map, node_data);
  3.1738 +            markPertinentPath(wnode, order_map, node_data, arc_lists,
  3.1739 +                              embed_arc, merge_roots);
  3.1740 +            markInternalPath(ipath);
  3.1741 +            int xlp = markExternalPath(xnode, order_map, child_lists,
  3.1742 +                                       pred_map, ancestor_map, low_map);
  3.1743 +            int zlp = markExternalPath(znode, order_map, child_lists,
  3.1744 +                                       pred_map, ancestor_map, low_map);
  3.1745 +            markPredPath(root, order_list[xlp < zlp ? xlp : zlp], pred_map);
  3.1746 +          } else {
  3.1747 +            markFacePath(ynode, root, order_map, node_data);
  3.1748 +            markPertinentPath(wnode, order_map, node_data, arc_lists,
  3.1749 +                              embed_arc, merge_roots);
  3.1750 +            markInternalPath(ipath);
  3.1751 +            int ylp = markExternalPath(ynode, order_map, child_lists,
  3.1752 +                                       pred_map, ancestor_map, low_map);
  3.1753 +            int zlp = markExternalPath(znode, order_map, child_lists,
  3.1754 +                                       pred_map, ancestor_map, low_map);
  3.1755 +            markPredPath(root, order_list[ylp < zlp ? ylp : zlp], pred_map);
  3.1756 +          }
  3.1757 +          return;
  3.1758 +        }
  3.1759 +
  3.1760 +        int xlp = markExternalPath(xnode, order_map, child_lists,
  3.1761 +                                   pred_map, ancestor_map, low_map);
  3.1762 +        int ylp = markExternalPath(ynode, order_map, child_lists,
  3.1763 +                                   pred_map, ancestor_map, low_map);
  3.1764 +        int wlp = markExternalPath(wnode, order_map, child_lists,
  3.1765 +                                   pred_map, ancestor_map, low_map);
  3.1766 +
  3.1767 +        if (wlp > xlp && wlp > ylp) {
  3.1768 +          markFacePath(root, root, order_map, node_data);
  3.1769 +          markPredPath(root, order_list[xlp < ylp ? xlp : ylp], pred_map);
  3.1770 +          return;
  3.1771 +        }
  3.1772 +
  3.1773 +        markInternalPath(ipath);
  3.1774 +        markPertinentPath(wnode, order_map, node_data, arc_lists,
  3.1775 +                          embed_arc, merge_roots);
  3.1776 +
  3.1777 +        if (xlp > ylp && xlp > wlp) {
  3.1778 +          markFacePath(root, pynode, order_map, node_data);
  3.1779 +          markFacePath(wnode, xnode, order_map, node_data);
  3.1780 +          markPredPath(root, order_list[ylp < wlp ? ylp : wlp], pred_map);
  3.1781 +          return;
  3.1782 +        }
  3.1783 +
  3.1784 +        if (ylp > xlp && ylp > wlp) {
  3.1785 +          markFacePath(pxnode, root, order_map, node_data);
  3.1786 +          markFacePath(ynode, wnode, order_map, node_data);
  3.1787 +          markPredPath(root, order_list[xlp < wlp ? xlp : wlp], pred_map);
  3.1788 +          return;
  3.1789 +        }
  3.1790 +
  3.1791 +        if (pynode != ynode) {
  3.1792 +          markFacePath(pxnode, wnode, order_map, node_data);
  3.1793 +
  3.1794 +          int minlp = xlp < ylp ? xlp : ylp;
  3.1795 +          if (wlp < minlp) minlp = wlp;
  3.1796 +
  3.1797 +          int maxlp = xlp > ylp ? xlp : ylp;
  3.1798 +          if (wlp > maxlp) maxlp = wlp;
  3.1799 +
  3.1800 +          markPredPath(order_list[maxlp], order_list[minlp], pred_map);
  3.1801 +          return;
  3.1802 +        }
  3.1803 +
  3.1804 +        if (pxnode != xnode) {
  3.1805 +          markFacePath(wnode, pynode, order_map, node_data);
  3.1806 +
  3.1807 +          int minlp = xlp < ylp ? xlp : ylp;
  3.1808 +          if (wlp < minlp) minlp = wlp;
  3.1809 +
  3.1810 +          int maxlp = xlp > ylp ? xlp : ylp;
  3.1811 +          if (wlp > maxlp) maxlp = wlp;
  3.1812 +
  3.1813 +          markPredPath(order_list[maxlp], order_list[minlp], pred_map);
  3.1814 +          return;
  3.1815 +        }
  3.1816 +
  3.1817 +        markFacePath(root, root, order_map, node_data);
  3.1818 +        int minlp = xlp < ylp ? xlp : ylp;
  3.1819 +        if (wlp < minlp) minlp = wlp;
  3.1820 +        markPredPath(root, order_list[minlp], pred_map);
  3.1821 +        return;
  3.1822 +      }
  3.1823 +
  3.1824 +    }
  3.1825 +
  3.1826 +  };
  3.1827 +
  3.1828 +  namespace _planarity_bits {
  3.1829 +
  3.1830 +    template <typename Graph, typename EmbeddingMap>
  3.1831 +    void makeConnected(Graph& graph, EmbeddingMap& embedding) {
  3.1832 +      DfsVisitor<Graph> null_visitor;
  3.1833 +      DfsVisit<Graph, DfsVisitor<Graph> > dfs(graph, null_visitor);
  3.1834 +      dfs.init();
  3.1835 +
  3.1836 +      typename Graph::Node u = INVALID;
  3.1837 +      for (typename Graph::NodeIt n(graph); n != INVALID; ++n) {
  3.1838 +        if (!dfs.reached(n)) {
  3.1839 +          dfs.addSource(n);
  3.1840 +          dfs.start();
  3.1841 +          if (u == INVALID) {
  3.1842 +            u = n;
  3.1843 +          } else {
  3.1844 +            typename Graph::Node v = n;
  3.1845 +
  3.1846 +            typename Graph::Arc ue = typename Graph::OutArcIt(graph, u);
  3.1847 +            typename Graph::Arc ve = typename Graph::OutArcIt(graph, v);
  3.1848 +
  3.1849 +            typename Graph::Arc e = graph.direct(graph.addEdge(u, v), true);
  3.1850 +
  3.1851 +            if (ue != INVALID) {
  3.1852 +              embedding[e] = embedding[ue];
  3.1853 +              embedding[ue] = e;
  3.1854 +            } else {
  3.1855 +              embedding[e] = e;
  3.1856 +            }
  3.1857 +
  3.1858 +            if (ve != INVALID) {
  3.1859 +              embedding[graph.oppositeArc(e)] = embedding[ve];
  3.1860 +              embedding[ve] = graph.oppositeArc(e);
  3.1861 +            } else {
  3.1862 +              embedding[graph.oppositeArc(e)] = graph.oppositeArc(e);
  3.1863 +            }
  3.1864 +          }
  3.1865 +        }
  3.1866 +      }
  3.1867 +    }
  3.1868 +
  3.1869 +    template <typename Graph, typename EmbeddingMap>
  3.1870 +    void makeBiNodeConnected(Graph& graph, EmbeddingMap& embedding) {
  3.1871 +      typename Graph::template ArcMap<bool> processed(graph);
  3.1872 +
  3.1873 +      std::vector<typename Graph::Arc> arcs;
  3.1874 +      for (typename Graph::ArcIt e(graph); e != INVALID; ++e) {
  3.1875 +        arcs.push_back(e);
  3.1876 +      }
  3.1877 +
  3.1878 +      IterableBoolMap<Graph, typename Graph::Node> visited(graph, false);
  3.1879 +
  3.1880 +      for (int i = 0; i < int(arcs.size()); ++i) {
  3.1881 +        typename Graph::Arc pp = arcs[i];
  3.1882 +        if (processed[pp]) continue;
  3.1883 +
  3.1884 +        typename Graph::Arc e = embedding[graph.oppositeArc(pp)];
  3.1885 +        processed[e] = true;
  3.1886 +        visited.set(graph.source(e), true);
  3.1887 +
  3.1888 +        typename Graph::Arc p = e, l = e;
  3.1889 +        e = embedding[graph.oppositeArc(e)];
  3.1890 +
  3.1891 +        while (e != l) {
  3.1892 +          processed[e] = true;
  3.1893 +
  3.1894 +          if (visited[graph.source(e)]) {
  3.1895 +
  3.1896 +            typename Graph::Arc n =
  3.1897 +              graph.direct(graph.addEdge(graph.source(p),
  3.1898 +                                           graph.target(e)), true);
  3.1899 +            embedding[n] = p;
  3.1900 +            embedding[graph.oppositeArc(pp)] = n;
  3.1901 +
  3.1902 +            embedding[graph.oppositeArc(n)] =
  3.1903 +              embedding[graph.oppositeArc(e)];
  3.1904 +            embedding[graph.oppositeArc(e)] =
  3.1905 +              graph.oppositeArc(n);
  3.1906 +
  3.1907 +            p = n;
  3.1908 +            e = embedding[graph.oppositeArc(n)];
  3.1909 +          } else {
  3.1910 +            visited.set(graph.source(e), true);
  3.1911 +            pp = p;
  3.1912 +            p = e;
  3.1913 +            e = embedding[graph.oppositeArc(e)];
  3.1914 +          }
  3.1915 +        }
  3.1916 +        visited.setAll(false);
  3.1917 +      }
  3.1918 +    }
  3.1919 +
  3.1920 +
  3.1921 +    template <typename Graph, typename EmbeddingMap>
  3.1922 +    void makeMaxPlanar(Graph& graph, EmbeddingMap& embedding) {
  3.1923 +
  3.1924 +      typename Graph::template NodeMap<int> degree(graph);
  3.1925 +
  3.1926 +      for (typename Graph::NodeIt n(graph); n != INVALID; ++n) {
  3.1927 +        degree[n] = countIncEdges(graph, n);
  3.1928 +      }
  3.1929 +
  3.1930 +      typename Graph::template ArcMap<bool> processed(graph);
  3.1931 +      IterableBoolMap<Graph, typename Graph::Node> visited(graph, false);
  3.1932 +
  3.1933 +      std::vector<typename Graph::Arc> arcs;
  3.1934 +      for (typename Graph::ArcIt e(graph); e != INVALID; ++e) {
  3.1935 +        arcs.push_back(e);
  3.1936 +      }
  3.1937 +
  3.1938 +      for (int i = 0; i < int(arcs.size()); ++i) {
  3.1939 +        typename Graph::Arc e = arcs[i];
  3.1940 +
  3.1941 +        if (processed[e]) continue;
  3.1942 +        processed[e] = true;
  3.1943 +
  3.1944 +        typename Graph::Arc mine = e;
  3.1945 +        int mind = degree[graph.source(e)];
  3.1946 +
  3.1947 +        int face_size = 1;
  3.1948 +
  3.1949 +        typename Graph::Arc l = e;
  3.1950 +        e = embedding[graph.oppositeArc(e)];
  3.1951 +        while (l != e) {
  3.1952 +          processed[e] = true;
  3.1953 +
  3.1954 +          ++face_size;
  3.1955 +
  3.1956 +          if (degree[graph.source(e)] < mind) {
  3.1957 +            mine = e;
  3.1958 +            mind = degree[graph.source(e)];
  3.1959 +          }
  3.1960 +
  3.1961 +          e = embedding[graph.oppositeArc(e)];
  3.1962 +        }
  3.1963 +
  3.1964 +        if (face_size < 4) {
  3.1965 +          continue;
  3.1966 +        }
  3.1967 +
  3.1968 +        typename Graph::Node s = graph.source(mine);
  3.1969 +        for (typename Graph::OutArcIt e(graph, s); e != INVALID; ++e) {
  3.1970 +          visited.set(graph.target(e), true);
  3.1971 +        }
  3.1972 +
  3.1973 +        typename Graph::Arc oppe = INVALID;
  3.1974 +
  3.1975 +        e = embedding[graph.oppositeArc(mine)];
  3.1976 +        e = embedding[graph.oppositeArc(e)];
  3.1977 +        while (graph.target(e) != s) {
  3.1978 +          if (visited[graph.source(e)]) {
  3.1979 +            oppe = e;
  3.1980 +            break;
  3.1981 +          }
  3.1982 +          e = embedding[graph.oppositeArc(e)];
  3.1983 +        }
  3.1984 +        visited.setAll(false);
  3.1985 +
  3.1986 +        if (oppe == INVALID) {
  3.1987 +
  3.1988 +          e = embedding[graph.oppositeArc(mine)];
  3.1989 +          typename Graph::Arc pn = mine, p = e;
  3.1990 +
  3.1991 +          e = embedding[graph.oppositeArc(e)];
  3.1992 +          while (graph.target(e) != s) {
  3.1993 +            typename Graph::Arc n =
  3.1994 +              graph.direct(graph.addEdge(s, graph.source(e)), true);
  3.1995 +
  3.1996 +            embedding[n] = pn;
  3.1997 +            embedding[graph.oppositeArc(n)] = e;
  3.1998 +            embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);
  3.1999 +
  3.2000 +            pn = n;
  3.2001 +
  3.2002 +            p = e;
  3.2003 +            e = embedding[graph.oppositeArc(e)];
  3.2004 +          }
  3.2005 +
  3.2006 +          embedding[graph.oppositeArc(e)] = pn;
  3.2007 +
  3.2008 +        } else {
  3.2009 +
  3.2010 +          mine = embedding[graph.oppositeArc(mine)];
  3.2011 +          s = graph.source(mine);
  3.2012 +          oppe = embedding[graph.oppositeArc(oppe)];
  3.2013 +          typename Graph::Node t = graph.source(oppe);
  3.2014 +
  3.2015 +          typename Graph::Arc ce = graph.direct(graph.addEdge(s, t), true);
  3.2016 +          embedding[ce] = mine;
  3.2017 +          embedding[graph.oppositeArc(ce)] = oppe;
  3.2018 +
  3.2019 +          typename Graph::Arc pn = ce, p = oppe;
  3.2020 +          e = embedding[graph.oppositeArc(oppe)];
  3.2021 +          while (graph.target(e) != s) {
  3.2022 +            typename Graph::Arc n =
  3.2023 +              graph.direct(graph.addEdge(s, graph.source(e)), true);
  3.2024 +
  3.2025 +            embedding[n] = pn;
  3.2026 +            embedding[graph.oppositeArc(n)] = e;
  3.2027 +            embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);
  3.2028 +
  3.2029 +            pn = n;
  3.2030 +
  3.2031 +            p = e;
  3.2032 +            e = embedding[graph.oppositeArc(e)];
  3.2033 +
  3.2034 +          }
  3.2035 +          embedding[graph.oppositeArc(e)] = pn;
  3.2036 +
  3.2037 +          pn = graph.oppositeArc(ce), p = mine;
  3.2038 +          e = embedding[graph.oppositeArc(mine)];
  3.2039 +          while (graph.target(e) != t) {
  3.2040 +            typename Graph::Arc n =
  3.2041 +              graph.direct(graph.addEdge(t, graph.source(e)), true);
  3.2042 +
  3.2043 +            embedding[n] = pn;
  3.2044 +            embedding[graph.oppositeArc(n)] = e;
  3.2045 +            embedding[graph.oppositeArc(p)] = graph.oppositeArc(n);
  3.2046 +
  3.2047 +            pn = n;
  3.2048 +
  3.2049 +            p = e;
  3.2050 +            e = embedding[graph.oppositeArc(e)];
  3.2051 +
  3.2052 +          }
  3.2053 +          embedding[graph.oppositeArc(e)] = pn;
  3.2054 +        }
  3.2055 +      }
  3.2056 +    }
  3.2057 +
  3.2058 +  }
  3.2059 +
  3.2060 +  /// \ingroup planar
  3.2061 +  ///
  3.2062 +  /// \brief Schnyder's planar drawing algorithm
  3.2063 +  ///
  3.2064 +  /// The planar drawing algorithm calculates positions for the nodes
  3.2065 +  /// in the plane which coordinates satisfy that if the arcs are
  3.2066 +  /// represented with straight lines then they will not intersect
  3.2067 +  /// each other.
  3.2068 +  ///
  3.2069 +  /// Scnyder's algorithm embeds the graph on \c (n-2,n-2) size grid,
  3.2070 +  /// i.e. each node will be located in the \c [0,n-2]x[0,n-2] square.
  3.2071 +  /// The time complexity of the algorithm is O(n).
  3.2072 +  template <typename Graph>
  3.2073 +  class PlanarDrawing {
  3.2074 +  public:
  3.2075 +
  3.2076 +    TEMPLATE_GRAPH_TYPEDEFS(Graph);
  3.2077 +
  3.2078 +    /// \brief The point type for store coordinates
  3.2079 +    typedef dim2::Point<int> Point;
  3.2080 +    /// \brief The map type for store coordinates
  3.2081 +    typedef typename Graph::template NodeMap<Point> PointMap;
  3.2082 +
  3.2083 +
  3.2084 +    /// \brief Constructor
  3.2085 +    ///
  3.2086 +    /// Constructor
  3.2087 +    /// \pre The graph should be simple, i.e. loop and parallel arc free.
  3.2088 +    PlanarDrawing(const Graph& graph)
  3.2089 +      : _graph(graph), _point_map(graph) {}
  3.2090 +
  3.2091 +  private:
  3.2092 +
  3.2093 +    template <typename AuxGraph, typename AuxEmbeddingMap>
  3.2094 +    void drawing(const AuxGraph& graph,
  3.2095 +                 const AuxEmbeddingMap& next,
  3.2096 +                 PointMap& point_map) {
  3.2097 +      TEMPLATE_GRAPH_TYPEDEFS(AuxGraph);
  3.2098 +
  3.2099 +      typename AuxGraph::template ArcMap<Arc> prev(graph);
  3.2100 +
  3.2101 +      for (NodeIt n(graph); n != INVALID; ++n) {
  3.2102 +        Arc e = OutArcIt(graph, n);
  3.2103 +
  3.2104 +        Arc p = e, l = e;
  3.2105 +
  3.2106 +        e = next[e];
  3.2107 +        while (e != l) {
  3.2108 +          prev[e] = p;
  3.2109 +          p = e;
  3.2110 +          e = next[e];
  3.2111 +        }
  3.2112 +        prev[e] = p;
  3.2113 +      }
  3.2114 +
  3.2115 +      Node anode, bnode, cnode;
  3.2116 +
  3.2117 +      {
  3.2118 +        Arc e = ArcIt(graph);
  3.2119 +        anode = graph.source(e);
  3.2120 +        bnode = graph.target(e);
  3.2121 +        cnode = graph.target(next[graph.oppositeArc(e)]);
  3.2122 +      }
  3.2123 +
  3.2124 +      IterableBoolMap<AuxGraph, Node> proper(graph, false);
  3.2125 +      typename AuxGraph::template NodeMap<int> conn(graph, -1);
  3.2126 +
  3.2127 +      conn[anode] = conn[bnode] = -2;
  3.2128 +      {
  3.2129 +        for (OutArcIt e(graph, anode); e != INVALID; ++e) {
  3.2130 +          Node m = graph.target(e);
  3.2131 +          if (conn[m] == -1) {
  3.2132 +            conn[m] = 1;
  3.2133 +          }
  3.2134 +        }
  3.2135 +        conn[cnode] = 2;
  3.2136 +
  3.2137 +        for (OutArcIt e(graph, bnode); e != INVALID; ++e) {
  3.2138 +          Node m = graph.target(e);
  3.2139 +          if (conn[m] == -1) {
  3.2140 +            conn[m] = 1;
  3.2141 +          } else if (conn[m] != -2) {
  3.2142 +            conn[m] += 1;
  3.2143 +            Arc pe = graph.oppositeArc(e);
  3.2144 +            if (conn[graph.target(next[pe])] == -2) {
  3.2145 +              conn[m] -= 1;
  3.2146 +            }
  3.2147 +            if (conn[graph.target(prev[pe])] == -2) {
  3.2148 +              conn[m] -= 1;
  3.2149 +            }
  3.2150 +
  3.2151 +            proper.set(m, conn[m] == 1);
  3.2152 +          }
  3.2153 +        }
  3.2154 +      }
  3.2155 +
  3.2156 +
  3.2157 +      typename AuxGraph::template ArcMap<int> angle(graph, -1);
  3.2158 +
  3.2159 +      while (proper.trueNum() != 0) {
  3.2160 +        Node n = typename IterableBoolMap<AuxGraph, Node>::TrueIt(proper);
  3.2161 +        proper.set(n, false);
  3.2162 +        conn[n] = -2;
  3.2163 +
  3.2164 +        for (OutArcIt e(graph, n); e != INVALID; ++e) {
  3.2165 +          Node m = graph.target(e);
  3.2166 +          if (conn[m] == -1) {
  3.2167 +            conn[m] = 1;
  3.2168 +          } else if (conn[m] != -2) {
  3.2169 +            conn[m] += 1;
  3.2170 +            Arc pe = graph.oppositeArc(e);
  3.2171 +            if (conn[graph.target(next[pe])] == -2) {
  3.2172 +              conn[m] -= 1;
  3.2173 +            }
  3.2174 +            if (conn[graph.target(prev[pe])] == -2) {
  3.2175 +              conn[m] -= 1;
  3.2176 +            }
  3.2177 +
  3.2178 +            proper.set(m, conn[m] == 1);
  3.2179 +          }
  3.2180 +        }
  3.2181 +
  3.2182 +        {
  3.2183 +          Arc e = OutArcIt(graph, n);
  3.2184 +          Arc p = e, l = e;
  3.2185 +
  3.2186 +          e = next[e];
  3.2187 +          while (e != l) {
  3.2188 +
  3.2189 +            if (conn[graph.target(e)] == -2 && conn[graph.target(p)] == -2) {
  3.2190 +              Arc f = e;
  3.2191 +              angle[f] = 0;
  3.2192 +              f = next[graph.oppositeArc(f)];
  3.2193 +              angle[f] = 1;
  3.2194 +              f = next[graph.oppositeArc(f)];
  3.2195 +              angle[f] = 2;
  3.2196 +            }
  3.2197 +
  3.2198 +            p = e;
  3.2199 +            e = next[e];
  3.2200 +          }
  3.2201 +
  3.2202 +          if (conn[graph.target(e)] == -2 && conn[graph.target(p)] == -2) {
  3.2203 +            Arc f = e;
  3.2204 +            angle[f] = 0;
  3.2205 +            f = next[graph.oppositeArc(f)];
  3.2206 +            angle[f] = 1;
  3.2207 +            f = next[graph.oppositeArc(f)];
  3.2208 +            angle[f] = 2;
  3.2209 +          }
  3.2210 +        }
  3.2211 +      }
  3.2212 +
  3.2213 +      typename AuxGraph::template NodeMap<Node> apred(graph, INVALID);
  3.2214 +      typename AuxGraph::template NodeMap<Node> bpred(graph, INVALID);
  3.2215 +      typename AuxGraph::template NodeMap<Node> cpred(graph, INVALID);
  3.2216 +
  3.2217 +      typename AuxGraph::template NodeMap<int> apredid(graph, -1);
  3.2218 +      typename AuxGraph::template NodeMap<int> bpredid(graph, -1);
  3.2219 +      typename AuxGraph::template NodeMap<int> cpredid(graph, -1);
  3.2220 +
  3.2221 +      for (ArcIt e(graph); e != INVALID; ++e) {
  3.2222 +        if (angle[e] == angle[next[e]]) {
  3.2223 +          switch (angle[e]) {
  3.2224 +          case 2:
  3.2225 +            apred[graph.target(e)] = graph.source(e);
  3.2226 +            apredid[graph.target(e)] = graph.id(graph.source(e));
  3.2227 +            break;
  3.2228 +          case 1:
  3.2229 +            bpred[graph.target(e)] = graph.source(e);
  3.2230 +            bpredid[graph.target(e)] = graph.id(graph.source(e));
  3.2231 +            break;
  3.2232 +          case 0:
  3.2233 +            cpred[graph.target(e)] = graph.source(e);
  3.2234 +            cpredid[graph.target(e)] = graph.id(graph.source(e));
  3.2235 +            break;
  3.2236 +          }
  3.2237 +        }
  3.2238 +      }
  3.2239 +
  3.2240 +      cpred[anode] = INVALID;
  3.2241 +      cpred[bnode] = INVALID;
  3.2242 +
  3.2243 +      std::vector<Node> aorder, border, corder;
  3.2244 +
  3.2245 +      {
  3.2246 +        typename AuxGraph::template NodeMap<bool> processed(graph, false);
  3.2247 +        std::vector<Node> st;
  3.2248 +        for (NodeIt n(graph); n != INVALID; ++n) {
  3.2249 +          if (!processed[n] && n != bnode && n != cnode) {
  3.2250 +            st.push_back(n);
  3.2251 +            processed[n] = true;
  3.2252 +            Node m = apred[n];
  3.2253 +            while (m != INVALID && !processed[m]) {
  3.2254 +              st.push_back(m);
  3.2255 +              processed[m] = true;
  3.2256 +              m = apred[m];
  3.2257 +            }
  3.2258 +            while (!st.empty()) {
  3.2259 +              aorder.push_back(st.back());
  3.2260 +              st.pop_back();
  3.2261 +            }
  3.2262 +          }
  3.2263 +        }
  3.2264 +      }
  3.2265 +
  3.2266 +      {
  3.2267 +        typename AuxGraph::template NodeMap<bool> processed(graph, false);
  3.2268 +        std::vector<Node> st;
  3.2269 +        for (NodeIt n(graph); n != INVALID; ++n) {
  3.2270 +          if (!processed[n] && n != cnode && n != anode) {
  3.2271 +            st.push_back(n);
  3.2272 +            processed[n] = true;
  3.2273 +            Node m = bpred[n];
  3.2274 +            while (m != INVALID && !processed[m]) {
  3.2275 +              st.push_back(m);
  3.2276 +              processed[m] = true;
  3.2277 +              m = bpred[m];
  3.2278 +            }
  3.2279 +            while (!st.empty()) {
  3.2280 +              border.push_back(st.back());
  3.2281 +              st.pop_back();
  3.2282 +            }
  3.2283 +          }
  3.2284 +        }
  3.2285 +      }
  3.2286 +
  3.2287 +      {
  3.2288 +        typename AuxGraph::template NodeMap<bool> processed(graph, false);
  3.2289 +        std::vector<Node> st;
  3.2290 +        for (NodeIt n(graph); n != INVALID; ++n) {
  3.2291 +          if (!processed[n] && n != anode && n != bnode) {
  3.2292 +            st.push_back(n);
  3.2293 +            processed[n] = true;
  3.2294 +            Node m = cpred[n];
  3.2295 +            while (m != INVALID && !processed[m]) {
  3.2296 +              st.push_back(m);
  3.2297 +              processed[m] = true;
  3.2298 +              m = cpred[m];
  3.2299 +            }
  3.2300 +            while (!st.empty()) {
  3.2301 +              corder.push_back(st.back());
  3.2302 +              st.pop_back();
  3.2303 +            }
  3.2304 +          }
  3.2305 +        }
  3.2306 +      }
  3.2307 +
  3.2308 +      typename AuxGraph::template NodeMap<int> atree(graph, 0);
  3.2309 +      for (int i = aorder.size() - 1; i >= 0; --i) {
  3.2310 +        Node n = aorder[i];
  3.2311 +        atree[n] = 1;
  3.2312 +        for (OutArcIt e(graph, n); e != INVALID; ++e) {
  3.2313 +          if (apred[graph.target(e)] == n) {
  3.2314 +            atree[n] += atree[graph.target(e)];
  3.2315 +          }
  3.2316 +        }
  3.2317 +      }
  3.2318 +
  3.2319 +      typename AuxGraph::template NodeMap<int> btree(graph, 0);
  3.2320 +      for (int i = border.size() - 1; i >= 0; --i) {
  3.2321 +        Node n = border[i];
  3.2322 +        btree[n] = 1;
  3.2323 +        for (OutArcIt e(graph, n); e != INVALID; ++e) {
  3.2324 +          if (bpred[graph.target(e)] == n) {
  3.2325 +            btree[n] += btree[graph.target(e)];
  3.2326 +          }
  3.2327 +        }
  3.2328 +      }
  3.2329 +
  3.2330 +      typename AuxGraph::template NodeMap<int> apath(graph, 0);
  3.2331 +      apath[bnode] = apath[cnode] = 1;
  3.2332 +      typename AuxGraph::template NodeMap<int> apath_btree(graph, 0);
  3.2333 +      apath_btree[bnode] = btree[bnode];
  3.2334 +      for (int i = 1; i < int(aorder.size()); ++i) {
  3.2335 +        Node n = aorder[i];
  3.2336 +        apath[n] = apath[apred[n]] + 1;
  3.2337 +        apath_btree[n] = btree[n] + apath_btree[apred[n]];
  3.2338 +      }
  3.2339 +
  3.2340 +      typename AuxGraph::template NodeMap<int> bpath_atree(graph, 0);
  3.2341 +      bpath_atree[anode] = atree[anode];
  3.2342 +      for (int i = 1; i < int(border.size()); ++i) {
  3.2343 +        Node n = border[i];
  3.2344 +        bpath_atree[n] = atree[n] + bpath_atree[bpred[n]];
  3.2345 +      }
  3.2346 +
  3.2347 +      typename AuxGraph::template NodeMap<int> cpath(graph, 0);
  3.2348 +      cpath[anode] = cpath[bnode] = 1;
  3.2349 +      typename AuxGraph::template NodeMap<int> cpath_atree(graph, 0);
  3.2350 +      cpath_atree[anode] = atree[anode];
  3.2351 +      typename AuxGraph::template NodeMap<int> cpath_btree(graph, 0);
  3.2352 +      cpath_btree[bnode] = btree[bnode];
  3.2353 +      for (int i = 1; i < int(corder.size()); ++i) {
  3.2354 +        Node n = corder[i];
  3.2355 +        cpath[n] = cpath[cpred[n]] + 1;
  3.2356 +        cpath_atree[n] = atree[n] + cpath_atree[cpred[n]];
  3.2357 +        cpath_btree[n] = btree[n] + cpath_btree[cpred[n]];
  3.2358 +      }
  3.2359 +
  3.2360 +      typename AuxGraph::template NodeMap<int> third(graph);
  3.2361 +      for (NodeIt n(graph); n != INVALID; ++n) {
  3.2362 +        point_map[n].x =
  3.2363 +          bpath_atree[n] + cpath_atree[n] - atree[n] - cpath[n] + 1;
  3.2364 +        point_map[n].y =
  3.2365 +          cpath_btree[n] + apath_btree[n] - btree[n] - apath[n] + 1;
  3.2366 +      }
  3.2367 +
  3.2368 +    }
  3.2369 +
  3.2370 +  public:
  3.2371 +
  3.2372 +    /// \brief Calculates the node positions
  3.2373 +    ///
  3.2374 +    /// This function calculates the node positions.
  3.2375 +    /// \return %True if the graph is planar.
  3.2376 +    bool run() {
  3.2377 +      PlanarEmbedding<Graph> pe(_graph);
  3.2378 +      if (!pe.run()) return false;
  3.2379 +
  3.2380 +      run(pe);
  3.2381 +      return true;
  3.2382 +    }
  3.2383 +
  3.2384 +    /// \brief Calculates the node positions according to a
  3.2385 +    /// combinatorical embedding
  3.2386 +    ///
  3.2387 +    /// This function calculates the node locations. The \c embedding
  3.2388 +    /// parameter should contain a valid combinatorical embedding, i.e.
  3.2389 +    /// a valid cyclic order of the arcs.
  3.2390 +    template <typename EmbeddingMap>
  3.2391 +    void run(const EmbeddingMap& embedding) {
  3.2392 +      typedef SmartEdgeSet<Graph> AuxGraph;
  3.2393 +
  3.2394 +      if (3 * countNodes(_graph) - 6 == countEdges(_graph)) {
  3.2395 +        drawing(_graph, embedding, _point_map);
  3.2396 +        return;
  3.2397 +      }
  3.2398 +
  3.2399 +      AuxGraph aux_graph(_graph);
  3.2400 +      typename AuxGraph::template ArcMap<typename AuxGraph::Arc>
  3.2401 +        aux_embedding(aux_graph);
  3.2402 +
  3.2403 +      {
  3.2404 +
  3.2405 +        typename Graph::template EdgeMap<typename AuxGraph::Edge>
  3.2406 +          ref(_graph);
  3.2407 +
  3.2408 +        for (EdgeIt e(_graph); e != INVALID; ++e) {
  3.2409 +          ref[e] = aux_graph.addEdge(_graph.u(e), _graph.v(e));
  3.2410 +        }
  3.2411 +
  3.2412 +        for (EdgeIt e(_graph); e != INVALID; ++e) {
  3.2413 +          Arc ee = embedding[_graph.direct(e, true)];
  3.2414 +          aux_embedding[aux_graph.direct(ref[e], true)] =
  3.2415 +            aux_graph.direct(ref[ee], _graph.direction(ee));
  3.2416 +          ee = embedding[_graph.direct(e, false)];
  3.2417 +          aux_embedding[aux_graph.direct(ref[e], false)] =
  3.2418 +            aux_graph.direct(ref[ee], _graph.direction(ee));
  3.2419 +        }
  3.2420 +      }
  3.2421 +      _planarity_bits::makeConnected(aux_graph, aux_embedding);
  3.2422 +      _planarity_bits::makeBiNodeConnected(aux_graph, aux_embedding);
  3.2423 +      _planarity_bits::makeMaxPlanar(aux_graph, aux_embedding);
  3.2424 +      drawing(aux_graph, aux_embedding, _point_map);
  3.2425 +    }
  3.2426 +
  3.2427 +    /// \brief The coordinate of the given node
  3.2428 +    ///
  3.2429 +    /// The coordinate of the given node.
  3.2430 +    Point operator[](const Node& node) const {
  3.2431 +      return _point_map[node];
  3.2432 +    }
  3.2433 +
  3.2434 +    /// \brief Returns the grid embedding in a \e NodeMap.
  3.2435 +    ///
  3.2436 +    /// Returns the grid embedding in a \e NodeMap of \c dim2::Point<int> .
  3.2437 +    const PointMap& coords() const {
  3.2438 +      return _point_map;
  3.2439 +    }
  3.2440 +
  3.2441 +  private:
  3.2442 +
  3.2443 +    const Graph& _graph;
  3.2444 +    PointMap _point_map;
  3.2445 +
  3.2446 +  };
  3.2447 +
  3.2448 +  namespace _planarity_bits {
  3.2449 +
  3.2450 +    template <typename ColorMap>
  3.2451 +    class KempeFilter {
  3.2452 +    public:
  3.2453 +      typedef typename ColorMap::Key Key;
  3.2454 +      typedef bool Value;
  3.2455 +
  3.2456 +      KempeFilter(const ColorMap& color_map,
  3.2457 +                  const typename ColorMap::Value& first,
  3.2458 +                  const typename ColorMap::Value& second)
  3.2459 +        : _color_map(color_map), _first(first), _second(second) {}
  3.2460 +
  3.2461 +      Value operator[](const Key& key) const {
  3.2462 +        return _color_map[key] == _first || _color_map[key] == _second;
  3.2463 +      }
  3.2464 +
  3.2465 +    private:
  3.2466 +      const ColorMap& _color_map;
  3.2467 +      typename ColorMap::Value _first, _second;
  3.2468 +    };
  3.2469 +  }
  3.2470 +
  3.2471 +  /// \ingroup planar
  3.2472 +  ///
  3.2473 +  /// \brief Coloring planar graphs
  3.2474 +  ///
  3.2475 +  /// The graph coloring problem is the coloring of the graph nodes
  3.2476 +  /// that there are not adjacent nodes with the same color. The
  3.2477 +  /// planar graphs can be always colored with four colors, it is
  3.2478 +  /// proved by Appel and Haken and their proofs provide a quadratic
  3.2479 +  /// time algorithm for four coloring, but it could not be used to
  3.2480 +  /// implement efficient algorithm. The five and six coloring can be
  3.2481 +  /// made in linear time, but in this class the five coloring has
  3.2482 +  /// quadratic worst case time complexity. The two coloring (if
  3.2483 +  /// possible) is solvable with a graph search algorithm and it is
  3.2484 +  /// implemented in \ref bipartitePartitions() function in LEMON. To
  3.2485 +  /// decide whether the planar graph is three colorable is
  3.2486 +  /// NP-complete.
  3.2487 +  ///
  3.2488 +  /// This class contains member functions for calculate colorings
  3.2489 +  /// with five and six colors. The six coloring algorithm is a simple
  3.2490 +  /// greedy coloring on the backward minimum outgoing order of nodes.
  3.2491 +  /// This order can be computed as in each phase the node with least
  3.2492 +  /// outgoing arcs to unprocessed nodes is chosen. This order
  3.2493 +  /// guarantees that when a node is chosen for coloring it has at
  3.2494 +  /// most five already colored adjacents. The five coloring algorithm
  3.2495 +  /// use the same method, but if the greedy approach fails to color
  3.2496 +  /// with five colors, i.e. the node has five already different
  3.2497 +  /// colored neighbours, it swaps the colors in one of the connected
  3.2498 +  /// two colored sets with the Kempe recoloring method.
  3.2499 +  template <typename Graph>
  3.2500 +  class PlanarColoring {
  3.2501 +  public:
  3.2502 +
  3.2503 +    TEMPLATE_GRAPH_TYPEDEFS(Graph);
  3.2504 +
  3.2505 +    /// \brief The map type for store color indexes
  3.2506 +    typedef typename Graph::template NodeMap<int> IndexMap;
  3.2507 +    /// \brief The map type for store colors
  3.2508 +    typedef ComposeMap<Palette, IndexMap> ColorMap;
  3.2509 +
  3.2510 +    /// \brief Constructor
  3.2511 +    ///
  3.2512 +    /// Constructor
  3.2513 +    /// \pre The graph should be simple, i.e. loop and parallel arc free.
  3.2514 +    PlanarColoring(const Graph& graph)
  3.2515 +      : _graph(graph), _color_map(graph), _palette(0) {
  3.2516 +      _palette.add(Color(1,0,0));
  3.2517 +      _palette.add(Color(0,1,0));
  3.2518 +      _palette.add(Color(0,0,1));
  3.2519 +      _palette.add(Color(1,1,0));
  3.2520 +      _palette.add(Color(1,0,1));
  3.2521 +      _palette.add(Color(0,1,1));
  3.2522 +    }
  3.2523 +
  3.2524 +    /// \brief Returns the \e NodeMap of color indexes
  3.2525 +    ///
  3.2526 +    /// Returns the \e NodeMap of color indexes. The values are in the
  3.2527 +    /// range \c [0..4] or \c [0..5] according to the coloring method.
  3.2528 +    IndexMap colorIndexMap() const {
  3.2529 +      return _color_map;
  3.2530 +    }
  3.2531 +
  3.2532 +    /// \brief Returns the \e NodeMap of colors
  3.2533 +    ///
  3.2534 +    /// Returns the \e NodeMap of colors. The values are five or six
  3.2535 +    /// distinct \ref lemon::Color "colors".
  3.2536 +    ColorMap colorMap() const {
  3.2537 +      return composeMap(_palette, _color_map);
  3.2538 +    }
  3.2539 +
  3.2540 +    /// \brief Returns the color index of the node
  3.2541 +    ///
  3.2542 +    /// Returns the color index of the node. The values are in the
  3.2543 +    /// range \c [0..4] or \c [0..5] according to the coloring method.
  3.2544 +    int colorIndex(const Node& node) const {
  3.2545 +      return _color_map[node];
  3.2546 +    }
  3.2547 +
  3.2548 +    /// \brief Returns the color of the node
  3.2549 +    ///
  3.2550 +    /// Returns the color of the node. The values are five or six
  3.2551 +    /// distinct \ref lemon::Color "colors".
  3.2552 +    Color color(const Node& node) const {
  3.2553 +      return _palette[_color_map[node]];
  3.2554 +    }
  3.2555 +
  3.2556 +
  3.2557 +    /// \brief Calculates a coloring with at most six colors
  3.2558 +    ///
  3.2559 +    /// This function calculates a coloring with at most six colors. The time
  3.2560 +    /// complexity of this variant is linear in the size of the graph.
  3.2561 +    /// \return %True when the algorithm could color the graph with six color.
  3.2562 +    /// If the algorithm fails, then the graph could not be planar.
  3.2563 +    /// \note This function can return true if the graph is not
  3.2564 +    /// planar but it can be colored with 6 colors.
  3.2565 +    bool runSixColoring() {
  3.2566 +
  3.2567 +      typename Graph::template NodeMap<int> heap_index(_graph, -1);
  3.2568 +      BucketHeap<typename Graph::template NodeMap<int> > heap(heap_index);
  3.2569 +
  3.2570 +      for (NodeIt n(_graph); n != INVALID; ++n) {
  3.2571 +        _color_map[n] = -2;
  3.2572 +        heap.push(n, countOutArcs(_graph, n));
  3.2573 +      }
  3.2574 +
  3.2575 +      std::vector<Node> order;
  3.2576 +
  3.2577 +      while (!heap.empty()) {
  3.2578 +        Node n = heap.top();
  3.2579 +        heap.pop();
  3.2580 +        _color_map[n] = -1;
  3.2581 +        order.push_back(n);
  3.2582 +        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
  3.2583 +          Node t = _graph.runningNode(e);
  3.2584 +          if (_color_map[t] == -2) {
  3.2585 +            heap.decrease(t, heap[t] - 1);
  3.2586 +          }
  3.2587 +        }
  3.2588 +      }
  3.2589 +
  3.2590 +      for (int i = order.size() - 1; i >= 0; --i) {
  3.2591 +        std::vector<bool> forbidden(6, false);
  3.2592 +        for (OutArcIt e(_graph, order[i]); e != INVALID; ++e) {
  3.2593 +          Node t = _graph.runningNode(e);
  3.2594 +          if (_color_map[t] != -1) {
  3.2595 +            forbidden[_color_map[t]] = true;
  3.2596 +          }
  3.2597 +        }
  3.2598 +               for (int k = 0; k < 6; ++k) {
  3.2599 +          if (!forbidden[k]) {
  3.2600 +            _color_map[order[i]] = k;
  3.2601 +            break;
  3.2602 +          }
  3.2603 +        }
  3.2604 +        if (_color_map[order[i]] == -1) {
  3.2605 +          return false;
  3.2606 +        }
  3.2607 +      }
  3.2608 +      return true;
  3.2609 +    }
  3.2610 +
  3.2611 +  private:
  3.2612 +
  3.2613 +    bool recolor(const Node& u, const Node& v) {
  3.2614 +      int ucolor = _color_map[u];
  3.2615 +      int vcolor = _color_map[v];
  3.2616 +      typedef _planarity_bits::KempeFilter<IndexMap> KempeFilter;
  3.2617 +      KempeFilter filter(_color_map, ucolor, vcolor);
  3.2618 +
  3.2619 +      typedef FilterNodes<const Graph, const KempeFilter> KempeGraph;
  3.2620 +      KempeGraph kempe_graph(_graph, filter);
  3.2621 +
  3.2622 +      std::vector<Node> comp;
  3.2623 +      Bfs<KempeGraph> bfs(kempe_graph);
  3.2624 +      bfs.init();
  3.2625 +      bfs.addSource(u);
  3.2626 +      while (!bfs.emptyQueue()) {
  3.2627 +        Node n = bfs.nextNode();
  3.2628 +        if (n == v) return false;
  3.2629 +        comp.push_back(n);
  3.2630 +        bfs.processNextNode();
  3.2631 +      }
  3.2632 +
  3.2633 +      int scolor = ucolor + vcolor;
  3.2634 +      for (int i = 0; i < static_cast<int>(comp.size()); ++i) {
  3.2635 +        _color_map[comp[i]] = scolor - _color_map[comp[i]];
  3.2636 +      }
  3.2637 +
  3.2638 +      return true;
  3.2639 +    }
  3.2640 +
  3.2641 +    template <typename EmbeddingMap>
  3.2642 +    void kempeRecoloring(const Node& node, const EmbeddingMap& embedding) {
  3.2643 +      std::vector<Node> nodes;
  3.2644 +      nodes.reserve(4);
  3.2645 +
  3.2646 +      for (Arc e = OutArcIt(_graph, node); e != INVALID; e = embedding[e]) {
  3.2647 +        Node t = _graph.target(e);
  3.2648 +        if (_color_map[t] != -1) {
  3.2649 +          nodes.push_back(t);
  3.2650 +          if (nodes.size() == 4) break;
  3.2651 +        }
  3.2652 +      }
  3.2653 +
  3.2654 +      int color = _color_map[nodes[0]];
  3.2655 +      if (recolor(nodes[0], nodes[2])) {
  3.2656 +        _color_map[node] = color;
  3.2657 +      } else {
  3.2658 +        color = _color_map[nodes[1]];
  3.2659 +        recolor(nodes[1], nodes[3]);
  3.2660 +        _color_map[node] = color;
  3.2661 +      }
  3.2662 +    }
  3.2663 +
  3.2664 +  public:
  3.2665 +
  3.2666 +    /// \brief Calculates a coloring with at most five colors
  3.2667 +    ///
  3.2668 +    /// This function calculates a coloring with at most five
  3.2669 +    /// colors. The worst case time complexity of this variant is
  3.2670 +    /// quadratic in the size of the graph.
  3.2671 +    template <typename EmbeddingMap>
  3.2672 +    void runFiveColoring(const EmbeddingMap& embedding) {
  3.2673 +
  3.2674 +      typename Graph::template NodeMap<int> heap_index(_graph, -1);
  3.2675 +      BucketHeap<typename Graph::template NodeMap<int> > heap(heap_index);
  3.2676 +
  3.2677 +      for (NodeIt n(_graph); n != INVALID; ++n) {
  3.2678 +        _color_map[n] = -2;
  3.2679 +        heap.push(n, countOutArcs(_graph, n));
  3.2680 +      }
  3.2681 +
  3.2682 +      std::vector<Node> order;
  3.2683 +
  3.2684 +      while (!heap.empty()) {
  3.2685 +        Node n = heap.top();
  3.2686 +        heap.pop();
  3.2687 +        _color_map[n] = -1;
  3.2688 +        order.push_back(n);
  3.2689 +        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
  3.2690 +          Node t = _graph.runningNode(e);
  3.2691 +          if (_color_map[t] == -2) {
  3.2692 +            heap.decrease(t, heap[t] - 1);
  3.2693 +          }
  3.2694 +        }
  3.2695 +      }
  3.2696 +
  3.2697 +      for (int i = order.size() - 1; i >= 0; --i) {
  3.2698 +        std::vector<bool> forbidden(5, false);
  3.2699 +        for (OutArcIt e(_graph, order[i]); e != INVALID; ++e) {
  3.2700 +          Node t = _graph.runningNode(e);
  3.2701 +          if (_color_map[t] != -1) {
  3.2702 +            forbidden[_color_map[t]] = true;
  3.2703 +          }
  3.2704 +        }
  3.2705 +        for (int k = 0; k < 5; ++k) {
  3.2706 +          if (!forbidden[k]) {
  3.2707 +            _color_map[order[i]] = k;
  3.2708 +            break;
  3.2709 +          }
  3.2710 +        }
  3.2711 +        if (_color_map[order[i]] == -1) {
  3.2712 +          kempeRecoloring(order[i], embedding);
  3.2713 +        }
  3.2714 +      }
  3.2715 +    }
  3.2716 +
  3.2717 +    /// \brief Calculates a coloring with at most five colors
  3.2718 +    ///
  3.2719 +    /// This function calculates a coloring with at most five
  3.2720 +    /// colors. The worst case time complexity of this variant is
  3.2721 +    /// quadratic in the size of the graph.
  3.2722 +    /// \return %True when the graph is planar.
  3.2723 +    bool runFiveColoring() {
  3.2724 +      PlanarEmbedding<Graph> pe(_graph);
  3.2725 +      if (!pe.run()) return false;
  3.2726 +
  3.2727 +      runFiveColoring(pe.embeddingMap());
  3.2728 +      return true;
  3.2729 +    }
  3.2730 +
  3.2731 +  private:
  3.2732 +
  3.2733 +    const Graph& _graph;
  3.2734 +    IndexMap _color_map;
  3.2735 +    Palette _palette;
  3.2736 +  };
  3.2737 +
  3.2738 +}
  3.2739 +
  3.2740 +#endif
     4.1 --- a/lemon/unionfind.h	Fri Nov 13 00:39:28 2009 +0100
     4.2 +++ b/lemon/unionfind.h	Mon Dec 14 06:07:52 2009 +0100
     4.3 @@ -739,7 +739,7 @@
     4.4      /// Erase each item from the data structure.
     4.5      void clear() {
     4.6        items.clear();
     4.7 -      classes.clear;
     4.8 +      classes.clear();
     4.9        firstClass = firstFreeClass = firstFreeItem = -1;
    4.10      }
    4.11  
     5.1 --- a/test/CMakeLists.txt	Fri Nov 13 00:39:28 2009 +0100
     5.2 +++ b/test/CMakeLists.txt	Mon Dec 14 06:07:52 2009 +0100
     5.3 @@ -34,6 +34,7 @@
     5.4    min_cost_flow_test
     5.5    min_mean_cycle_test
     5.6    path_test
     5.7 +  planarity_test
     5.8    preflow_test
     5.9    radix_sort_test
    5.10    random_test
     6.1 --- a/test/Makefile.am	Fri Nov 13 00:39:28 2009 +0100
     6.2 +++ b/test/Makefile.am	Mon Dec 14 06:07:52 2009 +0100
     6.3 @@ -36,6 +36,7 @@
     6.4  	test/min_cost_flow_test \
     6.5  	test/min_mean_cycle_test \
     6.6  	test/path_test \
     6.7 +	test/planarity_test \
     6.8  	test/preflow_test \
     6.9  	test/radix_sort_test \
    6.10  	test/random_test \
    6.11 @@ -85,6 +86,7 @@
    6.12  test_min_cost_flow_test_SOURCES = test/min_cost_flow_test.cc
    6.13  test_min_mean_cycle_test_SOURCES = test/min_mean_cycle_test.cc
    6.14  test_path_test_SOURCES = test/path_test.cc
    6.15 +test_planarity_test_SOURCES = test/planarity_test.cc
    6.16  test_preflow_test_SOURCES = test/preflow_test.cc
    6.17  test_radix_sort_test_SOURCES = test/radix_sort_test.cc
    6.18  test_suurballe_test_SOURCES = test/suurballe_test.cc
     7.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     7.2 +++ b/test/planarity_test.cc	Mon Dec 14 06:07:52 2009 +0100
     7.3 @@ -0,0 +1,262 @@
     7.4 +/* -*- mode: C++; indent-tabs-mode: nil; -*-
     7.5 + *
     7.6 + * This file is a part of LEMON, a generic C++ optimization library.
     7.7 + *
     7.8 + * Copyright (C) 2003-2009
     7.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    7.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    7.11 + *
    7.12 + * Permission to use, modify and distribute this software is granted
    7.13 + * provided that this copyright notice appears in all copies. For
    7.14 + * precise terms see the accompanying LICENSE file.
    7.15 + *
    7.16 + * This software is provided "AS IS" with no warranty of any kind,
    7.17 + * express or implied, and with no claim as to its suitability for any
    7.18 + * purpose.
    7.19 + *
    7.20 + */
    7.21 +
    7.22 +#include <iostream>
    7.23 +
    7.24 +#include <lemon/planarity.h>
    7.25 +
    7.26 +#include <lemon/smart_graph.h>
    7.27 +#include <lemon/lgf_reader.h>
    7.28 +#include <lemon/connectivity.h>
    7.29 +#include <lemon/dim2.h>
    7.30 +
    7.31 +#include "test_tools.h"
    7.32 +
    7.33 +using namespace lemon;
    7.34 +using namespace lemon::dim2;
    7.35 +
    7.36 +const int lgfn = 4;
    7.37 +const std::string lgf[lgfn] = {
    7.38 +  "@nodes\n"
    7.39 +  "label\n"
    7.40 +  "0\n"
    7.41 +  "1\n"
    7.42 +  "2\n"
    7.43 +  "3\n"
    7.44 +  "4\n"
    7.45 +  "@edges\n"
    7.46 +  "     label\n"
    7.47 +  "0 1  0\n"
    7.48 +  "0 2  0\n"
    7.49 +  "0 3  0\n"
    7.50 +  "0 4  0\n"
    7.51 +  "1 2  0\n"
    7.52 +  "1 3  0\n"
    7.53 +  "1 4  0\n"
    7.54 +  "2 3  0\n"
    7.55 +  "2 4  0\n"
    7.56 +  "3 4  0\n",
    7.57 +
    7.58 +  "@nodes\n"
    7.59 +  "label\n"
    7.60 +  "0\n"
    7.61 +  "1\n"
    7.62 +  "2\n"
    7.63 +  "3\n"
    7.64 +  "4\n"
    7.65 +  "@edges\n"
    7.66 +  "     label\n"
    7.67 +  "0 1  0\n"
    7.68 +  "0 2  0\n"
    7.69 +  "0 3  0\n"
    7.70 +  "0 4  0\n"
    7.71 +  "1 2  0\n"
    7.72 +  "1 3  0\n"
    7.73 +  "2 3  0\n"
    7.74 +  "2 4  0\n"
    7.75 +  "3 4  0\n",
    7.76 +
    7.77 +  "@nodes\n"
    7.78 +  "label\n"
    7.79 +  "0\n"
    7.80 +  "1\n"
    7.81 +  "2\n"
    7.82 +  "3\n"
    7.83 +  "4\n"
    7.84 +  "5\n"
    7.85 +  "@edges\n"
    7.86 +  "     label\n"
    7.87 +  "0 3  0\n"
    7.88 +  "0 4  0\n"
    7.89 +  "0 5  0\n"
    7.90 +  "1 3  0\n"
    7.91 +  "1 4  0\n"
    7.92 +  "1 5  0\n"
    7.93 +  "2 3  0\n"
    7.94 +  "2 4  0\n"
    7.95 +  "2 5  0\n",
    7.96 +
    7.97 +  "@nodes\n"
    7.98 +  "label\n"
    7.99 +  "0\n"
   7.100 +  "1\n"
   7.101 +  "2\n"
   7.102 +  "3\n"
   7.103 +  "4\n"
   7.104 +  "5\n"
   7.105 +  "@edges\n"
   7.106 +  "     label\n"
   7.107 +  "0 3  0\n"
   7.108 +  "0 4  0\n"
   7.109 +  "0 5  0\n"
   7.110 +  "1 3  0\n"
   7.111 +  "1 4  0\n"
   7.112 +  "1 5  0\n"
   7.113 +  "2 3  0\n"
   7.114 +  "2 5  0\n"
   7.115 +};
   7.116 +
   7.117 +
   7.118 +
   7.119 +typedef SmartGraph Graph;
   7.120 +GRAPH_TYPEDEFS(Graph);
   7.121 +
   7.122 +typedef PlanarEmbedding<SmartGraph> PE;
   7.123 +typedef PlanarDrawing<SmartGraph> PD;
   7.124 +typedef PlanarColoring<SmartGraph> PC;
   7.125 +
   7.126 +void checkEmbedding(const Graph& graph, PE& pe) {
   7.127 +  int face_num = 0;
   7.128 +
   7.129 +  Graph::ArcMap<int> face(graph, -1);
   7.130 +
   7.131 +  for (ArcIt a(graph); a != INVALID; ++a) {
   7.132 +    if (face[a] == -1) {
   7.133 +      Arc b = a;
   7.134 +      while (face[b] == -1) {
   7.135 +        face[b] = face_num;
   7.136 +        b = pe.next(graph.oppositeArc(b));
   7.137 +      }
   7.138 +      check(face[b] == face_num, "Wrong face");
   7.139 +      ++face_num;
   7.140 +    }
   7.141 +  }
   7.142 +  check(face_num + countNodes(graph) - countConnectedComponents(graph) ==
   7.143 +        countEdges(graph) + 1, "Euler test does not passed");
   7.144 +}
   7.145 +
   7.146 +void checkKuratowski(const Graph& graph, PE& pe) {
   7.147 +  std::map<int, int> degs;
   7.148 +  for (NodeIt n(graph); n != INVALID; ++n) {
   7.149 +    int deg = 0;
   7.150 +    for (IncEdgeIt e(graph, n); e != INVALID; ++e) {
   7.151 +      if (pe.kuratowski(e)) {
   7.152 +        ++deg;
   7.153 +      }
   7.154 +    }
   7.155 +    ++degs[deg];
   7.156 +  }
   7.157 +  for (std::map<int, int>::iterator it = degs.begin(); it != degs.end(); ++it) {
   7.158 +    check(it->first == 0 || it->first == 2 ||
   7.159 +          (it->first == 3 && it->second == 6) ||
   7.160 +          (it->first == 4 && it->second == 5),
   7.161 +          "Wrong degree in Kuratowski graph");
   7.162 +  }
   7.163 +
   7.164 +  // Not full test
   7.165 +  check((degs[3] == 0) != (degs[4] == 0), "Wrong Kuratowski graph");
   7.166 +}
   7.167 +
   7.168 +bool intersect(Point<int> e1, Point<int> e2, Point<int> f1, Point<int> f2) {
   7.169 +  int l, r;
   7.170 +  if (std::min(e1.x, e2.x) > std::max(f1.x, f2.x)) return false;
   7.171 +  if (std::max(e1.x, e2.x) < std::min(f1.x, f2.x)) return false;
   7.172 +  if (std::min(e1.y, e2.y) > std::max(f1.y, f2.y)) return false;
   7.173 +  if (std::max(e1.y, e2.y) < std::min(f1.y, f2.y)) return false;
   7.174 +
   7.175 +  l = (e2.x - e1.x) * (f1.y - e1.y) - (e2.y - e1.y) * (f1.x - e1.x);
   7.176 +  r = (e2.x - e1.x) * (f2.y - e1.y) - (e2.y - e1.y) * (f2.x - e1.x);
   7.177 +  if (!((l >= 0 && r <= 0) || (l <= 0 && r >= 0))) return false;
   7.178 +  l = (f2.x - f1.x) * (e1.y - f1.y) - (f2.y - f1.y) * (e1.x - f1.x);
   7.179 +  r = (f2.x - f1.x) * (e2.y - f1.y) - (f2.y - f1.y) * (e2.x - f1.x);
   7.180 +  if (!((l >= 0 && r <= 0) || (l <= 0 && r >= 0))) return false;
   7.181 +  return true;
   7.182 +}
   7.183 +
   7.184 +bool collinear(Point<int> p, Point<int> q, Point<int> r) {
   7.185 +  int v;
   7.186 +  v = (q.x - p.x) * (r.y - p.y) - (q.y - p.y) * (r.x - p.x);
   7.187 +  if (v != 0) return false;
   7.188 +  v = (q.x - p.x) * (r.x - p.x) + (q.y - p.y) * (r.y - p.y);
   7.189 +  if (v < 0) return false;
   7.190 +  return true;
   7.191 +}
   7.192 +
   7.193 +void checkDrawing(const Graph& graph, PD& pd) {
   7.194 +  for (Graph::NodeIt n(graph); n != INVALID; ++n) {
   7.195 +    Graph::NodeIt m(n);
   7.196 +    for (++m; m != INVALID; ++m) {
   7.197 +      check(pd[m] != pd[n], "Two nodes with identical coordinates");
   7.198 +    }
   7.199 +  }
   7.200 +
   7.201 +  for (Graph::EdgeIt e(graph); e != INVALID; ++e) {
   7.202 +    for (Graph::EdgeIt f(e); f != e; ++f) {
   7.203 +      Point<int> e1 = pd[graph.u(e)];
   7.204 +      Point<int> e2 = pd[graph.v(e)];
   7.205 +      Point<int> f1 = pd[graph.u(f)];
   7.206 +      Point<int> f2 = pd[graph.v(f)];
   7.207 +
   7.208 +      if (graph.u(e) == graph.u(f)) {
   7.209 +        check(!collinear(e1, e2, f2), "Wrong drawing");
   7.210 +      } else if (graph.u(e) == graph.v(f)) {
   7.211 +        check(!collinear(e1, e2, f1), "Wrong drawing");
   7.212 +      } else if (graph.v(e) == graph.u(f)) {
   7.213 +        check(!collinear(e2, e1, f2), "Wrong drawing");
   7.214 +      } else if (graph.v(e) == graph.v(f)) {
   7.215 +        check(!collinear(e2, e1, f1), "Wrong drawing");
   7.216 +      } else {
   7.217 +        check(!intersect(e1, e2, f1, f2), "Wrong drawing");
   7.218 +      }
   7.219 +    }
   7.220 +  }
   7.221 +}
   7.222 +
   7.223 +void checkColoring(const Graph& graph, PC& pc, int num) {
   7.224 +  for (NodeIt n(graph); n != INVALID; ++n) {
   7.225 +    check(pc.colorIndex(n) >= 0 && pc.colorIndex(n) < num,
   7.226 +          "Wrong coloring");
   7.227 +  }
   7.228 +  for (EdgeIt e(graph); e != INVALID; ++e) {
   7.229 +    check(pc.colorIndex(graph.u(e)) != pc.colorIndex(graph.v(e)),
   7.230 +          "Wrong coloring");
   7.231 +  }
   7.232 +}
   7.233 +
   7.234 +int main() {
   7.235 +
   7.236 +  for (int i = 0; i < lgfn; ++i) {
   7.237 +    std::istringstream lgfs(lgf[i]);
   7.238 +
   7.239 +    SmartGraph graph;
   7.240 +    graphReader(graph, lgfs).run();
   7.241 +
   7.242 +    check(simpleGraph(graph), "Test graphs must be simple");
   7.243 +
   7.244 +    PE pe(graph);
   7.245 +    bool planar = pe.run();
   7.246 +    check(checkPlanarity(graph) == planar, "Planarity checking failed");
   7.247 +
   7.248 +    if (planar) {
   7.249 +      checkEmbedding(graph, pe);
   7.250 +
   7.251 +      PlanarDrawing<Graph> pd(graph);
   7.252 +      pd.run(pe.embeddingMap());
   7.253 +      checkDrawing(graph, pd);
   7.254 +
   7.255 +      PlanarColoring<Graph> pc(graph);
   7.256 +      pc.runFiveColoring(pe.embeddingMap());
   7.257 +      checkColoring(graph, pc, 5);
   7.258 +
   7.259 +    } else {
   7.260 +      checkKuratowski(graph, pe);
   7.261 +    }
   7.262 +  }
   7.263 +
   7.264 +  return 0;
   7.265 +}