# HG changeset patch
# User gyorokp
# Date 1269025969 -3600
# Node ID eb9a0cb25eea358d39a0d2cecc4722c919ff1b7f
# Parent  4266cccd150b8bac2ae31fe1a9b862e73786e95f
Got rid of PlanarDigraph
Tidied code, added doxygen comments

diff -r 4266cccd150b -r eb9a0cb25eea lemon/planar_graph.h
--- a/lemon/planar_graph.h	Thu Mar 18 21:56:09 2010 +0100
+++ b/lemon/planar_graph.h	Fri Mar 19 20:12:49 2010 +0100
@@ -21,7 +21,7 @@
 
 ///\ingroup graphs
 ///\file
-///\brief ListDigraph and ListGraph classes.
+///\brief PlanarGraph classes. UNDER CONSTRUCTION.
 
 #include <lemon/core.h>
 #include <lemon/error.h>
@@ -36,2741 +36,1495 @@
 
 namespace lemon {
 
-    class PlanarDigraph;
+  class PlanarGraphBase {
 
-    class PlanarDigraphBase {
+  protected:
 
+    struct NodeT {
+      int first_out, last_out;
+      int prev, next;
+      int component;
+    };
+
+    struct ArcT {
+      int target;
+      int left_face;
+      int prev_out, next_out;
+    };
+
+    struct FaceT {
+      int first_arc;
+      int prev, next;
+    };
+
+    std::vector<NodeT> nodes;
+
+    int first_node;
+
+    int first_free_node;
+
+    std::vector<ArcT> arcs;
+
+    int first_free_arc;
+
+    std::vector<FaceT> faces;
+
+    int first_face;
+
+    int first_free_face;
+
+    int component_id;
+
+  public:
+
+    typedef PlanarGraphBase Graph;
+
+    class Node {
+      friend class PlanarGraphBase;
     protected:
-        struct NodeT {
-            int first_in, first_out;
-            int first_ccw, last_ccw;
-            int prev, next;
-            int component;
-        };
 
-        struct ArcT {
-            int target, source;
-            int prev_in, prev_out;
-            int next_in, next_out;
+      int id;
+      explicit Node(int pid) {
+        id = pid;
+      }
 
-            int prev_s, prev_t;
-            int next_s, next_t;
-            int left_face, right_face;
+    public:
+      Node() {}
+      Node (Invalid) {
+        id = -1;
+      }
+      bool operator==(const Node& node) const {
+        return id == node.id;
+      }
+      bool operator!=(const Node& node) const {
+        return id != node.id;
+      }
+      bool operator<(const Node& node) const {
+        return id < node.id;
+      }
+    };
 
-            int target_at(bool d) const {
-                return d?target:source;
+    class Edge {
+      friend class PlanarGraphBase;
+    protected:
+
+      int id;
+      explicit Edge(int pid) {
+        id = pid;
+      }
+
+    public:
+      Edge() {}
+      Edge (Invalid) {
+        id = -1;
+      }
+      bool operator==(const Edge& edge) const {
+        return id == edge.id;
+      }
+      bool operator!=(const Edge& edge) const {
+        return id != edge.id;
+      }
+      bool operator<(const Edge& edge) const {
+        return id < edge.id;
+      }
+    };
+
+    class Arc {
+      friend class PlanarGraphBase;
+    protected:
+
+      int id;
+      explicit Arc(int pid) {
+        id = pid;
+      }
+
+    public:
+      operator Edge() const {
+        return id != -1 ? edgeFromId(id / 2) : INVALID;
+      }
+
+      Arc() {}
+      Arc (Invalid) {
+        id = -1;
+      }
+      bool operator==(const Arc& arc) const {
+        return id == arc.id;
+      }
+      bool operator!=(const Arc& arc) const {
+        return id != arc.id;
+      }
+      bool operator<(const Arc& arc) const {
+        return id < arc.id;
+      }
+    };
+
+    class Face {
+      friend class PlanarGraphBase;
+    protected:
+
+      int id;
+      explicit Face(int pid) {
+        id = pid;
+      }
+
+    public:
+      Face() {}
+      Face (Invalid) {
+        id = -1;
+      }
+      bool operator==(const Face& face) const {
+        return id == face.id;
+      }
+      bool operator!=(const Face& face) const {
+        return id != face.id;
+      }
+      bool operator<(const Face& face) const {
+        return id < face.id;
+      }
+    };
+
+    PlanarGraphBase()
+        : nodes(), first_node(-1), first_free_node(-1),
+        arcs(), first_free_arc(-1),
+        faces(), first_face(-1), first_free_face(-1),
+        component_id(0) {
+    }
+
+
+    int maxNodeId() const {
+      return nodes.size()-1;
+    }
+    int maxEdgeId() const {
+      return arcs.size() / 2 - 1;
+    }
+    int maxArcId() const {
+      return arcs.size()-1;
+    }
+    int maxFaceId() const {
+      return faces.size()-1;
+    }
+
+    Node source(Arc e) const {
+      return Node(arcs[e.id ^ 1].target);
+    }
+    Node target(Arc e) const {
+      return Node(arcs[e.id].target);
+    }
+    Face leftFace(Arc e) const {
+      return Face(arcs[e.id].left_face);
+    }
+    Face rightFace(Arc e) const {
+      return Face(arcs[e.id ^ 1].left_face);
+    }
+
+    Node u(Edge e) const {
+      return Node(arcs[2 * e.id].target);
+    }
+    Node v(Edge e) const {
+      return Node(arcs[2 * e.id + 1].target);
+    }
+    Face w1(Edge e) const {
+      return Face(arcs[2 * e.id].left_face);
+    }
+    Face w2(Edge e) const {
+      return Face(arcs[2 * e.id + 1].left_face);
+    }
+
+    static bool direction(Arc e) {
+      return (e.id & 1) == 1;
+    }
+
+    static Arc direct(Edge e, bool d) {
+      return Arc(e.id * 2 + (d ? 1 : 0));
+    }
+
+    //Primitives to use in iterators
+    void first(Node& node) const {
+      node.id = first_node;
+    }
+
+    void next(Node& node) const {
+      node.id = nodes[node.id].next;
+    }
+
+    void first(Arc& e) const {
+      int n = first_node;
+      while (n != -1 && nodes[n].first_out == -1) {
+        n = nodes[n].next;
+      }
+      e.id = (n == -1) ? -1 : nodes[n].first_out;
+    }
+
+    void next(Arc& e) const {
+      if (arcs[e.id].next_out != -1) {
+        e.id = arcs[e.id].next_out;
+      } else {
+        int n = nodes[arcs[e.id ^ 1].target].next;
+        while (n != -1 && nodes[n].first_out == -1) {
+          n = nodes[n].next;
+        }
+        e.id = (n == -1) ? -1 : nodes[n].first_out;
+      }
+    }
+
+    void first(Edge& e) const {
+      int n = first_node;
+      while (n != -1) {
+        e.id = nodes[n].first_out;
+        while ((e.id & 1) != 1) {
+          e.id = arcs[e.id].next_out;
+        }
+        if (e.id != -1) {
+          e.id /= 2;
+          return;
+        }
+        n = nodes[n].next;
+      }
+      e.id = -1;
+    }
+
+    void next(Edge& e) const {
+      int n = arcs[e.id * 2].target;
+      e.id = arcs[(e.id * 2) | 1].next_out;
+      while ((e.id & 1) != 1) {
+        e.id = arcs[e.id].next_out;
+      }
+      if (e.id != -1) {
+        e.id /= 2;
+        return;
+      }
+      n = nodes[n].next;
+      while (n != -1) {
+        e.id = nodes[n].first_out;
+        while ((e.id & 1) != 1) {
+          e.id = arcs[e.id].next_out;
+        }
+        if (e.id != -1) {
+          e.id /= 2;
+          return;
+        }
+        n = nodes[n].next;
+      }
+      e.id = -1;
+    }
+
+    void firstOut(Arc &e, const Node& v) const {
+      e.id = nodes[v.id].first_out;
+    }
+    void nextOut(Arc &e) const {
+      e.id = arcs[e.id].next_out;
+    }
+
+    void firstIn(Arc &e, const Node& v) const {
+      e.id = ((nodes[v.id].first_out) ^ 1);
+      if (e.id == -2) e.id = -1;
+    }
+    void nextIn(Arc &e) const {
+      e.id = ((arcs[e.id ^ 1].next_out) ^ 1);
+      if (e.id == -2) e.id = -1;
+    }
+    void lastIn(Arc &e, const Node& v) const {
+      e.id = ((nodes[v.id].last_out) ^ 1);
+      if (e.id == -2) e.id = -1;
+    }
+    void prevIn(Arc &e) const {
+      e.id = ((arcs[e.id ^ 1].prev_out) ^ 1);
+      if (e.id == -2) e.id = -1;
+    }
+
+    void firstCwF(Arc &e, const Face &f) const {
+      e.id = faces[f.id].first_arc;
+    }
+    void nextCwF(Arc &e) const {
+      turnLeft(e);
+      setToOpposite(e);
+      if (e.id == faces[arcs[e.id].left_face].first_arc) e = INVALID;
+    }
+
+    void firstInF(Arc &e, const Face &f) const {
+      e.id = faces[f.id].first_arc;
+    }
+    void nextInF(Arc &e) const {
+      setToOpposite(e);
+      turnRight(e);
+      if (e.id == faces[arcs[e.id].left_face].first_arc) e = INVALID;
+    }
+    void lastInF(Arc &e, const Face &f) const {
+      e.id = faces[f.id].first_arc;
+      setToOpposite(e);
+      turnRight(e);
+    }
+    void prevInF(Arc &e) const {
+      if (e.id == faces[arcs[e.id].left_face].first_arc) {
+        e = INVALID;
+        return;
+      }
+      setToOpposite(e);
+      turnRight(e);
+    }
+
+    void firstOutF(Arc &e, const Face &f) const {
+      e.id = faces[e.id].first_arc ^ 1;
+    }
+    void nextOutF(Arc &e) const {
+      turnRight(e);
+      setToOpposite(e);
+      if (e.id == faces[arcs[e.id ^ 1].left_face].first_arc) e = INVALID;
+    }
+
+    void first(Arc &arc, const Face& face) const {
+      arc.id = faces[face.id].first_arc;
+    }
+
+    void turnLeftF(Arc &e) const {
+      setToOpposite(e);
+      turnRight(e);
+    }
+
+    void turnRightF(Arc &e) const {
+      turnLeft(e);
+      setToOpposite(e);
+    }
+
+    void turnLeft(Arc &e) const {
+      if (arcs[e.id].next_out > -1) {
+        e.id = arcs[e.id].next_out;
+      } else {
+        e.id = nodes[arcs[e.id ^ 1].target].first_out;
+      }
+    }
+    void turnRight(Arc &e) const {
+      if (arcs[e.id].prev_out > -1) {
+        e.id = arcs[e.id].prev_out;
+      } else {
+        e.id = nodes[arcs[e.id ^ 1].target].last_out;
+      }
+    }
+    void setToOpposite(Arc &a) const {
+      if (a.id != -1) a.id ^= 1;
+    }
+
+    void firstInc(Edge &e, bool& d, const Node& v) const {
+      int a = nodes[v.id].first_out;
+      if (a != -1 ) {
+        e.id = a / 2;
+        d = ((a & 1) == 1);
+      } else {
+        e.id = -1;
+        d = true;
+      }
+    }
+    void nextInc(Edge &e, bool& d) const {
+      int a = (arcs[(e.id * 2) | (d ? 1 : 0)].next_out);
+      if (a != -1 ) {
+        e.id = a / 2;
+        d = ((a & 1) == 1);
+      } else {
+        e.id = -1;
+        d = true;
+      }
+    }
+
+    void first(Face& face) const {
+      face.id = first_face;
+    }
+
+    void next(Face& face) const {
+      face.id = faces[face.id].next;
+    }
+
+    Arc arcAt(Node n, Edge e) {
+      if (e.id == -1) return INVALID;
+      return Arc((e.id*2) | (arcs[e.id*2].target == n.id?1:0));
+    }
+
+    static int id(Node v) {
+      return v.id;
+    }
+    static int id(Arc e) {
+      return e.id;
+    }
+    static int id(Edge e) {
+      return e.id;
+    }
+    static int id(Face f) {
+      return f.id;
+    }
+
+    static Node nodeFromId(int id) {
+      return Node(id);
+    }
+    static Arc arcFromId(int id) {
+      return Arc(id);
+    }
+    static Edge edgeFromId(int id) {
+      return Edge(id);
+    }
+    static Face faceFromId(int id) {
+      return Face(id);
+    }
+
+    bool valid(Node n) const {
+      return n.id >= 0 && n.id < static_cast<int>(nodes.size()) &&
+             nodes[n.id].prev != -2;
+    }
+
+    bool valid(Arc a) const {
+      return a.id >= 0 && a.id < static_cast<int>(arcs.size()) &&
+             arcs[a.id].prev_out != -2;
+    }
+
+    bool valid(Edge e) const {
+      return e.id >= 0 && 2 * e.id < static_cast<int>(arcs.size()) &&
+             arcs[2 * e.id].prev_out != -2;
+    }
+
+    bool valid(Face f) const {
+      return f.id >= 0 && f.id < static_cast<int>(faces.size()) &&
+             faces[f.id].prev != -2;
+    }
+
+    Node addNode() {
+      int n;
+
+      if (first_free_node==-1) {
+        n = nodes.size();
+        nodes.push_back(NodeT());
+      } else {
+        n = first_free_node;
+        first_free_node = nodes[n].next;
+      }
+
+      nodes[n].next = first_node;
+      nodes[n].component = component_id++;
+      if (first_node != -1) nodes[first_node].prev = n;
+      first_node = n;
+      nodes[n].prev = -1;
+
+      nodes[n].first_out = nodes[n].last_out = -1;
+
+      return Node(n);
+    }
+
+    Edge addEdge(Node u, Node v, Edge e_u, Edge e_v) {
+
+      Arc p_u = arcAt(u,e_u);
+      Arc p_v = arcAt(v,e_v);
+
+      if (p_u.id > -1 && p_v.id > -1 && arcs[p_u.id].left_face != arcs[p_v.id].
+          left_face && nodes[u.id].component == nodes[v.id].component) return
+          INVALID;
+
+      int n = addBlankEdge();
+
+      arcs[n].target = u.id;
+      arcs[n | 1].target = v.id;
+
+      arcs[n].prev_out = p_v.id;
+      if (p_v.id > -1) {
+        arcs[n].next_out = arcs[p_v.id].next_out;
+        arcs[p_v.id].next_out = n;
+      } else {
+        arcs[n].next_out = nodes[v.id].first_out;
+        nodes[v.id].first_out = n;
+      }
+      if (arcs[n].next_out > -1) {
+        arcs[arcs[n].next_out].prev_out = n;
+      } else {
+        nodes[v.id].last_out = n;
+      }
+
+      arcs[n | 1].prev_out = p_u.id;
+      if (p_u.id > -1) {
+        arcs[n | 1].next_out = arcs[p_u.id].next_out;
+        arcs[p_u.id].next_out = n | 1;
+      } else {
+        arcs[n | 1].next_out = nodes[u.id].first_out;
+        nodes[u.id].first_out = n | 1;
+      }
+      if (arcs[n | 1].next_out > -1) {
+        arcs[arcs[n | 1].next_out].prev_out = n | 1;
+      } else {
+        nodes[u.id].last_out = n | 1;
+      }
+
+      //Add the extra face, if needed
+      if (p_u.id > -1 & p_v.id > -1) {
+        int oldf = arcs[p_u.id].left_face;
+        int oldfb = arcs[p_v.id].left_face;
+        arcs[n].left_face = arcs[n | 1].left_face = oldf;
+        Face f = addFace();
+        faces[f.id].first_arc = n | 1;
+        faces[oldf].first_arc = n;
+        Arc arc(n | 1);
+        wall_paint(arc,f.id,arc);
+        if (nodes[v.id].component != nodes[u.id].component) {
+          erase(Face(oldf));
+          erase(Face(oldfb));
+          int ca = nodes[u.id].component;
+          int cb = nodes[v.id].component;
+          int k = first_node;
+          while (k != -1) {
+            if (nodes[k].component == cb)
+              nodes[k].component = ca;
+            k = nodes[k].next;
+          }
+        }
+      } else if (p_u.id > -1) {
+        arcs[n].left_face = arcs[n | 1].left_face = arcs[p_u.id].left_face;
+        faces[arcs[n].left_face].first_arc = n | 1;
+        nodes[v.id].component = nodes[u.id].component;
+      } else if (p_v.id > -1) {
+        arcs[n].left_face = arcs[n | 1].left_face = arcs[p_v.id].left_face;
+        faces[arcs[n].left_face].first_arc = n | 1;
+        nodes[u.id].component = nodes[v.id].component;
+      } else {    //both prevs are INVALID
+        Face f = addFace();
+        arcs[n].left_face = arcs[n | 1].left_face = f.id;
+        faces[f.id].first_arc = n | 1;
+        nodes[v.id].component = nodes[u.id].component;
+      }
+
+      return Edge(n / 2);
+    }
+
+    void erase(const Node& node) {
+      int n = node.id;
+
+      if (nodes[n].next != -1) {
+        nodes[nodes[n].next].prev = nodes[n].prev;
+      }
+
+      if (nodes[n].prev != -1) {
+        nodes[nodes[n].prev].next = nodes[n].next;
+      } else {
+        first_node = nodes[n].next;
+      }
+
+      nodes[n].next = first_free_node;
+      first_free_node = n;
+      nodes[n].prev = -2;
+    }
+
+    void erase(const Edge& edge) {
+      int n = edge.id * 2;
+
+      //"retreat" the incident faces' first arcs
+      int fl = arcs[n].left_face;
+      if ((faces[fl].first_arc | 1) == (n | 1)) {
+        Arc e(faces[fl].first_arc);
+        turnRightF(e);
+        if ((e.id | 1) == (n | 1)) turnRightF(e);
+        faces[fl].first_arc = e.id;
+      }
+
+      int fr = arcs[n | 1].left_face;
+
+      bool comp_split = false;
+      if (fr != fl) {
+        Arc arc(faces[fr].first_arc);
+        wall_paint(arc,fl,arc);
+        erase(Face(fr));
+      } else if ((arcs[n].next_out > -1 || arcs[n].prev_out > -1) &&
+                 (arcs[n | 1].next_out > -1 || arcs[n | 1].prev_out > -1)) {
+        comp_split = true;
+        Arc arc(n);
+        Arc ed = arc;
+        ed.id ^= 1;
+        turnRightF(arc);
+        Face f = addFace();
+        wall_paint(arc,f.id,ed);
+        faces[f.id].first_arc = arc.id;
+      }
+
+      if (arcs[n].next_out != -1) {
+        arcs[arcs[n].next_out].prev_out = arcs[n].prev_out;
+      } else {
+        nodes[arcs[n].target].last_out = arcs[n].prev_out;
+      }
+
+      if (arcs[n].prev_out != -1) {
+        arcs[arcs[n].prev_out].next_out = arcs[n].next_out;
+      } else {
+        nodes[arcs[n | 1].target].first_out = arcs[n].next_out;
+      }
+
+      if (arcs[n | 1].next_out != -1) {
+        arcs[arcs[n | 1].next_out].prev_out = arcs[n | 1].prev_out;
+      } else {
+        nodes[arcs[n | 1].target].last_out = arcs[n | 1].prev_out;
+      }
+
+      if (arcs[n | 1].prev_out != -1) {
+        arcs[arcs[n | 1].prev_out].next_out = arcs[n | 1].next_out;
+      } else {
+        nodes[arcs[n].target].first_out = arcs[n | 1].next_out;
+      }
+
+      arcs[n].next_out = first_free_arc;
+      first_free_arc = n;
+      arcs[n].prev_out = -2;
+      arcs[n | 1].prev_out = -2;
+
+      if (comp_split) component_relabel(Node(arcs[n | 1].target),
+        component_id++);
+
+    }
+
+    void clear() {
+      arcs.clear();
+      nodes.clear();
+      faces.clear();
+      first_node = first_free_node = first_free_arc = first_face =
+        first_free_face = -1;
+    }
+
+    Node split(Edge e) {
+      Node v = addNode();
+      Arc a(e.id*2);
+      int b = addBlankEdge();
+
+      nodes[v.id].component = nodes[arcs[a.id].target].component;
+      nodes[v.id].first_out = a.id;
+      nodes[v.id].last_out = b | 1;
+
+      arcs[b] = arcs[a.id];
+      arcs[b].target = v.id;
+      if (arcs[a.id].next_out > -1)
+        arcs[arcs[a.id].next_out].prev_out = b;
+      else
+        nodes[arcs[a.id | 1].target].last_out = b;
+      if (arcs[a.id].prev_out > -1)
+        arcs[arcs[a.id].prev_out].next_out = b;
+      else
+        nodes[arcs[a.id | 1].target].first_out = b;
+
+      arcs[b | 1] = arcs[a.id | 1];
+      arcs[b | 1].next_out = -1;
+      arcs[b | 1].prev_out = a.id;
+
+      arcs[a.id].next_out = b | 1;
+      arcs[a.id].prev_out = -1;
+      arcs[a.id | 1].target = v.id;
+
+
+      return v;
+    }
+
+  protected:
+
+    void wall_paint(Arc arc, int f_id, Arc ed) {
+      do {
+        arcs[arc.id].left_face = f_id;
+        turnRightF(arc);
+      } while (arc != ed);
+    }
+
+    void component_relabel(Node node, int comp_id) {
+      std::vector<int> ns(nodes.size());
+      std::list<int> q;
+      q.push_back(node.id);
+      ns[node.id] = 1;
+      while (!q.empty()) {
+        int n = q.front();
+        ns[n] = 2;
+        nodes[n].component = comp_id;
+        q.pop_front();
+        Arc arc;
+        firstOut(arc,Node(n));
+        while (arc.id > -1) {
+          int m = arcs[arc.id].target;
+          if (ns[m] == 0) {
+            ns[m] = 1;
+            q.push_back(m);
+          }
+          nextOut(arc);
+        }
+      }
+    }
+
+    Face addFace() {
+      int n;
+
+      if (first_free_face==-1) {
+        n = faces.size();
+        faces.push_back(FaceT());
+      } else {
+        n = first_free_face;
+        first_free_face = faces[n].next;
+      }
+
+      faces[n].next = first_face;
+      if (first_face != -1) faces[first_face].prev = n;
+      first_face = n;
+      faces[n].prev = -1;
+
+      faces[n].first_arc = -1;
+
+      return Face(n);
+    }
+
+    void erase(const Face& face) {
+      int n = face.id;
+
+      if (faces[n].next != -1) {
+        faces[faces[n].next].prev = faces[n].prev;
+      }
+
+      if (faces[n].prev != -1) {
+        faces[faces[n].prev].next = faces[n].next;
+      } else {
+        first_face = faces[n].next;
+      }
+
+      faces[n].next = first_free_face;
+      first_free_face = n;
+      faces[n].prev = -2;
+    }
+
+    int addBlankEdge() {
+      int n;
+      if (first_free_arc == -1) {
+        n = arcs.size();
+        arcs.push_back(ArcT());
+        arcs.push_back(ArcT());
+      } else {
+        n = first_free_arc;
+        first_free_arc = arcs[n].next_out;
+      }
+      return n;
+    }
+
+#ifdef REMOVE_BEFORE_RELEASE
+  public:
+    void print() {
+      std::cout << "Nodes: " << std::endl;
+      for (int i=0; i<nodes.size(); ++i)
+        std::cout << i << ":"
+        << " fo=" << nodes[i].first_out
+        << " pr=" << nodes[i].prev
+        << " nx=" << nodes[i].next
+        << " co=" << nodes[i].component
+        <<std::endl;
+      std::cout << "Arcs: " << std::endl;
+      for (int i=0; i<arcs.size(); ++i) {
+        if (arcs[i].next_out > -2) {
+          std::cout << i << ":"
+          << " tg=" << arcs[i].target
+          << " po=" << arcs[i].prev_out
+          << " no=" << arcs[i].next_out
+          << " lf=" << arcs[i].left_face
+          <<std::endl;
+        } else std::cout << i << ": (deleted)" << std::endl;
+      }
+      std::cout << "Faces: " << std::endl;
+      for (int i=0; i<faces.size(); ++i)
+        std::cout << i
+        << " pr=" << faces[i].prev
+        << " nx=" << faces[i].next
+        << " fa=" << faces[i].first_arc
+        <<std::endl;
+    }
+#endif
+
+  };
+
+  template<typename Base>
+  class PlanarGraphExtender : public Base{
+
+    typedef Base Parent;
+
+  public:
+    typedef PlanarGraphExtender Graph;
+
+    PlanarGraphExtender() {}
+
+    typedef typename Parent::Node Node;
+    typedef typename Parent::Arc Arc;
+    typedef typename Parent::Edge Edge;
+    typedef typename Parent::Face Face;
+
+  /// Iterator class for the faces.
+
+  /// This iterator goes through the faces of a planar graph.
+  class FaceIt : public Face {
+      const Graph* _graph;
+    public:
+
+      FaceIt() {}
+
+      FaceIt(Invalid i) : Face(i) { }
+
+      explicit FaceIt(const Graph& graph) : _graph(&graph) {
+        _graph->first(static_cast<Face&>(*this));
+      }
+
+      FaceIt(const Graph& graph, const Face& face)
+          : Face(face), _graph(&graph) {}
+
+      FaceIt& operator++() {
+        _graph->next(*this);
+        return *this;
+      }
+
+    };
+
+
+  /// Iterator class for the arcs on the boundary of a face.
+
+  /// This iterator goes through the arcs on the boundary of a face clockwise.
+  class CwBoundaryArcIt : public Arc {
+      const Graph* _graph;
+      Face _face;
+      Arc f_arc;
+    public:
+
+      CwBoundaryArcIt() { }
+
+      CwBoundaryArcIt(Invalid i) : Arc(i) { }
+
+      CwBoundaryArcIt(const Graph& graph, const Face& face)
+          : _graph(&graph), _face(face) {
+        _graph->firstCwF(static_cast<Arc&>(*this),face);
+        f_arc = *this;
+      }
+
+      CwBoundaryArcIt(const Graph& graph, const Arc& arc)
+          : Arc(arc), _graph(&graph) {}
+
+      CwBoundaryArcIt& operator++() {
+        _graph->nextCwF(*this);
+        return *this;
+      }
+
+    };
+
+  /// Iterator class for the arcs around a node.
+
+  /// This iterator goes through the arcs around a node anticlockwise.
+  class CcwArcIt : public Arc {
+      const Graph* _graph;
+      const Node _node;
+    public:
+
+      CcwArcIt() { }
+
+      CcwArcIt(Invalid i) : Arc(i) { }
+
+      CcwArcIt(const Graph& graph, const Node& node)
+          : _graph(&graph), _node(node) {
+        _graph->firstCcw(*this, node);
+      }
+
+      CcwArcIt(const Graph& graph, const Arc& arc)
+          : Arc(arc), _graph(&graph) {}
+
+      CcwArcIt& operator++() {
+        _graph->nextCcw(*this, _node);
+        return *this;
+      }
+
+    };
+
+  };
+
+  typedef PlanarGraphExtender<GraphExtender<PlanarGraphBase> >
+    ExtendedPlanarGraphBase;
+
+
+  /// \addtogroup graphs
+  /// @{
+
+  ///An undirected planar graph structure.
+
+  ///\ref PlanarGraph is a versatile and fast undirected planar graph
+  ///implementation based on linked lists that are stored in
+  ///\c std::vector structures. It maintains a topology of nodes, edges
+  ///and faces.
+  ///
+  ///This type fully conforms to the \ref concepts::Graph "Graph concept"
+  ///and it also provides several useful additional functionalities, including
+  ///those specific to planar graphs.
+  ///Most of its member functions and nested classes are documented
+  ///only in the concept class.
+  ///
+  ///This class provides only linear time counting for nodes, edges, arcs and
+  ///faces.
+  ///
+  ///A disconnected planar graph has have an outer face for each of its
+  ///components, effectively turning them into separate graphs. Each component
+  ///has a corresponding component in the dual.
+  ///\sa concepts::Graph
+  ///\sa PlaneGraph
+  class PlanarGraph : public ExtendedPlanarGraphBase {
+    typedef ExtendedPlanarGraphBase Parent;
+
+  public:
+    /// Graphs are \e not copy constructible. Use GraphCopy instead.
+    PlanarGraph(const PlanarGraph &) = delete;
+    /// \brief Assignment of a graph to another one is \e not allowed.
+    /// Use GraphCopy instead.
+    void operator=(const PlanarGraph &) = delete;
+
+    /// \brief Constructor
+
+    /// Constructor.
+    ///
+    PlanarGraph() {}
+
+    typedef Parent::OutArcIt IncEdgeIt;
+
+    /// \brief Add a new node to the graph.
+    ///
+    /// This function adds a new node to the graph. A new outer face is created
+    /// for the node.
+    /// \return The new node.
+    Node addNode() {
+      return Parent::addNode();
+    }
+
+    /// \brief Add a new edge to the graph.
+    ///
+    /// This function adds a new edge to the graph between nodes
+    /// \c u and \c v with inherent orientation from node \c u to
+    /// node \c v. \c p_u and \c p_v are the edges that directly precede the new
+    /// edge in anticlockwise order at the nodes \c u and \c v, respectively.
+    /// INVALID should be passed as \c p_u or \c p_v if and only if the
+    /// respective node is isolated.
+    ///
+    /// If \c u and \c v are in the same component, \c p_u and \c p_v must share
+    /// the same left face (when looking from \c u or \c v). This face will be
+    /// split in two. If \c u and \c v are in different components, one of the
+    /// outer faces will be deleted.
+    /// \note It can be costly to join components unless one of them is an
+    /// isolated node.
+    /// \return The new edge, or INVALID if the parameters don't meet the
+    /// preconditions.
+    Edge addEdge(Node u, Node v, Edge p_u, Edge p_v) {
+      return PlanarGraphBase::addEdge(u, v, p_u, p_v);
+    }
+
+    ///\brief Erase a node from the graph.
+    ///
+    /// This function erases the given node along with its incident arcs
+    /// from the graph.
+    ///
+    /// \note All iterators referencing the removed node or the incident
+    /// edges are invalidated, of course.
+    void erase(Node n) {
+      Parent::erase(n);
+    }
+
+    ///\brief Erase an edge from the graph.
+    ///
+    /// This function erases the given edge from the graph. The faces on the two
+    /// sides are merged into one, unless the edge holds two components
+    /// together. In the latter case, a new outer face is added and the
+    /// component is split in two.
+    ///
+    /// \note It can be costly to split a component, unless one of the resulting
+    /// components is an isolated node.
+    /// \note All iterators referencing the removed edge are invalidated,
+    /// of course.
+    void erase(Edge e) {
+      Parent::erase(e);
+    }
+    /// Node validity check
+
+    /// This function gives back \c true if the given node is valid,
+    /// i.e. it is a real node of the graph.
+    ///
+    /// \warning A removed node could become valid again if new nodes are
+    /// added to the graph.
+    bool valid(Node n) const {
+      return Parent::valid(n);
+    }
+    /// Edge validity check
+
+    /// This function gives back \c true if the given edge is valid,
+    /// i.e. it is a real edge of the graph.
+    ///
+    /// \warning A removed edge could become valid again if new edges are
+    /// added to the graph.
+    bool valid(Edge e) const {
+      return Parent::valid(e);
+    }
+    /// Arc validity check
+
+    /// This function gives back \c true if the given arc is valid,
+    /// i.e. it is a real arc of the graph.
+    ///
+    /// \warning A removed arc could become valid again if new edges are
+    /// added to the graph.
+    bool valid(Arc a) const {
+      return Parent::valid(a);
+    }
+
+    /// Face validity check
+
+    /// This function gives back \c true if the given face is valid,
+    /// i.e. it is a real face of the graph.
+    ///
+    /// \warning A removed face could become valid again if new edges are
+    /// added to the graph.
+    bool valid(Face f) const {
+      return Parent::valid(f);
+    }
+
+    /// \brief Change the first node of an edge.
+    ///
+    /// Planar graphs don't support changing the endpoints of edges.
+    void changeU(Edge e, Node n) = delete;
+    /// \brief Change the second node of an edge.
+    ///
+    /// Planar graphs don't support changing the endpoints of edges.
+    void changeV(Edge e, Node n) = delete;
+
+    /// \brief Contract two nodes.
+    ///
+    /// This function contracts the given two nodes.
+    /// Node \c b is removed, but instead of deleting
+    /// its incident edges, they are joined to node \c a.
+    /// The two nodes must have exactly one edge between them, otherwise the
+    /// function will fail. The joining edge will be deleted.
+    ///
+    /// \note All edge and arc iterators whose base node is
+    /// \c b are invalidated.
+    /// Moreover all iterators referencing node \c b or the removed
+    /// edge are also invalidated. Other iterators remain valid.
+    ///
+    ///\warning This functionality cannot be used together with the
+    ///Snapshot feature.
+    /*TODO: rewrite this function
+            void contract(Node a, Node b, bool r = true) {
+                for(IncEdgeIt e(*this, b); e!=INVALID;) {
+                    IncEdgeIt f = e; ++f;
+                    if (r && runningNode(e) == a) {
+                        erase(e);
+                    } else if (u(e) == b) {
+                        changeU(e, a);
+                    } else {
+                        changeV(e, a);
+                    }
+                    e = f;
+                }
+                erase(b);
             }
-            int &target_at(bool d) {
-                return d?target:source;
-            }
-            int source_at(bool d) const {
-                return d?source:target;
-            }
-            int &source_at(bool d) {
-                return d?source:target;
-            }
-            int prev_at(bool d) const {
-                return d?prev_s:prev_t;
-            }
-            int &prev_at(bool d) {
-                return d?prev_s:prev_t;
-            }
-            int next_at(bool d) const {
-                return d?next_s:next_t;
-            }
-            int &next_at(bool d) {
-                return d?next_s:next_t;
-            }
-            int left_face_at(bool d) const {
-                return d?left_face:right_face;
-            }
-            int &left_face_at(bool d) {
-                return d?left_face:right_face;
-            }
-            int right_face_at(bool d) const {
-                return d?right_face:left_face;
-            }
-            int &right_face_at(bool d) {
-                return d?right_face:left_face;
-            }
-        };
+    */
 
-        struct FaceT {
-            int prev, next;
-            int first_arc;
-            bool first_arc_dir;
-        };
+    ///Clear the graph.
 
-        std::vector<NodeT> nodes;
+    ///This function erases all nodes and arcs from the graph.
+    ///
+    ///\note All iterators of the graph are invalidated, of course.
+    void clear() {
+      Parent::clear();
+    }
 
-        int first_node;
+    /// Reserve memory for nodes.
 
-        int first_free_node;
+    /// Using this function, it is possible to avoid superfluous memory
+    /// allocation: if you know that the graph you want to build will
+    /// be large (e.g. it will contain millions of nodes and/or edges),
+    /// then it is worth reserving space for this amount before starting
+    /// to build the graph.
+    /// \sa reserveEdge()
+    void reserveNode(int n) {
+      nodes.reserve(n);
+    };
 
-        std::vector<ArcT> arcs;
+    /// Reserve memory for edges.
 
-        int first_free_arc;
+    /// Using this function, it is possible to avoid superfluous memory
+    /// allocation: if you know that the graph you want to build will
+    /// be large (e.g. it will contain millions of nodes and/or edges),
+    /// then it is worth reserving space for this amount before starting
+    /// to build the graph.
+    /// \sa reserveNode()
+    void reserveEdge(int m) {
+      arcs.reserve(2 * m);
+    };
 
-        std::vector<FaceT> faces;
+    /// Reserve memory for faces.
 
-        int first_face;
+    /// Using this function, it is possible to avoid superfluous memory
+    /// allocation: if you know that the graph you want to build will
+    /// be large (e.g. it will contain millions of nodes and/or edges),
+    /// then it is worth reserving space for this amount before starting
+    /// to build the graph.
+    /// \sa reserveFace()
+    void reserveFace(int n) {
+      faces.reserve(n);
+    };
 
-        int first_free_face;
+    class DualBase {
+      const Graph *_graph;
+    protected:
+      void initialize(const Graph *graph) {
+        _graph = graph;
+      }
+    public:
 
-        int component_id;
+      typedef PlanarGraph::Face Node;
+      typedef PlanarGraph::Arc Arc;
+      typedef PlanarGraph::Edge Edge;
+      typedef PlanarGraph::Node Face;
+
+      int maxNodeId() const {
+        return _graph->maxFaceId();
+      }
+      int maxArcId() const {
+        return _graph->maxArcId();
+      }
+      int maxFaceId() const {
+        return _graph->maxNodeId();
+      }
+
+      Node source(Arc e) const {
+        return _graph->leftFace(e);
+      }
+      Node target(Arc e) const {
+        return _graph->rightFace(e);
+      }
+      Face leftFace(Arc e) const {
+        return _graph->target(e);
+      }
+      Face rightFace(Arc e) const {
+        return _graph->source(e);
+      }
+      Arc direct(const Edge &edge, const Node &node) const {
+        return _graph->direct(edge, _graph->w1(edge) == node);
+      }
+
+      void first(Node &i) const {
+        _graph->first(i);
+      }
+      void next(Node &i) const {
+        _graph->next(i);
+      }
+      void first(Arc &i) const {
+        _graph->first(i);
+      }
+      void next(Arc &i) const {
+        _graph->next(i);
+      }
+      void firstCcw(Arc& i, const Node& n) const {
+        _graph->lastInF(i, n);
+      }
+      void nextCcw(Arc& i, const Node &n) const {
+        _graph->prevInF(i);
+      }
+      void firstIn(Arc& i, const Node& n) const {
+        _graph->firstInF(i, n);
+      }
+      void nextIn(Arc& i) const {
+        _graph->nextInF(i);
+      }
+      void firstCwF(Arc& i, const Face& n) const {
+        _graph->lastIn(i, n);
+      }
+      void nextCwF(Arc& i) const {
+        _graph->prevIn(i);
+      }
+      void firstOut(Arc& i, const Node& n ) const {
+        _graph->firstOutF(i, n);
+      }
+      void nextOut(Arc& i) const {
+        _graph->nextOutF(i);
+      }
+      void first(Face &i) const {
+        _graph->first(i);
+      }
+      void next(Face &i) const {
+        _graph->next(i);
+      }
+
+      static int id(Node v) {
+        return PlanarGraph::id(v);
+      }
+      static int id(Arc e) {
+        return PlanarGraph::id(e);
+      }
+      static int id(Face f) {
+        return PlanarGraph::id(f);
+      }
+      static Node nodeFromId(int id) {
+        return PlanarGraph::faceFromId(id);
+      }
+      static Arc arcFromId(int id) {
+        return PlanarGraph::arcFromId(id);
+      }
+      static Face faceFromId(int id) {
+        return PlanarGraph::nodeFromId(id);
+      }
+
+      bool valid(Node n) const {
+        return _graph->valid(n);
+      }
+      bool valid(Arc n) const {
+        return _graph->valid(n);
+      }
+      bool valid(Face n) const {
+        return _graph->valid(n);
+      }
+
+    };
+
+    typedef PlanarGraphExtender<GraphExtender<DualBase> > ExtendedDualBase;
+
+  /// Adaptor class for the dual of a planar graph.
+
+  /// This is an adaptor class for the dual of a planar graph.
+  class Dual : public ExtendedDualBase {
+    public:
+      Dual(const PlanarGraph &graph) {
+        initialize(&graph);
+      }
+
+    };
+    /// \brief Class to make a snapshot of the graph and restore
+    /// it later.
+    ///
+    /// Class to make a snapshot of the graph and restore it later.
+    ///
+    /// The newly added nodes and edges can be removed
+    /// using the restore() function.
+    ///
+    /// \note After a state is restored, you cannot restore a later state,
+    /// i.e. you cannot add the removed nodes and edges again using
+    /// another Snapshot instance.
+    ///
+    /// \warning Node and edge deletions and other modifications
+    /// (e.g. changing the end-nodes of edges or contracting nodes)
+    /// cannot be restored. These events invalidate the snapshot.
+    /// However, the edges and nodes that were added to the graph after
+    /// making the current snapshot can be removed without invalidating it.
+    class Snapshot {
+    protected:
+
+      typedef Parent::NodeNotifier NodeNotifier;
+
+    class NodeObserverProxy : public NodeNotifier::ObserverBase {
+      public:
+
+        NodeObserverProxy(Snapshot& _snapshot)
+            : snapshot(_snapshot) {}
+
+        using NodeNotifier::ObserverBase::attach;
+        using NodeNotifier::ObserverBase::detach;
+        using NodeNotifier::ObserverBase::attached;
+
+      protected:
+
+        virtual void add(const Node& node) {
+          snapshot.addNode(node);
+        }
+        virtual void add(const std::vector<Node>& nodes) {
+          for (int i = nodes.size() - 1; i >= 0; ++i) {
+            snapshot.addNode(nodes[i]);
+          }
+        }
+        virtual void erase(const Node& node) {
+          snapshot.eraseNode(node);
+        }
+        virtual void erase(const std::vector<Node>& nodes) {
+          for (int i = 0; i < int(nodes.size()); ++i) {
+            snapshot.eraseNode(nodes[i]);
+          }
+        }
+        virtual void build() {
+          Node node;
+          std::vector<Node> nodes;
+          for (notifier()->first(node); node != INVALID;
+               notifier()->next(node)) {
+            nodes.push_back(node);
+          }
+          for (int i = nodes.size() - 1; i >= 0; --i) {
+            snapshot.addNode(nodes[i]);
+          }
+        }
+        virtual void clear() {
+          Node node;
+          for (notifier()->first(node); node != INVALID;
+               notifier()->next(node)) {
+            snapshot.eraseNode(node);
+          }
+        }
+
+        Snapshot& snapshot;
+      };
+
+    class EdgeObserverProxy : public EdgeNotifier::ObserverBase {
+      public:
+
+        EdgeObserverProxy(Snapshot& _snapshot)
+            : snapshot(_snapshot) {}
+
+        using EdgeNotifier::ObserverBase::attach;
+        using EdgeNotifier::ObserverBase::detach;
+        using EdgeNotifier::ObserverBase::attached;
+
+      protected:
+
+        virtual void add(const Edge& edge) {
+          snapshot.addEdge(edge);
+        }
+        virtual void add(const std::vector<Edge>& edges) {
+          for (int i = edges.size() - 1; i >= 0; ++i) {
+            snapshot.addEdge(edges[i]);
+          }
+        }
+        virtual void erase(const Edge& edge) {
+          snapshot.eraseEdge(edge);
+        }
+        virtual void erase(const std::vector<Edge>& edges) {
+          for (int i = 0; i < int(edges.size()); ++i) {
+            snapshot.eraseEdge(edges[i]);
+          }
+        }
+        virtual void build() {
+          Edge edge;
+          std::vector<Edge> edges;
+          for (notifier()->first(edge); edge != INVALID;
+               notifier()->next(edge)) {
+            edges.push_back(edge);
+          }
+          for (int i = edges.size() - 1; i >= 0; --i) {
+            snapshot.addEdge(edges[i]);
+          }
+        }
+        virtual void clear() {
+          Edge edge;
+          for (notifier()->first(edge); edge != INVALID;
+               notifier()->next(edge)) {
+            snapshot.eraseEdge(edge);
+          }
+        }
+
+        Snapshot& snapshot;
+      };
+
+      PlanarGraph *graph;
+
+      NodeObserverProxy node_observer_proxy;
+      EdgeObserverProxy edge_observer_proxy;
+
+      std::list<Node> added_nodes;
+      std::list<Edge> added_edges;
+
+
+      void addNode(const Node& node) {
+        added_nodes.push_front(node);
+      }
+      void eraseNode(const Node& node) {
+        std::list<Node>::iterator it =
+          std::find(added_nodes.begin(), added_nodes.end(), node);
+        if (it == added_nodes.end()) {
+          clear();
+          edge_observer_proxy.detach();
+          throw NodeNotifier::ImmediateDetach();
+        } else {
+          added_nodes.erase(it);
+        }
+      }
+
+      void addEdge(const Edge& edge) {
+        added_edges.push_front(edge);
+      }
+      void eraseEdge(const Edge& edge) {
+        std::list<Edge>::iterator it =
+          std::find(added_edges.begin(), added_edges.end(), edge);
+        if (it == added_edges.end()) {
+          clear();
+          node_observer_proxy.detach();
+          throw EdgeNotifier::ImmediateDetach();
+        } else {
+          added_edges.erase(it);
+        }
+      }
+
+      void attach(PlanarGraph &_graph) {
+        graph = &_graph;
+        node_observer_proxy.attach(graph->notifier(Node()));
+        edge_observer_proxy.attach(graph->notifier(Edge()));
+      }
+
+      void detach() {
+        node_observer_proxy.detach();
+        edge_observer_proxy.detach();
+      }
+
+      bool attached() const {
+        return node_observer_proxy.attached();
+      }
+
+      void clear() {
+        added_nodes.clear();
+        added_edges.clear();
+      }
 
     public:
 
-        typedef PlanarDigraphBase Digraph;
+      /// \brief Default constructor.
+      ///
+      /// Default constructor.
+      /// You have to call save() to actually make a snapshot.
+      Snapshot()
+          : graph(0), node_observer_proxy(*this),
+          edge_observer_proxy(*this) {}
 
-        class Node {
-            friend class PlanarDigraphBase;
-            friend class PlanarDigraph;
-        protected:
+      /// \brief Constructor that immediately makes a snapshot.
+      ///
+      /// This constructor immediately makes a snapshot of the given graph.
+      Snapshot(PlanarGraph &gr)
+          : node_observer_proxy(*this),
+          edge_observer_proxy(*this) {
+        attach(gr);
+      }
 
-            int id;
-            explicit Node(int pid) { id = pid;}
+      /// \brief Make a snapshot.
+      ///
+      /// This function makes a snapshot of the given graph.
+      /// It can be called more than once. In case of a repeated
+      /// call, the previous snapshot gets lost.
+      void save(PlanarGraph &gr) {
+        if (attached()) {
+          detach();
+          clear();
+        }
+        attach(gr);
+      }
 
-        public:
-            Node() {}
-            Node (Invalid) { id = -1; }
-            bool operator==(const Node& node) const {return id == node.id;}
-            bool operator!=(const Node& node) const {return id != node.id;}
-            bool operator<(const Node& node) const {return id < node.id;}
-        };
+      /// \brief Undo the changes until the last snapshot.
+      ///
+      /// This function undos the changes until the last snapshot
+      /// created by save() or Snapshot(PlanarGraph&).
+      ///
+      /// \warning This method invalidates the snapshot, i.e. repeated
+      /// restoring is not supported unless you call save() again.
+      void restore() {
+        detach();
+        for (std::list<Edge>::iterator it = added_edges.begin();
+             it != added_edges.end(); ++it) {
+          graph->erase(*it);
+        }
+        for (std::list<Node>::iterator it = added_nodes.begin();
+             it != added_nodes.end(); ++it) {
+          graph->erase(*it);
+        }
+        clear();
+      }
 
-        class Arc {
-            friend class PlanarDigraphBase;
-            friend class PlanarDigraph;
-        protected:
+      /// \brief Returns \c true if the snapshot is valid.
+      ///
+      /// This function returns \c true if the snapshot is valid.
+      bool valid() const {
+        return attached();
+      }
+    };
+  };
 
-            int id;
-            bool dir;   //travel direction during iteration,
-                        //true means we are travelling from source to target
-            explicit Arc(int pid) { id = pid; dir = true; }
-
-        public:
-            Arc() {}
-            Arc (Invalid) { id = -1; }
-            bool operator==(const Arc& arc) const {return id == arc.id;}
-            bool operator!=(const Arc& arc) const {return id != arc.id;}
-            bool operator<(const Arc& arc) const {return id < arc.id;}
-        };
-
-        class Face {
-            friend class PlanarDigraphBase;
-            friend class PlanarDigraph;
-        protected:
-
-            int id;
-            explicit Face(int pid) { id = pid;}
-
-        public:
-            Face() {}
-            Face (Invalid) { id = -1; }
-            bool operator==(const Face& face) const {return id == face.id;}
-            bool operator!=(const Face& face) const {return id != face.id;}
-            bool operator<(const Face& face) const {return id < face.id;}
-        };
-
-
-
-        PlanarDigraphBase()
-            : nodes(), first_node(-1), first_free_node(-1),
-              arcs(), first_free_arc(-1),
-              faces(), first_face(-1), first_free_face(-1),
-              component_id(0) {
-        }
-
-
-        int maxNodeId() const { return nodes.size()-1; }
-        int maxArcId() const { return arcs.size()-1; }
-        int maxFaceId() const { return faces.size()-1; }
-
-        Node source(Arc e) const { return Node(arcs[e.id].source); }
-        Node target(Arc e) const { return Node(arcs[e.id].target); }
-        Face leftFace(Arc e) const { return Face(arcs[e.id].left_face); }
-        Face rightFace(Arc e) const { return Face(arcs[e.id].right_face); }
-
-
-        void first(Node& node) const {
-            node.id = first_node;
-        }
-
-        void next(Node& node) const {
-            node.id = nodes[node.id].next;
-        }
-
-
-        void first(Arc& arc) const {
-            int n;
-            for(n = first_node;
-                n != -1 && nodes[n].first_out == -1;
-                n = nodes[n].next) {}
-            arc.id = (n == -1) ? -1 : nodes[n].first_out;
-            arc.dir = true;
-        }
-
-        void next(Arc& arc) const {
-            if (arcs[arc.id].next_out != -1) {
-                arc.id = arcs[arc.id].next_out;
-            } else {
-                int n;
-                for(n = nodes[arcs[arc.id].source].next;
-                    n != -1 && nodes[n].first_out == -1;
-                    n = nodes[n].next) {}
-                arc.id = (n == -1) ? -1 : nodes[n].first_out;
-            }
-        }
-
-        void firstOut(Arc &e, const Node& v) const {
-            e.id = nodes[v.id].first_out;
-            e.dir = true;
-        }
-        void nextOut(Arc &e) const {
-            e.id=arcs[e.id].next_out;
-        }
-
-        void firstIn(Arc &e, const Node& v) const {
-            e.id = nodes[v.id].first_in;
-            e.dir = false;
-        }
-        void nextIn(Arc &e) const {
-            e.id=arcs[e.id].next_in;
-        }
-
-
-        void firstCcw(Arc &e, const Node &v) const {
-            e.id = nodes[v.id].first_ccw;
-            if (e.id > -1) e.dir = v.id == arcs[e.id].source;
-        }
-        void nextCcw(Arc &e) const {
-            int n = arcs[e.id].source_at(e.dir);
-            e.id = arcs[e.id].next_at(e.dir);
-            if (e.id > -1) e.dir = arcs[e.id].source == n;
-        }
-
-        void firstCw(Arc &e, const Node &v) const {
-            e.id = nodes[v.id].last_ccw;
-            if (e.id > -1) e.dir = v.id == arcs[e.id].source;
-        }
-
-        void nextCw(Arc &e) const {
-            int n = arcs[e.id].source_at(e.dir);
-            e.id = arcs[e.id].prev_at(e.dir);
-            if (e.id > -1) e.dir = arcs[e.id].source == n;
-        }
-
-        void turnLeft(Arc &e) const {
-            int n = arcs[e.id].source_at(e.dir);
-            e.id = arcs[e.id].next_at(e.dir);
-            if (e.id == -1) e.id = nodes[n].first_ccw;
-            if (e.id > -1) e.dir = arcs[e.id].source == n;
-        }
-        void turnRight(Arc &e) const {
-            int n = arcs[e.id].source_at(e.dir);
-            e.id = arcs[e.id].prev_at(e.dir);
-            if (e.id == -1) e.id = nodes[n].last_ccw;
-            if (e.id > -1) e.dir = arcs[e.id].source == n;
-        }
-
-        void first(Face& face) const {
-            face.id = first_face;
-        }
-
-        void next(Face& face) const {
-            face.id = faces[face.id].next;
-        }
-
-        void firstCcwF(Arc &e, const Face &f) const {
-            e.id = faces[f.id].first_arc;
-            e.dir = faces[f.id].first_arc_dir;
-        }
-        void nextCcwF(Arc &e) const {
-            int f = arcs[e.id].left_face_at(e.dir);
-            turnLeftF(e);
-            if (e.id == faces[f].first_arc && e.dir == faces[f].first_arc_dir) e.id = -1;
-        }
-
-        void firstCwF(Arc &e, const Face &f) const {
-            e.id = faces[f.id].first_arc;
-            e.dir = faces[f.id].first_arc_dir;
-        }
-        void nextCwF(Arc &e) const {
-            int f = arcs[e.id].left_face_at(e.dir);
-            turnRightF(e);
-            if (e.id == faces[f].first_arc && e.dir == faces[f].first_arc_dir) e = INVALID;
-        }
-
-        void turnLeftF(Arc &e) const {
-            e.dir = !e.dir;
-            turnRight(e);
-        }
-
-        void turnRightF(Arc &e) const {
-            turnLeft(e);
-            e.dir = !e.dir;
-        }
-
-        void firstInF(Arc &e, const Face &f) const {
-            e.id = faces[f.id].first_arc;
-            e.dir = faces[f.id].first_arc_dir;
-            while (e.id > -1 && arcs[e.id].right_face != f.id) {
-                turnLeft(e);
-            }
-        }
-        void nextInF(Arc &e) const {
-            int f = arcs[e.id].right_face_at(e.dir);
-            do {
-                turnLeft(e);
-                if (e.id == faces[f].first_arc && e.dir == faces[f].first_arc_dir) e.id = -1;
-            } while (e.id > -1 && arcs[e.id].right_face != f);
-        }
-
-        void firstOutF(Arc &e, const Face &f) const {
-            e.id = faces[f.id].first_arc;
-            e.dir = faces[f.id].first_arc_dir;
-            while (e.id > -1 && arcs[e.id].left_face != f.id) {
-                turnRight(e);
-            }
-        }
-        void nextOutF(Arc &e) const {
-            int f = arcs[e.id].left_face_at(e.dir);
-            do {
-                turnLeft(e);
-                if (e.id == faces[f].first_arc && e.dir == faces[f].first_arc_dir) e.id = -1;
-            } while (e.id > -1 && arcs[e.id].left_face != f);
-        }
-
-
-        static int id(Node v) { return v.id; }
-        static int id(Arc e) { return e.id; }
-        static int id(Face f) { return f.id; }
-
-        static Node nodeFromId(int id) { return Node(id);}
-        static Arc arcFromId(int id) { return Arc(id);}
-        static Face faceFromId(int id) { return Face(id);}
-
-        bool valid(Node n) const {
-            return n.id >= 0 && n.id < static_cast<int>(nodes.size()) &&
-                    nodes[n.id].prev != -2;
-        }
-
-        bool valid(Arc a) const {
-            return a.id >= 0 && a.id < static_cast<int>(arcs.size()) &&
-                    arcs[a.id].prev_in != -2;
-        }
-
-        bool valid(Face f) const {
-            return f.id >= 0 && f.id < static_cast<int>(faces.size()) &&
-                    faces[f.id].prev != -2;
-        }
-
-        Node addNode() {
-            int n;
-
-            if(first_free_node==-1) {
-                n = nodes.size();
-                nodes.push_back(NodeT());
-            } else {
-                n = first_free_node;
-                first_free_node = nodes[n].next;
-            }
-
-            nodes[n].next = first_node;
-            nodes[n].component = component_id++;
-            if(first_node != -1) nodes[first_node].prev = n;
-            first_node = n;
-            nodes[n].prev = -1;
-
-            nodes[n].first_in = nodes[n].first_out = -1;
-            nodes[n].first_ccw = -1;
-            nodes[n].last_ccw = -1;
-
-            return Node(n);
-        }
-
-        Arc addArc(Node u, Node v, Arc p_u, Arc p_v) {
-            bool ud = p_u.id > -1 && u.id == arcs[p_u.id].source;
-            bool vd = p_v.id > -1 && v.id == arcs[p_v.id].source;
-
-            if (p_u.id > -1 && p_v.id > -1 && arcs[p_u.id].left_face_at(ud) != arcs[p_v.id].left_face_at(vd)
-                && nodes[u.id].component == nodes[v.id].component) return INVALID;
-            int n = addBlankArc();
-
-            arcs[n].source = u.id;
-            arcs[n].target = v.id;
-
-            arcs[n].next_out = nodes[u.id].first_out;
-            if(nodes[u.id].first_out != -1) {
-                arcs[nodes[u.id].first_out].prev_out = n;
-            }
-
-            arcs[n].next_in = nodes[v.id].first_in;
-            if(nodes[v.id].first_in != -1) {
-                arcs[nodes[v.id].first_in].prev_in = n;
-            }
-
-            arcs[n].prev_in = arcs[n].prev_out = -1;
-
-            nodes[u.id].first_out = nodes[v.id].first_in = n;
-
-            arcs[n].prev_s = p_u.id;
-            if (p_u.id == -1) {
-                arcs[n].next_s = nodes[u.id].first_ccw;
-                nodes[u.id].first_ccw = n;
-            } else {
-                arcs[n].next_s = arcs[p_u.id].next_at(ud);
-                arcs[p_u.id].next_at(ud) = n;
-            }
-            if (arcs[n].next_s != -1) {
-                bool pvd = u.id == arcs[arcs[n].next_s].source;
-                arcs[n].prev_s = arcs[arcs[n].next_s].prev_at(pvd);
-                arcs[arcs[n].next_s].prev_at(pvd) = n;
-            } else {
-                nodes[u.id].last_ccw = n;
-            }
-
-            arcs[n].prev_t = p_v.id;
-            if (p_v.id == -1) {
-                arcs[n].next_t = nodes[v.id].first_ccw;
-                nodes[v.id].first_ccw = n;
-            } else {
-                arcs[n].next_t = arcs[p_v.id].next_at(vd);
-                arcs[p_v.id].next_at(vd) = n;
-            }
-            if (arcs[n].next_t != -1) {
-                bool nvd = v.id == arcs[arcs[n].next_t].source;
-                arcs[n].prev_t = arcs[arcs[n].next_t].prev_at(nvd);
-                arcs[arcs[n].next_t].prev_at(nvd) = n;
-            } else {
-                nodes[v.id].last_ccw = n;
-            }
-
-            //Add the extra face, if needed
-            if (p_u.id > -1 && p_v.id > -1) {
-                int oldf = arcs[p_u.id].left_face_at(ud);
-                int oldfb = arcs[p_v.id].left_face_at(vd);
-                arcs[n].left_face = arcs[n].right_face = oldf;
-                Face f = addFace();
-                faces[f.id].first_arc = n;
-                faces[f.id].first_arc_dir = true;
-                faces[oldf].first_arc = n;
-                faces[oldf].first_arc_dir = false;
-                Arc arc(n);
-                wall_paint(arc,f.id,arc);
-                if (nodes[v.id].component != nodes[u.id].component) {
-                    erase(Face(oldf));
-                    erase(Face(oldfb));
-                    int ca = nodes[u.id].component;
-                    int cb = nodes[v.id].component;
-                    int k = first_node;
-                    while (k != -1) {
-                        if (nodes[k].component == cb)
-                            nodes[k].component = ca;
-                        k = nodes[k].next;
-                    }
-                }
-            } else if (p_u.id > -1) {
-                arcs[n].left_face = arcs[n].right_face = arcs[p_u.id].left_face_at(ud);
-                faces[arcs[n].left_face].first_arc = n;
-                faces[arcs[n].left_face].first_arc_dir = true;
-                nodes[v.id].component = nodes[u.id].component;
-            } else if (p_v.id > -1) {
-                arcs[n].left_face = arcs[n].right_face = arcs[p_v.id].left_face_at(vd);
-                faces[arcs[n].left_face].first_arc = n;
-                faces[arcs[n].left_face].first_arc_dir = true;
-                nodes[u.id].component = nodes[v.id].component;
-            } else {    //both prevs are INVALID
-                Face f = addFace();
-                arcs[n].left_face = arcs[n].right_face = f.id;
-                faces[f.id].first_arc = n;
-                faces[f.id].first_arc_dir = true;
-                nodes[v.id].component = nodes[u.id].component;
-            }
-
-            return Arc(n);
-        }
-
-    public:
-
-        void erase(const Node& node) {  //breaks graph unless node is isolated
-            int n = node.id;
-
-            if(nodes[n].next != -1) {
-                nodes[nodes[n].next].prev = nodes[n].prev;
-            }
-
-            if(nodes[n].prev != -1) {
-                nodes[nodes[n].prev].next = nodes[n].next;
-            } else {
-                first_node = nodes[n].next;
-            }
-
-            nodes[n].next = first_free_node;
-            first_free_node = n;
-            nodes[n].prev = -2;
-
-        }
-
-        void erase(const Arc& arc) {
-            int n = arc.id;
-
-            //"retreat" the incident faces' first arcs
-            int fl = arcs[n].left_face;
-            if (faces[fl].first_arc == n) {
-                Arc e(faces[fl].first_arc);
-                e.dir = faces[fl].first_arc;
-                turnLeftF(e);
-                if (e.id == n) turnLeftF(e);
-                faces[fl].first_arc = e.id;
-                faces[fl].first_arc_dir = e.dir;
-            }
-
-            int fr = arcs[n].right_face;
-
-            bool comp_split = false;
-            if (fr != fl) {
-                Arc arc(faces[fr].first_arc);
-                arc.dir = faces[fr].first_arc_dir;
-                wall_paint(arc,fl,arc);
-                erase(Face(fr));
-            } else if ((arcs[n].next_s > -1 || arcs[n].prev_s > -1) && (arcs[n].next_t > -1 || arcs[n].prev_t > -1)) {
-                comp_split = true;
-                Arc arc(n);
-                Arc ed = arc;
-                ed.dir = false;
-                turnRightF(arc);
-                Face f = addFace();
-                wall_paint(arc,f.id,ed);
-                faces[f.id].first_arc = arc.id;
-                faces[f.id].first_arc_dir = arc.dir;
-            }
-
-            unlink_arc_inout(n);
-
-            if (arcs[n].next_s != -1) {
-                arcs[arcs[n].next_s].prev_at(arcs[n].source == arcs[arcs[n].next_s].source) = arcs[n].prev_s;
-            } else {
-                nodes[arcs[n].source].last_ccw = arcs[n].prev_s;
-            }
-            if (arcs[n].prev_s != -1) {
-                arcs[arcs[n].prev_s].next_at(arcs[n].source == arcs[arcs[n].prev_s].source) = arcs[n].next_s;
-            } else {
-                nodes[arcs[n].source].first_ccw = arcs[n].next_s;
-            }
-
-            if (arcs[n].next_t != -1) {
-                arcs[arcs[n].next_t].prev_at(arcs[n].target == arcs[arcs[n].next_t].source) = arcs[n].prev_t;
-            } else {
-                nodes[arcs[n].target].last_ccw = arcs[n].prev_t;
-            }
-            if (arcs[n].prev_t != -1) {
-                arcs[arcs[n].prev_t].next_at(arcs[n].target == arcs[arcs[n].prev_t].source) = arcs[n].next_t;
-            } else {
-                nodes[arcs[n].target].first_ccw = arcs[n].next_t;
-            }
-
-            unlink_arc_id(n);
-
-            if (comp_split) component_relabel(Node(arcs[n].target), component_id++);
-        }
-
-
-        void clear() {
-            arcs.clear();
-            nodes.clear();
-            first_node = first_free_node = first_free_arc = first_face = first_free_face = -1;
-            faces.clear();
-        }
-
-        void reverseArc(Arc a) {
-            int n = a.id;
-            int tmp;
-            if (faces[arcs[n].left_face].first_arc == n)
-                faces[arcs[n].left_face].first_arc_dir = !faces[arcs[n].left_face].first_arc_dir;
-            if (faces[arcs[n].right_face].first_arc == n)
-                faces[arcs[n].right_face].first_arc_dir = !faces[arcs[n].right_face].first_arc_dir;
-
-            if (arcs[n].prev_out > -1)
-                arcs[arcs[n].prev_out].next_out = arcs[n].next_out;
-            else
-                nodes[arcs[n].source].first_out = arcs[n].next_out;
-            if (arcs[n].next_out > -1)
-                arcs[arcs[n].next_out].prev_out = arcs[n].prev_out;
-
-            if (arcs[n].prev_in > -1)
-                arcs[arcs[n].prev_in].next_in = arcs[n].next_in;
-            else
-                nodes[arcs[n].source].first_in = arcs[n].next_in;
-            if (arcs[n].next_in > -1)
-                arcs[arcs[n].next_in].prev_in = arcs[n].prev_in;
-
-            tmp = arcs[n].left_face; arcs[n].left_face = arcs[n].right_face; arcs[n].right_face = tmp;
-            tmp = arcs[n].prev_s; arcs[n].prev_s = arcs[n].prev_t; arcs[n].prev_t = tmp;
-            tmp = arcs[n].next_s; arcs[n].next_s = arcs[n].next_t; arcs[n].next_t = tmp;
-            tmp = arcs[n].source; arcs[n].source = arcs[n].target; arcs[n].target = tmp;
-
-            arcs[n].prev_out = -1;
-            arcs[n].next_out = nodes[arcs[n].source].first_out;
-            nodes[arcs[n].source].first_out = n;
-
-            arcs[n].prev_in = -1;
-            arcs[n].next_in = nodes[arcs[n].target].first_in;
-            nodes[arcs[n].target].first_in = n;
-        }
-
-        Node split(Arc a) {
-            Node v = addNode();
-            nodes[v.id].component = nodes[arcs[a.id].target].component;
-            int b = addBlankArc();
-            arcs[b] = arcs[a.id];
-            arcs[b].source = v.id;
-            arcs[b].next_out = -1;
-            arcs[b].prev_out = -1;
-            arcs[b].next_s = -1;
-            arcs[b].prev_s = a.id;
-
-            if (arcs[b].next_in != -1) {
-                arcs[arcs[b].next_in].prev_in = b;
-            }
-            if (arcs[b].prev_in != -1) {
-                arcs[arcs[b].prev_in].next_in = b;
-            } else {
-                nodes[arcs[b].target].first_in = b;
-            }
-
-            if (arcs[b].next_t != -1) {
-                arcs[arcs[b].next_t].prev_at(arcs[b].target) = b;
-            } else {
-                nodes[arcs[b].target].last_ccw = b;
-            }
-            if (arcs[b].prev_t != -1) {
-                arcs[arcs[b].prev_t].next_at(arcs[b].target) = b;
-            } else {
-                nodes[arcs[b].target].first_ccw = b;
-            }
-
-            arcs[a.id].target = v.id;
-            arcs[a.id].next_in = -1;
-            arcs[a.id].prev_in = -1;
-            arcs[a.id].next_t = b;
-            arcs[a.id].prev_t = -1;
-
-            nodes[v.id].first_in = a.id;
-            nodes[v.id].first_out = b;
-            nodes[v.id].first_ccw = a.id;
-            nodes[v.id].last_ccw = b;
-
-            return v;
-        }
-
-        bool contract(Node n1, Node n2) {
-            Arc fs;
-            firstCcw(fs,n1);
-            while (fs != INVALID && arcs[fs.id].target_at(fs.dir) != n2.id)
-                nextCcw(fs);
-            if (fs == INVALID) return false;
-            Arc bs = fs;
-            turnRight(bs);
-            while (bs != fs && arcs[bs.id].target_at(bs.dir) == n2.id)
-                turnRight(bs);
-            if (bs == fs) {
-                while (nodes[n1.id].first_ccw != -1)
-                    erase(Arc(nodes[n1.id].first_ccw));
-                erase(n1);
-                return true;
-            }
-            turnLeft(bs);
-            while (bs != fs) {
-                Arc b2s = bs;
-                turnLeft(b2s);
-                erase(bs);
-                bs = b2s;
-            }
-
-            int fl = arcs[fs.id].left_face;
-            if (faces[fl].first_arc == fs.id) {
-                Arc e(faces[fl].first_arc);
-                e.dir = faces[fl].first_arc;
-                turnLeftF(e);
-                if (e.id == fs.id) turnLeftF(e);
-                faces[fl].first_arc = e.id;
-                faces[fl].first_arc_dir = e.dir;
-            }
-            int fr = arcs[fs.id].right_face;
-            if (faces[fr].first_arc == fs.id) {
-                Arc e(faces[fr].first_arc);
-                e.dir = faces[fr].first_arc;
-                turnLeftF(e);
-                if (e.id == fs.id) turnLeftF(e);
-                faces[fr].first_arc = e.id;
-                faces[fr].first_arc_dir = e.dir;
-            }
-
-            Arc a;
-            firstIn(a,n2);
-            while (arcs[a.id].next_in != -1) {
-                if (a.id != fs.id) arcs[a.id].target = n1.id;
-                nextIn(a);
-            }
-            if (a.id != fs.id) arcs[a.id].target = n1.id;
-            arcs[a.id].next_in = nodes[n1.id].first_in;
-            nodes[n1.id].first_in = nodes[n2.id].first_in;
-
-            firstOut(a,n2);
-            while (arcs[a.id].next_out != -1) {
-                if (a.id != fs.id) arcs[a.id].source = n1.id;
-                nextOut(a);
-            }
-            if (a.id != fs.id) arcs[a.id].source = n1.id;
-            arcs[a.id].next_out = nodes[n1.id].first_out;
-            nodes[n1.id].first_out = nodes[n2.id].first_out;
-            unlink_arc_inout(fs.id);
-
-            int la = nodes[n2.id].last_ccw;
-            arcs[la].next_at(arcs[la].source == n2.id) = nodes[n2.id].first_ccw;
-            la = nodes[n2.id].first_ccw;
-            arcs[la].prev_at(arcs[la].source == n2.id) = nodes[n2.id].last_ccw;
-            nodes[n2.id].first_ccw = -1;
-            nodes[n2.id].last_ccw = -1;
-
-            Arc ta = fs;
-            turnRight(ta);
-            if (nodes[n1.id].last_ccw == fs.id) nodes[n1.id].last_ccw = ta.id;
-            arcs[ta.id].next_at(ta.dir) = arcs[fs.id].next_at(!fs.dir);
-            ta = fs;
-            turnLeft(ta);
-            if (nodes[n1.id].first_ccw == fs.id) nodes[n1.id].first_ccw = ta.id;
-            arcs[ta.id].prev_at(ta.dir) = arcs[fs.id].prev_at(!fs.dir);
-            ta = fs;
-            ta.dir = !ta.dir;
-            turnLeft(ta);
-            arcs[ta.id].prev_at(ta.dir) = arcs[fs.id].next_at(fs.dir);
-            ta = fs;
-            ta.dir = !ta.dir;
-            turnRight(ta);
-            arcs[ta.id].next_at(ta.dir) = arcs[fs.id].prev_at(fs.dir);
-
-            unlink_arc_id(fs.id);
-            erase(n2);
-            return true;
-        }
-
-
-    protected:
-
-        //traverses a face clockwise, relabelling the arcs with the face id
-        void wall_paint(Arc arc, int f_id, Arc ed) {
-            do {
-                arcs[arc.id].left_face_at(arc.dir) = f_id;
-                turnRightF(arc);
-            } while (arc.id != ed.id || arc.dir != ed.dir);
-        }
-
-        void component_relabel(Node node, int comp_id) {
-            std::vector<int> ns(nodes.size());
-            std::list<int> q;
-            q.push_back(node.id);
-            ns[node.id] = 1;
-            while (!q.empty()) {
-                int n = q.front();
-                ns[n] = 2;
-                nodes[n].component = comp_id;
-                q.pop_front();
-                Arc arc;
-                firstCcw(arc,Node(n));
-                while (arc.id > -1) {
-                    int m = arcs[arc.id].target_at(arcs[arc.id].source == n);
-                    if (ns[m] == 0) {
-                        ns[m] = 1;
-                        q.push_back(m);
-                    }
-                    nextCcw(arc);
-                }
-            }
-        }
-
-        void unlink_arc_inout(int n) {
-            if(arcs[n].next_in!=-1) {
-                arcs[arcs[n].next_in].prev_in = arcs[n].prev_in;
-            }
-
-            if(arcs[n].prev_in!=-1) {
-                arcs[arcs[n].prev_in].next_in = arcs[n].next_in;
-            } else {
-                nodes[arcs[n].target].first_in = arcs[n].next_in;
-            }
-
-
-            if(arcs[n].next_out!=-1) {
-                arcs[arcs[n].next_out].prev_out = arcs[n].prev_out;
-            }
-
-            if(arcs[n].prev_out!=-1) {
-                arcs[arcs[n].prev_out].next_out = arcs[n].next_out;
-            } else {
-                nodes[arcs[n].source].first_out = arcs[n].next_out;
-            }
-
-        }
-
-        void unlink_arc_id(int n) {
-            arcs[n].next_in = first_free_arc;
-            first_free_arc = n;
-            arcs[n].prev_in = -2;
-        }
-
-        int addBlankArc() {
-            int n;
-            if (first_free_arc == -1) {
-                n = arcs.size();
-                arcs.push_back(ArcT());
-            } else {
-                n = first_free_arc;
-                first_free_arc = arcs[n].next_in;
-            }
-            return n;
-
-        }
-
-        void erase(const Face& face) {
-            int n = face.id;
-
-            if(faces[n].next != -1) {
-                faces[faces[n].next].prev = faces[n].prev;
-            }
-
-            if(faces[n].prev != -1) {
-                faces[faces[n].prev].next = faces[n].next;
-            } else {
-                first_face = faces[n].next;
-            }
-
-            faces[n].next = first_free_face;
-            first_free_face = n;
-            faces[n].prev = -2;
-
-        }
-
-
-        Face addFace() {
-            int n;
-
-            if(first_free_face==-1) {
-                n = faces.size();
-                faces.push_back(FaceT());
-            } else {
-                n = first_free_face;
-                first_free_face = faces[n].next;
-            }
-
-            faces[n].next = first_face;
-            if(first_face != -1) faces[first_face].prev = n;
-            first_face = n;
-            faces[n].prev = -1;
-
-            return Face(n);
-        }
-
-#ifdef REMOVE_BEFORE_RELEASE
-    public:
-        void print() {
-            std::cout << "Nodes: " << std::endl;
-            for (int i=0; i<nodes.size(); ++i)
-                std::cout << i << ": fi=" << nodes[i].first_in
-                          << " fo=" << nodes[i].first_out
-                          << " fc=" << nodes[i].first_ccw
-                          << " lc=" << nodes[i].last_ccw
-                          << " pr=" << nodes[i].prev
-                          << " nx=" << nodes[i].next
-                          << " co=" << nodes[i].component
-                          <<std::endl;
-            std::cout << "Arcs: " << std::endl;
-            for (int i=0; i<arcs.size(); ++i) {
-                if (arcs[i].prev_in > -2) {
-                    std::cout << i << ": sc=" << arcs[i].source
-                          << " tg=" << arcs[i].target
-                          << " pi=" << arcs[i].prev_in
-                          << " ni=" << arcs[i].next_in
-                          << " po=" << arcs[i].prev_out
-                          << " no=" << arcs[i].next_out
-                          << " ps=" << arcs[i].prev_s
-                          << " ns=" << arcs[i].next_s
-                          << " pt=" << arcs[i].prev_t
-                          << " nt=" << arcs[i].next_t
-                          << " lf=" << arcs[i].left_face
-                          << " rf=" << arcs[i].right_face
-                          <<std::endl;
-                } else std::cout << i << ": (deleted)" << std::endl;
-            }
-            std::cout << "Faces: " << std::endl;
-            for (int i=0; i<faces.size(); ++i)
-                std::cout << i
-                          << " pr=" << faces[i].prev
-                          << " nx=" << faces[i].next
-                          << " fa=" << faces[i].first_arc
-                          << " fd=" << faces[i].first_arc_dir
-                          << (faces[i].prev > -2? "" : " (deleted)")
-                          <<std::endl;
-        }
-#endif
-
-    };
-
-    template<typename Base>
-    class PlanarDigraphExtender : public Base{
-
-        typedef Base Parent;
-
-            public:
-        typedef PlanarDigraphExtender Digraph;
-
-        PlanarDigraphExtender() {}
-
-        typedef typename Parent::Node Node;
-        typedef typename Parent::Arc Arc;
-        typedef typename Parent::Face Face;
-
-        class FaceIt : public Face {
-            const Digraph* _digraph;
-        public:
-
-            FaceIt() {}
-
-            FaceIt(Invalid i) : Face(i) { }
-
-            explicit FaceIt(const Digraph& digraph) : _digraph(&digraph) {
-                _digraph->first(static_cast<Face&>(*this));
-            }
-
-            FaceIt(const Digraph& digraph, const Face& face)
-                : Face(face), _digraph(&digraph) {}
-
-            FaceIt& operator++() {
-                _digraph->next(*this);
-                return *this;
-            }
-
-        };
-
-
-        class CwPerimeterArcIt : public Arc {
-          const Digraph* _digraph;
-          Face _face;
-        public:
-
-          CwPerimeterArcIt() { }
-
-          CwPerimeterArcIt(Invalid i) : Arc(i) { }
-
-          CwPerimeterArcIt(const Digraph& digraph, const Face& face)
-            : _digraph(&digraph), _face(face) {
-            _digraph->firstCwF(*this, _face);
-          }
-
-          CwPerimeterArcIt(const Digraph& digraph, const Arc& arc)
-            : Arc(arc), _digraph(&digraph) {}
-
-          CwPerimeterArcIt& operator++() {
-            _digraph->nextCwF(*this);
-            return *this;
-          }
-
-        };
-
-        class CcwArcIt : public Arc {
-          const Digraph* _digraph;
-          const Node _node;
-        public:
-
-          CcwArcIt() { }
-
-          CcwArcIt(Invalid i) : Arc(i) { }
-
-          CcwArcIt(const Digraph& digraph, const Node& node)
-            : _digraph(&digraph), _node(node) {
-            _digraph->firstCcw(*this, node);
-          }
-
-          CcwArcIt(const Digraph& digraph, const Arc& arc)
-            : Arc(arc), _digraph(&digraph) {}
-
-          CcwArcIt& operator++() {
-            _digraph->nextCcw(*this);
-            return *this;
-          }
-
-        };
-
-    };
-
-    typedef PlanarDigraphExtender<DigraphExtender<PlanarDigraphBase> > ExtendedPlanarDigraphBase;
-
-    /// \addtogroup graphs
-    /// @{
-
-    ///A general directed graph structure.
-
-    ///\ref ListDigraph is a versatile and fast directed graph
-    ///implementation based on linked lists that are stored in
-    ///\c std::vector structures.
-    ///
-    ///This type fully conforms to the \ref concepts::Digraph "Digraph concept"
-    ///and it also provides several useful additional functionalities.
-    ///Most of its member functions and nested classes are documented
-    ///only in the concept class.
-    ///
-    ///This class provides only linear time counting for nodes and arcs.
-    ///
-    ///\sa concepts::Digraph
-    ///\sa ListGraph
-    class PlanarDigraph : public ExtendedPlanarDigraphBase {
-        typedef ExtendedPlanarDigraphBase Parent;
-
-    private:
-        /// Digraphs are \e not copy constructible. Use DigraphCopy instead.
-        PlanarDigraph(const PlanarDigraph &) :ExtendedPlanarDigraphBase() {};
-        /// \brief Assignment of a digraph to another one is \e not allowed.
-        /// Use DigraphCopy instead.
-        void operator=(const PlanarDigraph &) {}
-    public:
-
-        /// Constructor
-
-        /// Constructor.
-        ///
-        PlanarDigraph() {}
-
-        ///Add a new node to the digraph.
-
-        ///This function adds a new node to the digraph.
-        ///\return The new node.
-        Node addNode() { return Parent::addNode(); }
-
-        ///Add a new arc to the digraph.
-
-        ///This function adds a new arc to the digraph with source node \c s
-        ///and target node \c t.
-        ///\return The new arc.
-        Arc addArc(Node s, Node t, Arc p_u, Arc p_v) {
-            Arc arc = PlanarDigraphBase::addArc(s, t, p_u, p_v);
-            if (arc != INVALID)
-                notifier(Arc()).add(arc);
-            return arc;
-        }
-
-        ///\brief Erase a node from the digraph.
-        ///
-        ///This function erases the given node along with its outgoing and
-        ///incoming arcs from the digraph.
-        ///
-        ///\note All iterators referencing the removed node or the connected
-        ///arcs are invalidated, of course.
-        void erase(Node n) { Parent::erase(n); }
-
-        ///\brief Erase an arc from the digraph.
-        ///
-        ///This function erases the given arc from the digraph.
-        ///
-        ///\note All iterators referencing the removed arc are invalidated,
-        ///of course.
-        void erase(Arc a) { Parent::erase(a); }
-
-        /// Node validity check
-
-        /// This function gives back \c true if the given node is valid,
-        /// i.e. it is a real node of the digraph.
-        ///
-        /// \warning A removed node could become valid again if new nodes are
-        /// added to the digraph.
-        bool valid(Node n) const { return Parent::valid(n); }
-
-        /// Arc validity check
-
-        /// This function gives back \c true if the given arc is valid,
-        /// i.e. it is a real arc of the digraph.
-        ///
-        /// \warning A removed arc could become valid again if new arcs are
-        /// added to the digraph.
-        bool valid(Arc a) const { return Parent::valid(a); }
-
-        /// Change the target node of an arc
-
-        /// Planar graphs don't support changing endpoints of arcs.
-        void changeTarget(Arc a, Node n) = delete;
-        /// Change the source node of an arc
-
-        ///Planar graphs don't support changing endpoints of arcs.
-        void changeSource(Arc a, Node n) = delete;
-
-        /// Reverse the direction of an arc.
-
-        /// This function reverses the direction of the given arc.
-        ///\note \c ArcIt, \c OutArcIt and \c InArcIt iterators referencing
-        ///the changed arc are invalidated.
-        ///
-        ///\warning This functionality cannot be used together with the Snapshot
-        ///feature.
-        void reverseArc(Arc a) {
-            Parent::reverseArc(a);
-        }
-
-        ///Contract two nodes.
-
-        ///This function contracts the given two nodes.
-        ///Node \c v is removed, but instead of deleting its
-        ///incident arcs, they are joined to node \c u.
-        ///If the last parameter \c r is \c true (this is the default value),
-        ///then the newly created loops are removed.
-        ///
-        ///\note The moved arcs are joined to node \c u using changeSource()
-        ///or changeTarget(), thus \c ArcIt and \c OutArcIt iterators are
-        ///invalidated for the outgoing arcs of node \c v and \c InArcIt
-        ///iterators are invalidated for the incomming arcs of \c v.
-        ///Moreover all iterators referencing node \c v or the removed
-        ///loops are also invalidated. Other iterators remain valid.
-        ///
-        ///\warning This functionality cannot be used together with the Snapshot
-        ///feature.
-        bool contract(Node u, Node v)
-        {
-            return Parent::contract(u,v);
-        }
-
-        ///Split a node.
-
-        ///Planar digraphs don't support splitting of nodes.
-        Node split(Node n, bool connect = true) = delete;
-/*        Node split(Node n, bool connect = true) {
-            Node b = addNode();
-            nodes[b.id].first_out=nodes[n.id].first_out;
-            nodes[n.id].first_out=-1;
-            for(int i=nodes[b.id].first_out; i!=-1; i=arcs[i].next_out) {
-                arcs[i].source=b.id;
-            }
-            if (connect) addArc(n,b);
-            return b;
-        }*/
-
-        ///Split an arc.
-
-        ///This function splits the given arc. First, a new node \c v is
-        ///added to the digraph, then the target node of the original arc
-        ///is set to \c v. Finally, an arc from \c v to the original target
-        ///is added.
-        ///\return The newly created node.
-        ///
-        ///\note \c InArcIt iterators referencing the original arc are
-        ///invalidated. Other iterators remain valid.
-        ///
-        ///\warning This functionality cannot be used together with the
-        ///Snapshot feature.
-        Node split(Arc a) {
-            return Parent::split(a);
-        }
-
-        ///Clear the digraph.
-
-        ///This function erases all nodes and arcs from the digraph.
-        ///
-        ///\note All iterators of the digraph are invalidated, of course.
-        void clear() {
-            Parent::clear();
-        }
-
-        /// Reserve memory for nodes.
-
-        /// Using this function, it is possible to avoid superfluous memory
-        /// allocation: if you know that the digraph you want to build will
-        /// be large (e.g. it will contain millions of nodes and/or arcs),
-        /// then it is worth reserving space for this amount before starting
-        /// to build the digraph.
-        /// \sa reserveArc()
-        void reserveNode(int n) { nodes.reserve(n); };
-
-        /// Reserve memory for arcs.
-
-        /// Using this function, it is possible to avoid superfluous memory
-        /// allocation: if you know that the digraph you want to build will
-        /// be large (e.g. it will contain millions of nodes and/or arcs),
-        /// then it is worth reserving space for this amount before starting
-        /// to build the digraph.
-        /// \sa reserveNode()
-        void reserveArc(int m) { arcs.reserve(m); };
-
-        class DualBase {
-            const Digraph *_digraph;
-        protected:
-            void initialize(const Digraph *digraph) { _digraph = digraph; }
-        public:
-
-            typedef PlanarDigraph::Face Node;
-            typedef PlanarDigraph::Arc Arc;
-            typedef PlanarDigraph::Node Face;
-
-            int maxNodeId() const { return _digraph->maxFaceId(); }
-            int maxArcId() const { return _digraph->maxArcId(); }
-            int maxFaceId() const { return _digraph->maxNodeId(); }
-
-            Node source(Arc e) const { return _digraph->leftFace(e); }
-            Node target(Arc e) const { return _digraph->rightFace(e); }
-            Face leftFace(Arc e) const { return _digraph->target(e); }
-            Face rightFace(Arc e) const { return _digraph->source(e); }
-
-            void first(Node &i) const { _digraph->first(i); }
-            void next(Node &i) const { _digraph->next(i); }
-            void first(Arc &i) const { _digraph->first(i); }
-            void next(Arc &i) const { _digraph->next(i); }
-            void firstIn(Arc& i, const Node& n) const { _digraph->firstInF(i, n); }
-            void nextIn(Arc& i) const { _digraph->nextInF(i); }
-            void firstOut(Arc& i, const Node& n ) const { _digraph->firstOutF(i, n); }
-            void nextOut(Arc& i) const { _digraph->nextOutF(i); }
-            void firstCcw(Arc &e, const Node &v) const { _digraph->firstCcwF(e, v); }
-            void nextCcw(Arc &e) const { _digraph->nextCcwF(e); }
-            void turnLeft(Arc &e, const Node &v) const { _digraph->turnLeftF(e); }
-            void turnRight(Arc &e, const Node &v) const { _digraph->turnRightF(e); }
-            void first(Face &i) const { _digraph->first(i); }
-            void next(Face &i) const { _digraph->next(i); }
-            void firstCwF(Arc &arc, const Face& face) const { _digraph->firstCw(arc, face); }
-            void nextCwF(Arc &arc) const { _digraph->nextCw(arc); }
-
-            static int id(Node v) { return v.id; }
-            static int id(Arc e) { return e.id; }
-            static int id(Face f) { return f.id; }
-            static Node nodeFromId(int id) { return Node(id);}
-            static Arc arcFromId(int id) { return Arc(id);}
-            static Face faceFromId(int id) { return Face(id);}
-
-            bool valid(Node n) const { return _digraph->valid(n); }
-            bool valid(Arc n) const { return _digraph->valid(n); }
-            bool valid(Face n) const { return _digraph->valid(n); }
-
-        };
-
-        typedef PlanarDigraphExtender<DigraphExtender<DualBase> > ExtendedDualBase;
-        /// \brief Adaptor class for the dual of a planar graph.
-        ///
-        /// Adaptor class for the dual of a planar graph.
-        class Dual : public ExtendedDualBase {
-                public:
-            Dual(const PlanarDigraph &digraph) { initialize(&digraph); }
-
-        };
-
-        /// \brief Class to make a snapshot of the digraph and restore
-        /// it later.
-        ///
-        /// Class to make a snapshot of the digraph and restore it later.
-        ///
-        /// The newly added nodes and arcs can be removed using the
-        /// restore() function.
-        ///
-        /// \note After a state is restored, you cannot restore a later state,
-        /// i.e. you cannot add the removed nodes and arcs again using
-        /// another Snapshot instance.
-        ///
-        /// \warning Node and arc deletions and other modifications (e.g.
-        /// reversing, contracting, splitting arcs or nodes) cannot be
-        /// restored. These events invalidate the snapshot.
-        /// However, the arcs and nodes that were added to the digraph after
-        /// making the current snapshot can be removed without invalidating it.
-        class Snapshot {
-        protected:
-
-            typedef Parent::NodeNotifier NodeNotifier;
-
-            class NodeObserverProxy : public NodeNotifier::ObserverBase {
-            public:
-
-                NodeObserverProxy(Snapshot& _snapshot)
-                    : snapshot(_snapshot) {}
-
-                using NodeNotifier::ObserverBase::attach;
-                using NodeNotifier::ObserverBase::detach;
-                using NodeNotifier::ObserverBase::attached;
-
-            protected:
-
-                virtual void add(const Node& node) {
-                    snapshot.addNode(node);
-                }
-                virtual void add(const std::vector<Node>& nodes) {
-                    for (int i = nodes.size() - 1; i >= 0; ++i) {
-                        snapshot.addNode(nodes[i]);
-                    }
-                }
-                virtual void erase(const Node& node) {
-                    snapshot.eraseNode(node);
-                }
-                virtual void erase(const std::vector<Node>& nodes) {
-                    for (int i = 0; i < int(nodes.size()); ++i) {
-                        snapshot.eraseNode(nodes[i]);
-                    }
-                }
-                virtual void build() {
-                    Node node;
-                    std::vector<Node> nodes;
-                    for (notifier()->first(node); node != INVALID;
-                    notifier()->next(node)) {
-                        nodes.push_back(node);
-                    }
-                    for (int i = nodes.size() - 1; i >= 0; --i) {
-                        snapshot.addNode(nodes[i]);
-                    }
-                }
-                virtual void clear() {
-                    Node node;
-                    for (notifier()->first(node); node != INVALID;
-                    notifier()->next(node)) {
-                        snapshot.eraseNode(node);
-                    }
-                }
-
-                Snapshot& snapshot;
-            };
-
-            class ArcObserverProxy : public ArcNotifier::ObserverBase {
-            public:
-
-                ArcObserverProxy(Snapshot& _snapshot)
-                    : snapshot(_snapshot) {}
-
-                using ArcNotifier::ObserverBase::attach;
-                using ArcNotifier::ObserverBase::detach;
-                using ArcNotifier::ObserverBase::attached;
-
-            protected:
-
-                virtual void add(const Arc& arc) {
-                    snapshot.addArc(arc);
-                }
-                virtual void add(const std::vector<Arc>& arcs) {
-                    for (int i = arcs.size() - 1; i >= 0; ++i) {
-                        snapshot.addArc(arcs[i]);
-                    }
-                }
-                virtual void erase(const Arc& arc) {
-                    snapshot.eraseArc(arc);
-                }
-                virtual void erase(const std::vector<Arc>& arcs) {
-                    for (int i = 0; i < int(arcs.size()); ++i) {
-                        snapshot.eraseArc(arcs[i]);
-                    }
-                }
-                virtual void build() {
-                    Arc arc;
-                    std::vector<Arc> arcs;
-                    for (notifier()->first(arc); arc != INVALID;
-                    notifier()->next(arc)) {
-                        arcs.push_back(arc);
-                    }
-                    for (int i = arcs.size() - 1; i >= 0; --i) {
-                        snapshot.addArc(arcs[i]);
-                    }
-                }
-                virtual void clear() {
-                    Arc arc;
-                    for (notifier()->first(arc); arc != INVALID;
-                    notifier()->next(arc)) {
-                        snapshot.eraseArc(arc);
-                    }
-                }
-
-                Snapshot& snapshot;
-            };
-
-            PlanarDigraph *digraph;
-
-            NodeObserverProxy node_observer_proxy;
-            ArcObserverProxy arc_observer_proxy;
-
-            std::list<Node> added_nodes;
-            std::list<Arc> added_arcs;
-
-
-            void addNode(const Node& node) {
-                added_nodes.push_front(node);
-            }
-            void eraseNode(const Node& node) {
-                std::list<Node>::iterator it =
-                        std::find(added_nodes.begin(), added_nodes.end(), node);
-                if (it == added_nodes.end()) {
-                    clear();
-                    arc_observer_proxy.detach();
-                    throw NodeNotifier::ImmediateDetach();
-                } else {
-                    added_nodes.erase(it);
-                }
-            }
-
-            void addArc(const Arc& arc) {
-                added_arcs.push_front(arc);
-            }
-            void eraseArc(const Arc& arc) {
-                std::list<Arc>::iterator it =
-                        std::find(added_arcs.begin(), added_arcs.end(), arc);
-                if (it == added_arcs.end()) {
-                    clear();
-                    node_observer_proxy.detach();
-                    throw ArcNotifier::ImmediateDetach();
-                } else {
-                    added_arcs.erase(it);
-                }
-            }
-
-            void attach(PlanarDigraph &_digraph) {
-                digraph = &_digraph;
-                node_observer_proxy.attach(digraph->notifier(Node()));
-                arc_observer_proxy.attach(digraph->notifier(Arc()));
-            }
-
-            void detach() {
-                node_observer_proxy.detach();
-                arc_observer_proxy.detach();
-            }
-
-            bool attached() const {
-                return node_observer_proxy.attached();
-            }
-
-            void clear() {
-                added_nodes.clear();
-                added_arcs.clear();
-            }
-
-        public:
-
-            /// \brief Default constructor.
-            ///
-            /// Default constructor.
-            /// You have to call save() to actually make a snapshot.
-            Snapshot()
-                : digraph(0), node_observer_proxy(*this),
-                arc_observer_proxy(*this) {}
-
-            /// \brief Constructor that immediately makes a snapshot.
-            ///
-            /// This constructor immediately makes a snapshot of the given digraph.
-            Snapshot(PlanarDigraph &gr)
-                : node_observer_proxy(*this),
-                arc_observer_proxy(*this) {
-                attach(gr);
-            }
-
-            /// \brief Make a snapshot.
-            ///
-            /// This function makes a snapshot of the given digraph.
-            /// It can be called more than once. In case of a repeated
-            /// call, the previous snapshot gets lost.
-            void save(PlanarDigraph &gr) {
-                if (attached()) {
-                    detach();
-                    clear();
-                }
-                attach(gr);
-            }
-
-            /// \brief Undo the changes until the last snapshot.
-            ///
-            /// This function undos the changes until the last snapshot
-            /// created by save() or Snapshot(ListDigraph&).
-            ///
-            /// \warning This method invalidates the snapshot, i.e. repeated
-            /// restoring is not supported unless you call save() again.
-            void restore() {
-                detach();
-                for(std::list<Arc>::iterator it = added_arcs.begin();
-                it != added_arcs.end(); ++it) {
-                    digraph->erase(*it);
-                }
-                for(std::list<Node>::iterator it = added_nodes.begin();
-                it != added_nodes.end(); ++it) {
-                    digraph->erase(*it);
-                }
-                clear();
-            }
-
-            /// \brief Returns \c true if the snapshot is valid.
-            ///
-            /// This function returns \c true if the snapshot is valid.
-            bool valid() const {
-                return attached();
-            }
-        };
-
-    };
-
-    ///@}
-
-    class PlanarGraphBase {
-
-    protected:
-
-        struct NodeT {
-            int first_out, last_out;
-            int prev, next;
-            int component;
-        };
-
-        struct ArcT {
-            int target;
-            int left_face;
-            int prev_out, next_out;
-        };
-
-        struct FaceT {
-            int first_arc;
-            int prev, next;
-        };
-
-        std::vector<NodeT> nodes;
-
-        int first_node;
-
-        int first_free_node;
-
-        std::vector<ArcT> arcs;
-
-        int first_free_arc;
-
-        std::vector<FaceT> faces;
-
-        int first_face;
-
-        int first_free_face;
-
-        int component_id;
-
-    public:
-
-        typedef PlanarGraphBase Graph;
-
-        class Node {
-            friend class PlanarGraphBase;
-        protected:
-
-            int id;
-            explicit Node(int pid) { id = pid;}
-
-        public:
-            Node() {}
-            Node (Invalid) { id = -1; }
-            bool operator==(const Node& node) const {return id == node.id;}
-            bool operator!=(const Node& node) const {return id != node.id;}
-            bool operator<(const Node& node) const {return id < node.id;}
-        };
-
-        class Edge {
-            friend class PlanarGraphBase;
-        protected:
-
-            int id;
-            explicit Edge(int pid) { id = pid;}
-
-        public:
-            Edge() {}
-            Edge (Invalid) { id = -1; }
-            bool operator==(const Edge& edge) const {return id == edge.id;}
-            bool operator!=(const Edge& edge) const {return id != edge.id;}
-            bool operator<(const Edge& edge) const {return id < edge.id;}
-        };
-
-        class Arc {
-            friend class PlanarGraphBase;
-        protected:
-
-            int id;
-            explicit Arc(int pid) { id = pid;}
-
-        public:
-            operator Edge() const {
-                return id != -1 ? edgeFromId(id / 2) : INVALID;
-            }
-
-            Arc() {}
-            Arc (Invalid) { id = -1; }
-            bool operator==(const Arc& arc) const {return id == arc.id;}
-            bool operator!=(const Arc& arc) const {return id != arc.id;}
-            bool operator<(const Arc& arc) const {return id < arc.id;}
-        };
-
-        class Face {
-            friend class PlanarGraphBase;
-        protected:
-
-            int id;
-            explicit Face(int pid) { id = pid;}
-
-        public:
-            Face() {}
-            Face (Invalid) { id = -1; }
-            bool operator==(const Face& face) const {return id == face.id;}
-            bool operator!=(const Face& face) const {return id != face.id;}
-            bool operator<(const Face& face) const {return id < face.id;}
-        };
-
-        PlanarGraphBase()
-            : nodes(), first_node(-1), first_free_node(-1),
-              arcs(), first_free_arc(-1),
-              faces(), first_face(-1), first_free_face(-1),
-              component_id(0) {
-        }
-
-
-        int maxNodeId() const { return nodes.size()-1; }
-        int maxEdgeId() const { return arcs.size() / 2 - 1; }
-        int maxArcId() const { return arcs.size()-1; }
-        int maxFaceId() const { return faces.size()-1; }
-
-        Node source(Arc e) const { return Node(arcs[e.id ^ 1].target); }
-        Node target(Arc e) const { return Node(arcs[e.id].target); }
-        Face leftFace(Arc e) const { return Face(arcs[e.id].left_face); }
-        Face rightFace(Arc e) const { return Face(arcs[e.id ^ 1].left_face); }
-
-        Node u(Edge e) const { return Node(arcs[2 * e.id].target); }
-        Node v(Edge e) const { return Node(arcs[2 * e.id + 1].target); }
-        Face w1(Edge e) const { return Face(arcs[2 * e.id].left_face); }
-        Face w2(Edge e) const { return Face(arcs[2 * e.id + 1].left_face); }
-
-        static bool direction(Arc e) {
-            return (e.id & 1) == 1;
-        }
-
-        static Arc direct(Edge e, bool d) {
-            return Arc(e.id * 2 + (d ? 1 : 0));
-        }
-
-        //Primitives to use in iterators
-        void first(Node& node) const {
-            node.id = first_node;
-        }
-
-        void next(Node& node) const {
-            node.id = nodes[node.id].next;
-        }
-
-        void first(Arc& e) const {
-            int n = first_node;
-            while (n != -1 && nodes[n].first_out == -1) {
-                n = nodes[n].next;
-            }
-            e.id = (n == -1) ? -1 : nodes[n].first_out;
-        }
-
-        void next(Arc& e) const {
-            if (arcs[e.id].next_out != -1) {
-                e.id = arcs[e.id].next_out;
-            } else {
-                int n = nodes[arcs[e.id ^ 1].target].next;
-                while(n != -1 && nodes[n].first_out == -1) {
-                    n = nodes[n].next;
-                }
-                e.id = (n == -1) ? -1 : nodes[n].first_out;
-            }
-        }
-
-        void first(Edge& e) const {
-            int n = first_node;
-            while (n != -1) {
-                e.id = nodes[n].first_out;
-                while ((e.id & 1) != 1) {
-                    e.id = arcs[e.id].next_out;
-                }
-                if (e.id != -1) {
-                    e.id /= 2;
-                    return;
-                }
-                n = nodes[n].next;
-            }
-            e.id = -1;
-        }
-
-        void next(Edge& e) const {
-            int n = arcs[e.id * 2].target;
-            e.id = arcs[(e.id * 2) | 1].next_out;
-            while ((e.id & 1) != 1) {
-                e.id = arcs[e.id].next_out;
-            }
-            if (e.id != -1) {
-                e.id /= 2;
-                return;
-            }
-            n = nodes[n].next;
-            while (n != -1) {
-                e.id = nodes[n].first_out;
-                while ((e.id & 1) != 1) {
-                    e.id = arcs[e.id].next_out;
-                }
-                if (e.id != -1) {
-                    e.id /= 2;
-                    return;
-                }
-                n = nodes[n].next;
-            }
-            e.id = -1;
-        }
-
-        void firstOut(Arc &e, const Node& v) const {
-            e.id = nodes[v.id].first_out;
-        }
-        void nextOut(Arc &e) const {
-            e.id = arcs[e.id].next_out;
-        }
-
-        void firstIn(Arc &e, const Node& v) const {
-            e.id = ((nodes[v.id].first_out) ^ 1);
-            if (e.id == -2) e.id = -1;
-        }
-        void nextIn(Arc &e) const {
-            e.id = ((arcs[e.id ^ 1].next_out) ^ 1);
-            if (e.id == -2) e.id = -1;
-        }
-        void lastIn(Arc &e, const Node& v) const {
-            e.id = ((nodes[v.id].last_out) ^ 1);
-            if (e.id == -2) e.id = -1;
-        }
-        void prevIn(Arc &e) const {
-            e.id = ((arcs[e.id ^ 1].prev_out) ^ 1);
-            if (e.id == -2) e.id = -1;
-        }
-
-        void firstCwF(Arc &e, const Face &f) const {
-            e.id = faces[f.id].first_arc;
-        }
-        void nextCwF(Arc &e) const {
-            turnLeft(e);
-            setToOpposite(e);
-            if (e.id == faces[arcs[e.id].left_face].first_arc) e = INVALID;
-        }
-
-        void firstInF(Arc &e, const Face &f) const {
-            e.id = faces[f.id].first_arc;
-        }
-        void nextInF(Arc &e) const {
-            setToOpposite(e);
-            turnRight(e);
-            if (e.id == faces[arcs[e.id].left_face].first_arc) e = INVALID;
-        }
-        void lastInF(Arc &e, const Face &f) const {
-            e.id = faces[f.id].first_arc;
-            setToOpposite(e);
-            turnRight(e);
-        }
-        void prevInF(Arc &e) const {
-            if (e.id == faces[arcs[e.id].left_face].first_arc) {
-                e = INVALID;
-                return;
-            }
-            setToOpposite(e);
-            turnRight(e);
-        }
-
-        void firstOutF(Arc &e, const Face &f) const {
-            e.id = faces[e.id].first_arc ^ 1;
-        }
-        void nextOutF(Arc &e) const {
-            turnRight(e);
-            setToOpposite(e);
-            if (e.id == faces[arcs[e.id ^ 1].left_face].first_arc) e = INVALID;
-        }
-
-        void first(Arc &arc, const Face& face) const {
-            arc.id = faces[face.id].first_arc;
-        }
-
-        void turnLeftF(Arc &e) const {
-            setToOpposite(e);
-            turnRight(e);
-        }
-
-        void turnRightF(Arc &e) const {
-            turnLeft(e);
-            setToOpposite(e);
-        }
-
-        void turnLeft(Arc &e) const {
-            if (arcs[e.id].next_out > -1) {
-                e.id = arcs[e.id].next_out;
-            } else {
-                e.id = nodes[arcs[e.id ^ 1].target].first_out;
-            }
-        }
-        void turnRight(Arc &e) const {
-            if (arcs[e.id].prev_out > -1) {
-                e.id = arcs[e.id].prev_out;
-            } else {
-                e.id = nodes[arcs[e.id ^ 1].target].last_out;
-            }
-        }
-        void setToOpposite(Arc &a) const {
-            if (a.id != -1) a.id ^= 1;
-        }
-
-        void firstInc(Edge &e, bool& d, const Node& v) const {
-            int a = nodes[v.id].first_out;
-            if (a != -1 ) {
-                e.id = a / 2;
-                d = ((a & 1) == 1);
-            } else {
-                e.id = -1;
-                d = true;
-            }
-        }
-        void nextInc(Edge &e, bool& d) const {
-            int a = (arcs[(e.id * 2) | (d ? 1 : 0)].next_out);
-            if (a != -1 ) {
-                e.id = a / 2;
-                d = ((a & 1) == 1);
-            } else {
-                e.id = -1;
-                d = true;
-            }
-        }
-
-        void first(Face& face) const {
-            face.id = first_face;
-        }
-
-        void next(Face& face) const {
-            face.id = faces[face.id].next;
-        }
-
-        Arc arcAt(Node n, Edge e) {
-            if (e.id == -1) return INVALID;
-            return Arc((e.id*2) | (arcs[e.id*2].target == n.id?1:0));
-        }
-
-        static int id(Node v) { return v.id; }
-        static int id(Arc e) { return e.id; }
-        static int id(Edge e) { return e.id; }
-        static int id(Face f) { return f.id; }
-
-        static Node nodeFromId(int id) { return Node(id);}
-        static Arc arcFromId(int id) { return Arc(id);}
-        static Edge edgeFromId(int id) { return Edge(id);}
-        static Face faceFromId(int id) { return Face(id);}
-
-        bool valid(Node n) const {
-            return n.id >= 0 && n.id < static_cast<int>(nodes.size()) &&
-                    nodes[n.id].prev != -2;
-        }
-
-        bool valid(Arc a) const {
-            return a.id >= 0 && a.id < static_cast<int>(arcs.size()) &&
-                    arcs[a.id].prev_out != -2;
-        }
-
-        bool valid(Edge e) const {
-            return e.id >= 0 && 2 * e.id < static_cast<int>(arcs.size()) &&
-                    arcs[2 * e.id].prev_out != -2;
-        }
-
-        bool valid(Face f) const {
-            return f.id >= 0 && f.id < static_cast<int>(faces.size()) &&
-                    faces[f.id].prev != -2;
-        }
-
-        Node addNode() {
-            int n;
-
-            if(first_free_node==-1) {
-                n = nodes.size();
-                nodes.push_back(NodeT());
-            } else {
-                n = first_free_node;
-                first_free_node = nodes[n].next;
-            }
-
-            nodes[n].next = first_node;
-            nodes[n].component = component_id++;
-            if (first_node != -1) nodes[first_node].prev = n;
-            first_node = n;
-            nodes[n].prev = -1;
-
-            nodes[n].first_out = nodes[n].last_out = -1;
-
-            return Node(n);
-        }
-
-        Edge addEdge(Node u, Node v, Edge e_u, Edge e_v) {
-
-            Arc p_u = arcAt(u,e_u);
-            Arc p_v = arcAt(v,e_v);
-
-            if (p_u.id > -1 && p_v.id > -1 && arcs[p_u.id].left_face != arcs[p_v.id].left_face
-                && nodes[u.id].component == nodes[v.id].component) return INVALID;
-
-            int n = addBlankEdge();
-
-            arcs[n].target = u.id;
-            arcs[n | 1].target = v.id;
-
-            arcs[n].prev_out = p_v.id;
-            if (p_v.id > -1) {
-                arcs[n].next_out = arcs[p_v.id].next_out;
-                arcs[p_v.id].next_out = n;
-            } else {
-                arcs[n].next_out = nodes[v.id].first_out;
-                nodes[v.id].first_out = n;
-            }
-            if (arcs[n].next_out > -1) {
-                arcs[arcs[n].next_out].prev_out = n;
-            } else {
-                nodes[v.id].last_out = n;
-            }
-
-            arcs[n | 1].prev_out = p_u.id;
-            if (p_u.id > -1) {
-                arcs[n | 1].next_out = arcs[p_u.id].next_out;
-                arcs[p_u.id].next_out = n | 1;
-            } else {
-                arcs[n | 1].next_out = nodes[u.id].first_out;
-                nodes[u.id].first_out = n | 1;
-            }
-            if (arcs[n | 1].next_out > -1) {
-                arcs[arcs[n | 1].next_out].prev_out = n | 1;
-            } else {
-                nodes[u.id].last_out = n | 1;
-            }
-
-            //Add the extra face, if needed
-            if (p_u.id > -1 & p_v.id > -1) {
-                int oldf = arcs[p_u.id].left_face;
-                int oldfb = arcs[p_v.id].left_face;
-                arcs[n].left_face = arcs[n | 1].left_face = oldf;
-                Face f = addFace();
-                faces[f.id].first_arc = n | 1;
-                faces[oldf].first_arc = n;
-                Arc arc(n | 1);
-                wall_paint(arc,f.id,arc);
-                if (nodes[v.id].component != nodes[u.id].component) {
-                    erase(Face(oldf));
-                    erase(Face(oldfb));
-                    int ca = nodes[u.id].component;
-                    int cb = nodes[v.id].component;
-                    int k = first_node;
-                    while (k != -1) {
-                        if (nodes[k].component == cb)
-                            nodes[k].component = ca;
-                        k = nodes[k].next;
-                    }
-                }
-            } else if (p_u.id > -1) {
-                arcs[n].left_face = arcs[n | 1].left_face = arcs[p_u.id].left_face;
-                faces[arcs[n].left_face].first_arc = n | 1;
-                nodes[v.id].component = nodes[u.id].component;
-            } else if (p_v.id > -1) {
-                arcs[n].left_face = arcs[n | 1].left_face = arcs[p_v.id].left_face;
-                faces[arcs[n].left_face].first_arc = n | 1;
-                nodes[u.id].component = nodes[v.id].component;
-            } else {    //both prevs are INVALID
-                Face f = addFace();
-                arcs[n].left_face = arcs[n | 1].left_face = f.id;
-                faces[f.id].first_arc = n | 1;
-                nodes[v.id].component = nodes[u.id].component;
-            }
-
-            return Edge(n / 2);
-        }
-
-        void erase(const Node& node) {
-            int n = node.id;
-
-            if(nodes[n].next != -1) {
-                nodes[nodes[n].next].prev = nodes[n].prev;
-            }
-
-            if(nodes[n].prev != -1) {
-                nodes[nodes[n].prev].next = nodes[n].next;
-            } else {
-                first_node = nodes[n].next;
-            }
-
-            nodes[n].next = first_free_node;
-            first_free_node = n;
-            nodes[n].prev = -2;
-        }
-
-        void erase(const Edge& edge) {
-            int n = edge.id * 2;
-
-            //"retreat" the incident faces' first arcs
-            int fl = arcs[n].left_face;
-            if ((faces[fl].first_arc | 1) == (n | 1)) {
-                Arc e(faces[fl].first_arc);
-                turnRightF(e);
-                if ((e.id | 1) == (n | 1)) turnRightF(e);
-                faces[fl].first_arc = e.id;
-            }
-
-            int fr = arcs[n | 1].left_face;
-
-            bool comp_split = false;
-            if (fr != fl) {
-                Arc arc(faces[fr].first_arc);
-                wall_paint(arc,fl,arc);
-                erase(Face(fr));
-            } else if ((arcs[n].next_out > -1 || arcs[n].prev_out > -1) &&
-                (arcs[n | 1].next_out > -1 || arcs[n | 1].prev_out > -1)) {
-                comp_split = true;
-                Arc arc(n);
-                Arc ed = arc;
-                ed.id ^= 1;
-                turnRightF(arc);
-                Face f = addFace();
-                wall_paint(arc,f.id,ed);
-                faces[f.id].first_arc = arc.id;
-            }
-
-            if (arcs[n].next_out != -1) {
-                arcs[arcs[n].next_out].prev_out = arcs[n].prev_out;
-            } else {
-                nodes[arcs[n].target].last_out = arcs[n].prev_out;
-            }
-
-            if (arcs[n].prev_out != -1) {
-                arcs[arcs[n].prev_out].next_out = arcs[n].next_out;
-            } else {
-                nodes[arcs[n | 1].target].first_out = arcs[n].next_out;
-            }
-
-            if (arcs[n | 1].next_out != -1) {
-                arcs[arcs[n | 1].next_out].prev_out = arcs[n | 1].prev_out;
-            } else {
-                nodes[arcs[n | 1].target].last_out = arcs[n | 1].prev_out;
-            }
-
-            if (arcs[n | 1].prev_out != -1) {
-                arcs[arcs[n | 1].prev_out].next_out = arcs[n | 1].next_out;
-            } else {
-                nodes[arcs[n].target].first_out = arcs[n | 1].next_out;
-            }
-
-            arcs[n].next_out = first_free_arc;
-            first_free_arc = n;
-            arcs[n].prev_out = -2;
-            arcs[n | 1].prev_out = -2;
-
-            if (comp_split) component_relabel(Node(arcs[n | 1].target), component_id++);
-
-        }
-
-        void clear() {
-            arcs.clear();
-            nodes.clear();
-            faces.clear();
-            first_node = first_free_node = first_free_arc = first_face = first_free_face = -1;
-        }
-
-        Node split(Edge e) {
-            Node v = addNode();
-            Arc a(e.id*2);
-            int b = addBlankEdge();
-
-            nodes[v.id].component = nodes[arcs[a.id].target].component;
-            nodes[v.id].first_out = a.id;
-            nodes[v.id].last_out = b | 1;
-
-            arcs[b] = arcs[a.id];
-            arcs[b].target = v.id;
-            if (arcs[a.id].next_out > -1)
-                arcs[arcs[a.id].next_out].prev_out = b;
-            else
-                nodes[arcs[a.id | 1].target].last_out = b;
-            if (arcs[a.id].prev_out > -1)
-                arcs[arcs[a.id].prev_out].next_out = b;
-            else
-                nodes[arcs[a.id | 1].target].first_out = b;
-
-            arcs[b | 1] = arcs[a.id | 1];
-            arcs[b | 1].next_out = -1;
-            arcs[b | 1].prev_out = a.id;
-
-            arcs[a.id].next_out = b | 1;
-            arcs[a.id].prev_out = -1;
-            arcs[a.id | 1].target = v.id;
-
-
-            return v;
-        }
-
-    protected:
-
-        void wall_paint(Arc arc, int f_id, Arc ed) {
-            do {
-                arcs[arc.id].left_face = f_id;
-                turnRightF(arc);
-            } while (arc != ed);
-        }
-
-        void component_relabel(Node node, int comp_id) {
-            std::vector<int> ns(nodes.size());
-            std::list<int> q;
-            q.push_back(node.id);
-            ns[node.id] = 1;
-            while (!q.empty()) {
-                int n = q.front();
-                ns[n] = 2;
-                nodes[n].component = comp_id;
-                q.pop_front();
-                Arc arc;
-                firstOut(arc,Node(n));
-                while (arc.id > -1) {
-                    int m = arcs[arc.id].target;
-                    if (ns[m] == 0) {
-                        ns[m] = 1;
-                        q.push_back(m);
-                    }
-                    nextOut(arc);
-                }
-            }
-        }
-
-        Face addFace() {
-            int n;
-
-            if(first_free_face==-1) {
-                n = faces.size();
-                faces.push_back(FaceT());
-            } else {
-                n = first_free_face;
-                first_free_face = faces[n].next;
-            }
-
-            faces[n].next = first_face;
-            if (first_face != -1) faces[first_face].prev = n;
-            first_face = n;
-            faces[n].prev = -1;
-
-            faces[n].first_arc = -1;
-
-            return Face(n);
-        }
-
-        void erase(const Face& face) {
-            int n = face.id;
-
-            if(faces[n].next != -1) {
-                faces[faces[n].next].prev = faces[n].prev;
-            }
-
-            if(faces[n].prev != -1) {
-                faces[faces[n].prev].next = faces[n].next;
-            } else {
-                first_face = faces[n].next;
-            }
-
-            faces[n].next = first_free_face;
-            first_free_face = n;
-            faces[n].prev = -2;
-        }
-
-        int addBlankEdge() {
-            int n;
-            if (first_free_arc == -1) {
-                n = arcs.size();
-                arcs.push_back(ArcT());
-                arcs.push_back(ArcT());
-            } else {
-                n = first_free_arc;
-                first_free_arc = arcs[n].next_out;
-            }
-            return n;
-        }
-
-#ifdef REMOVE_BEFORE_RELEASE
-    public:
-        void print() {
-            std::cout << "Nodes: " << std::endl;
-            for (int i=0; i<nodes.size(); ++i)
-                std::cout << i << ":"
-                          << " fo=" << nodes[i].first_out
-                          << " pr=" << nodes[i].prev
-                          << " nx=" << nodes[i].next
-                          << " co=" << nodes[i].component
-                          <<std::endl;
-            std::cout << "Arcs: " << std::endl;
-            for (int i=0; i<arcs.size(); ++i) {
-                if (arcs[i].next_out > -2) {
-                    std::cout << i << ":"
-                          << " tg=" << arcs[i].target
-                          << " po=" << arcs[i].prev_out
-                          << " no=" << arcs[i].next_out
-                          << " lf=" << arcs[i].left_face
-                          <<std::endl;
-                } else std::cout << i << ": (deleted)" << std::endl;
-            }
-            std::cout << "Faces: " << std::endl;
-            for (int i=0; i<faces.size(); ++i)
-                std::cout << i
-                          << " pr=" << faces[i].prev
-                          << " nx=" << faces[i].next
-                          << " fa=" << faces[i].first_arc
-                          <<std::endl;
-        }
-#endif
-
-    };
-
-    template<typename Base>
-    class PlanarGraphExtender : public Base{
-
-        typedef Base Parent;
-
-            public:
-        typedef PlanarGraphExtender Graph;
-
-        PlanarGraphExtender() {}
-
-        typedef typename Parent::Node Node;
-        typedef typename Parent::Arc Arc;
-        typedef typename Parent::Edge Edge;
-        typedef typename Parent::Face Face;
-
-        class FaceIt : public Face {
-            const Graph* _graph;
-        public:
-
-            FaceIt() {}
-
-            FaceIt(Invalid i) : Face(i) { }
-
-            explicit FaceIt(const Graph& graph) : _graph(&graph) {
-                _graph->first(static_cast<Face&>(*this));
-            }
-
-            FaceIt(const Graph& graph, const Face& face)
-                : Face(face), _graph(&graph) {}
-
-            FaceIt& operator++() {
-                _graph->next(*this);
-                return *this;
-            }
-
-        };
-
-
-        class CwPerimeterArcIt : public Arc {
-          const Graph* _graph;
-          Face _face;
-          Arc f_arc;
-        public:
-
-          CwPerimeterArcIt() { }
-
-          CwPerimeterArcIt(Invalid i) : Arc(i) { }
-
-          CwPerimeterArcIt(const Graph& graph, const Face& face)
-            : _graph(&graph), _face(face) {
-            _graph->firstCwF(static_cast<Arc&>(*this),face);
-            f_arc = *this;
-          }
-
-          CwPerimeterArcIt(const Graph& graph, const Arc& arc)
-            : Arc(arc), _graph(&graph) {}
-
-          CwPerimeterArcIt& operator++() {
-            _graph->nextCwF(*this);
-            return *this;
-          }
-
-        };
-
-        class CcwArcIt : public Arc {
-          const Graph* _graph;
-          const Node _node;
-        public:
-
-          CcwArcIt() { }
-
-          CcwArcIt(Invalid i) : Arc(i) { }
-
-          CcwArcIt(const Graph& graph, const Node& node)
-            : _graph(&graph), _node(node) {
-            _graph->firstCcw(*this, node);
-          }
-
-          CcwArcIt(const Graph& graph, const Arc& arc)
-            : Arc(arc), _graph(&graph) {}
-
-          CcwArcIt& operator++() {
-            _graph->nextCcw(*this, _node);
-            return *this;
-          }
-
-        };
-
-    };
-
-    typedef PlanarGraphExtender<GraphExtender<PlanarGraphBase> > ExtendedPlanarGraphBase;
-
-
-    /// \addtogroup graphs
-    /// @{
-
-    ///A general undirected graph structure.
-
-    ///\ref ListGraph is a versatile and fast undirected graph
-    ///implementation based on linked lists that are stored in
-    ///\c std::vector structures.
-    ///
-    ///This type fully conforms to the \ref concepts::Graph "Graph concept"
-    ///and it also provides several useful additional functionalities.
-    ///Most of its member functions and nested classes are documented
-    ///only in the concept class.
-    ///
-    ///This class provides only linear time counting for nodes, edges and arcs.
-    ///
-    ///\sa concepts::Graph
-    ///\sa ListDigraph
-    class PlanarGraph : public ExtendedPlanarGraphBase {
-        typedef ExtendedPlanarGraphBase Parent;
-
-    private:
-        /// Graphs are \e not copy constructible. Use GraphCopy instead.
-        PlanarGraph(const PlanarGraph &) :ExtendedPlanarGraphBase()  {};
-        /// \brief Assignment of a graph to another one is \e not allowed.
-        /// Use GraphCopy instead.
-        void operator=(const PlanarGraph &) {}
-    public:
-        /// Constructor
-
-        /// Constructor.
-        ///
-        PlanarGraph() {}
-
-        typedef Parent::OutArcIt IncEdgeIt;
-
-        /// \brief Add a new node to the graph.
-        ///
-        /// This function adds a new node to the graph.
-        /// \return The new node.
-        Node addNode() { return Parent::addNode(); }
-
-        /// \brief Add a new edge to the graph.
-        ///
-        /// This function adds a new edge to the graph between nodes
-        /// \c u and \c v with inherent orientation from node \c u to
-        /// node \c v.
-        /// \return The new edge.
-        Edge addEdge(Node u, Node v, Edge p_u, Edge p_v) {
-            return PlanarGraphBase::addEdge(u, v, p_u, p_v);
-        }
-
-        ///\brief Erase a node from the graph.
-        ///
-        /// This function erases the given node along with its incident arcs
-        /// from the graph.
-        ///
-        /// \note All iterators referencing the removed node or the incident
-        /// edges are invalidated, of course.
-        void erase(Node n) { Parent::erase(n); }
-
-        ///\brief Erase an edge from the graph.
-        ///
-        /// This function erases the given edge from the graph.
-        ///
-        /// \note All iterators referencing the removed edge are invalidated,
-        /// of course.
-        void erase(Edge e) { Parent::erase(e); }
-        /// Node validity check
-
-        /// This function gives back \c true if the given node is valid,
-        /// i.e. it is a real node of the graph.
-        ///
-        /// \warning A removed node could become valid again if new nodes are
-        /// added to the graph.
-        bool valid(Node n) const { return Parent::valid(n); }
-        /// Edge validity check
-
-        /// This function gives back \c true if the given edge is valid,
-        /// i.e. it is a real edge of the graph.
-        ///
-        /// \warning A removed edge could become valid again if new edges are
-        /// added to the graph.
-        bool valid(Edge e) const { return Parent::valid(e); }
-        /// Arc validity check
-
-        /// This function gives back \c true if the given arc is valid,
-        /// i.e. it is a real arc of the graph.
-        ///
-        /// \warning A removed arc could become valid again if new edges are
-        /// added to the graph.
-        bool valid(Arc a) const { return Parent::valid(a); }
-
-        /// \brief Change the first node of an edge.
-        ///
-        /// Planar graphs don't support changing the endpoints of edges.
-        void changeU(Edge e, Node n) = delete;
-        /// \brief Change the second node of an edge.
-        ///
-        /// Planar graphs don't support changing the endpoints of edges.
-        void changeV(Edge e, Node n) = delete;
-
-        /// \brief Contract two nodes.
-        ///
-        /// This function contracts the given two nodes.
-        /// Node \c b is removed, but instead of deleting
-        /// its incident edges, they are joined to node \c a.
-        /// If the last parameter \c r is \c true (this is the default value),
-        /// then the newly created loops are removed.
-        ///
-        /// \note The moved edges are joined to node \c a using changeU()
-        /// or changeV(), thus all edge and arc iterators whose base node is
-        /// \c b are invalidated.
-        /// Moreover all iterators referencing node \c b or the removed
-        /// loops are also invalidated. Other iterators remain valid.
-        ///
-        ///\warning This functionality cannot be used together with the
-        ///Snapshot feature.
-/*TODO: rewrite this function
-        void contract(Node a, Node b, bool r = true) {
-            for(IncEdgeIt e(*this, b); e!=INVALID;) {
-                IncEdgeIt f = e; ++f;
-                if (r && runningNode(e) == a) {
-                    erase(e);
-                } else if (u(e) == b) {
-                    changeU(e, a);
-                } else {
-                    changeV(e, a);
-                }
-                e = f;
-            }
-            erase(b);
-        }
-*/
-
-        ///Clear the graph.
-
-        ///This function erases all nodes and arcs from the graph.
-        ///
-        ///\note All iterators of the graph are invalidated, of course.
-        void clear() {
-            Parent::clear();
-        }
-
-        /// Reserve memory for nodes.
-
-        /// Using this function, it is possible to avoid superfluous memory
-        /// allocation: if you know that the graph you want to build will
-        /// be large (e.g. it will contain millions of nodes and/or edges),
-        /// then it is worth reserving space for this amount before starting
-        /// to build the graph.
-        /// \sa reserveEdge()
-        void reserveNode(int n) { nodes.reserve(n); };
-
-        /// Reserve memory for edges.
-
-        /// Using this function, it is possible to avoid superfluous memory
-        /// allocation: if you know that the graph you want to build will
-        /// be large (e.g. it will contain millions of nodes and/or edges),
-        /// then it is worth reserving space for this amount before starting
-        /// to build the graph.
-        /// \sa reserveNode()
-        void reserveEdge(int m) { arcs.reserve(2 * m); };
-
-        class DualBase {
-            const Graph *_graph;
-        protected:
-            void initialize(const Graph *graph) { _graph = graph; }
-        public:
-
-            typedef PlanarGraph::Face Node;
-            typedef PlanarGraph::Arc Arc;
-            typedef PlanarGraph::Edge Edge;
-            typedef PlanarGraph::Node Face;
-
-            int maxNodeId() const { return _graph->maxFaceId(); }
-            int maxArcId() const { return _graph->maxArcId(); }
-            int maxFaceId() const { return _graph->maxNodeId(); }
-
-            Node source(Arc e) const { return _graph->leftFace(e); }
-            Node target(Arc e) const { return _graph->rightFace(e); }
-            Face leftFace(Arc e) const { return _graph->target(e); }
-            Face rightFace(Arc e) const { return _graph->source(e); }
-            Arc direct(const Edge &edge, const Node &node) const { return _graph->direct(edge, _graph->w1(edge) == node); }
-
-            void first(Node &i) const { _graph->first(i); }
-            void next(Node &i) const { _graph->next(i); }
-            void first(Arc &i) const { _graph->first(i); }
-            void next(Arc &i) const { _graph->next(i); }
-            void firstCcw(Arc& i, const Node& n) const { _graph->lastInF(i, n); }
-            void nextCcw(Arc& i, const Node &n) const { _graph->prevInF(i); }
-            void firstIn(Arc& i, const Node& n) const { _graph->firstInF(i, n); }
-            void nextIn(Arc& i) const { _graph->nextInF(i); }
-            void firstCwF(Arc& i, const Face& n) const { _graph->lastIn(i, n); }
-            void nextCwF(Arc& i) const { _graph->prevIn(i); }
-            void firstOut(Arc& i, const Node& n ) const { _graph->firstOutF(i, n); }
-            void nextOut(Arc& i) const { _graph->nextOutF(i); }
-            void first(Face &i) const { _graph->first(i); }
-            void next(Face &i) const { _graph->next(i); }
-
-            static int id(Node v) { return PlanarGraph::id(v); }
-            static int id(Arc e) { return PlanarGraph::id(e); }
-            static int id(Face f) { return PlanarGraph::id(f); }
-            static Node nodeFromId(int id) { return PlanarGraph::faceFromId(id);}
-            static Arc arcFromId(int id) { return PlanarGraph::arcFromId(id);}
-            static Face faceFromId(int id) { return PlanarGraph::nodeFromId(id);}
-
-            bool valid(Node n) const { return _graph->valid(n); }
-            bool valid(Arc n) const { return _graph->valid(n); }
-            bool valid(Face n) const { return _graph->valid(n); }
-
-        };
-
-        typedef PlanarGraphExtender<GraphExtender<DualBase> > ExtendedDualBase;
-        /// \brief Adaptor class for the dual of a planar graph.
-        ///
-        /// Adaptor class for the dual of a planar graph.
-        class Dual : public ExtendedDualBase {
-                public:
-            Dual(const PlanarGraph &graph) { initialize(&graph); }
-
-        };
-        /// \brief Class to make a snapshot of the graph and restore
-        /// it later.
-        ///
-        /// Class to make a snapshot of the graph and restore it later.
-        ///
-        /// The newly added nodes and edges can be removed
-        /// using the restore() function.
-        ///
-        /// \note After a state is restored, you cannot restore a later state,
-        /// i.e. you cannot add the removed nodes and edges again using
-        /// another Snapshot instance.
-        ///
-        /// \warning Node and edge deletions and other modifications
-        /// (e.g. changing the end-nodes of edges or contracting nodes)
-        /// cannot be restored. These events invalidate the snapshot.
-        /// However, the edges and nodes that were added to the graph after
-        /// making the current snapshot can be removed without invalidating it.
-        class Snapshot {
-        protected:
-
-            typedef Parent::NodeNotifier NodeNotifier;
-
-            class NodeObserverProxy : public NodeNotifier::ObserverBase {
-            public:
-
-                NodeObserverProxy(Snapshot& _snapshot)
-                    : snapshot(_snapshot) {}
-
-                using NodeNotifier::ObserverBase::attach;
-                using NodeNotifier::ObserverBase::detach;
-                using NodeNotifier::ObserverBase::attached;
-
-            protected:
-
-                virtual void add(const Node& node) {
-                    snapshot.addNode(node);
-                }
-                virtual void add(const std::vector<Node>& nodes) {
-                    for (int i = nodes.size() - 1; i >= 0; ++i) {
-                        snapshot.addNode(nodes[i]);
-                    }
-                }
-                virtual void erase(const Node& node) {
-                    snapshot.eraseNode(node);
-                }
-                virtual void erase(const std::vector<Node>& nodes) {
-                    for (int i = 0; i < int(nodes.size()); ++i) {
-                        snapshot.eraseNode(nodes[i]);
-                    }
-                }
-                virtual void build() {
-                    Node node;
-                    std::vector<Node> nodes;
-                    for (notifier()->first(node); node != INVALID;
-                    notifier()->next(node)) {
-                        nodes.push_back(node);
-                    }
-                    for (int i = nodes.size() - 1; i >= 0; --i) {
-                        snapshot.addNode(nodes[i]);
-                    }
-                }
-                virtual void clear() {
-                    Node node;
-                    for (notifier()->first(node); node != INVALID;
-                    notifier()->next(node)) {
-                        snapshot.eraseNode(node);
-                    }
-                }
-
-                Snapshot& snapshot;
-            };
-
-            class EdgeObserverProxy : public EdgeNotifier::ObserverBase {
-            public:
-
-                EdgeObserverProxy(Snapshot& _snapshot)
-                    : snapshot(_snapshot) {}
-
-                using EdgeNotifier::ObserverBase::attach;
-                using EdgeNotifier::ObserverBase::detach;
-                using EdgeNotifier::ObserverBase::attached;
-
-            protected:
-
-                virtual void add(const Edge& edge) {
-                    snapshot.addEdge(edge);
-                }
-                virtual void add(const std::vector<Edge>& edges) {
-                    for (int i = edges.size() - 1; i >= 0; ++i) {
-                        snapshot.addEdge(edges[i]);
-                    }
-                }
-                virtual void erase(const Edge& edge) {
-                    snapshot.eraseEdge(edge);
-                }
-                virtual void erase(const std::vector<Edge>& edges) {
-                    for (int i = 0; i < int(edges.size()); ++i) {
-                        snapshot.eraseEdge(edges[i]);
-                    }
-                }
-                virtual void build() {
-                    Edge edge;
-                    std::vector<Edge> edges;
-                    for (notifier()->first(edge); edge != INVALID;
-                    notifier()->next(edge)) {
-                        edges.push_back(edge);
-                    }
-                    for (int i = edges.size() - 1; i >= 0; --i) {
-                        snapshot.addEdge(edges[i]);
-                    }
-                }
-                virtual void clear() {
-                    Edge edge;
-                    for (notifier()->first(edge); edge != INVALID;
-                    notifier()->next(edge)) {
-                        snapshot.eraseEdge(edge);
-                    }
-                }
-
-                Snapshot& snapshot;
-            };
-
-            PlanarGraph *graph;
-
-            NodeObserverProxy node_observer_proxy;
-            EdgeObserverProxy edge_observer_proxy;
-
-            std::list<Node> added_nodes;
-            std::list<Edge> added_edges;
-
-
-            void addNode(const Node& node) {
-                added_nodes.push_front(node);
-            }
-            void eraseNode(const Node& node) {
-                std::list<Node>::iterator it =
-                        std::find(added_nodes.begin(), added_nodes.end(), node);
-                if (it == added_nodes.end()) {
-                    clear();
-                    edge_observer_proxy.detach();
-                    throw NodeNotifier::ImmediateDetach();
-                } else {
-                    added_nodes.erase(it);
-                }
-            }
-
-            void addEdge(const Edge& edge) {
-                added_edges.push_front(edge);
-            }
-            void eraseEdge(const Edge& edge) {
-                std::list<Edge>::iterator it =
-                        std::find(added_edges.begin(), added_edges.end(), edge);
-                if (it == added_edges.end()) {
-                    clear();
-                    node_observer_proxy.detach();
-                    throw EdgeNotifier::ImmediateDetach();
-                } else {
-                    added_edges.erase(it);
-                }
-            }
-
-            void attach(PlanarGraph &_graph) {
-                graph = &_graph;
-                node_observer_proxy.attach(graph->notifier(Node()));
-                edge_observer_proxy.attach(graph->notifier(Edge()));
-            }
-
-            void detach() {
-                node_observer_proxy.detach();
-                edge_observer_proxy.detach();
-            }
-
-            bool attached() const {
-                return node_observer_proxy.attached();
-            }
-
-            void clear() {
-                added_nodes.clear();
-                added_edges.clear();
-            }
-
-        public:
-
-            /// \brief Default constructor.
-            ///
-            /// Default constructor.
-            /// You have to call save() to actually make a snapshot.
-            Snapshot()
-                : graph(0), node_observer_proxy(*this),
-                edge_observer_proxy(*this) {}
-
-            /// \brief Constructor that immediately makes a snapshot.
-            ///
-            /// This constructor immediately makes a snapshot of the given graph.
-            Snapshot(PlanarGraph &gr)
-                : node_observer_proxy(*this),
-                edge_observer_proxy(*this) {
-                attach(gr);
-            }
-
-            /// \brief Make a snapshot.
-            ///
-            /// This function makes a snapshot of the given graph.
-            /// It can be called more than once. In case of a repeated
-            /// call, the previous snapshot gets lost.
-            void save(PlanarGraph &gr) {
-                if (attached()) {
-                    detach();
-                    clear();
-                }
-                attach(gr);
-            }
-
-            /// \brief Undo the changes until the last snapshot.
-            ///
-            /// This function undos the changes until the last snapshot
-            /// created by save() or Snapshot(ListGraph&).
-            ///
-            /// \warning This method invalidates the snapshot, i.e. repeated
-            /// restoring is not supported unless you call save() again.
-            void restore() {
-                detach();
-                for(std::list<Edge>::iterator it = added_edges.begin();
-                it != added_edges.end(); ++it) {
-                    graph->erase(*it);
-                }
-                for(std::list<Node>::iterator it = added_nodes.begin();
-                it != added_nodes.end(); ++it) {
-                    graph->erase(*it);
-                }
-                clear();
-            }
-
-            /// \brief Returns \c true if the snapshot is valid.
-            ///
-            /// This function returns \c true if the snapshot is valid.
-            bool valid() const {
-                return attached();
-            }
-        };
-    };
-
-    /// @}
+  /// @}
 } //namespace lemon
 
 
diff -r 4266cccd150b -r eb9a0cb25eea lemon/plane_graph.h
--- a/lemon/plane_graph.h	Thu Mar 18 21:56:09 2010 +0100
+++ b/lemon/plane_graph.h	Fri Mar 19 20:12:49 2010 +0100
@@ -21,7 +21,7 @@
 
 ///\ingroup graphs
 ///\file
-///\brief PlaneDigraph and PlaneGraph classes.
+///\brief PlaneGraph classes. UNDER CONSTRUCTION.
 
 #include <lemon/planar_graph.h>
 #include <lemon/dim2.h>
@@ -30,137 +30,112 @@
 
 namespace lemon {
 
-    class PlaneDigraph : public PlanarDigraph {
+  /// \addtogroup graphs
+  /// @{
 
-    private:
-        NodeMap<dim2::Point<double> > nodepos;
+  ///A data structure for a graph embedded in the plane.
 
-        static double angle(double x, double y) {
-            return atan2(y,x);
+  ///\ref PlaneGraph adds geometry features to PlanarGraph. While PlanarGraph
+  ///must be built from the topology, PlaneGraph requires coordinates for nodes
+  ///and determines the topology that fits the coordinates.
+  ///\sa PlanarGraph
+  class PlaneGraph : public PlanarGraph {
+
+  private:
+    NodeMap<dim2::Point<double> > nodepos;
+
+    static double angle(double x, double y) {
+      return atan2(y,x);
+    }
+  public:
+
+    ///\brief Constructor.
+
+    ///Constructor.
+    PlaneGraph() : nodepos(*this) {}
+
+    ///Add a new node to the graph.
+
+    ///PlaneGraph doesn't support this method. Use the overload that takes a
+    ///point parameter.
+    Node addNode() = delete;
+
+    ///Add a new edge to the graph.
+
+    ///PlaneGraph doesn't support this method. Use the overload that takes only
+    ///two node parameters.
+    Arc addEdge(Node, Node, Arc, Arc) = delete;
+
+    ///Add a new node to the graph.
+
+    ///Add a new node to the graph at the given coordinates. The coordinates for
+    ///different nodes should be distinct.
+    ///\return The new node.
+    Node addNode(const dim2::Point<double> &point) {
+      Node n = PlanarGraph::addNode();
+      nodepos[n] = point;
+      return n;
+    }
+
+    ///Add a new edge to the graph.
+
+    ///Add a new edge to the graph between the nodes \c n1 and \c n2. The
+    ///coordinates of \c n1 and \c n2 will be used to determine the correct
+    ///topology.
+    ///\return The new edge, or INVALID if the topology doesn't allow the
+    ///addition.
+    Edge addEdge(Node n1, Node n2) {
+      Arc p_u, p_v;
+      firstOut(p_u,n1);
+      if (p_u != INVALID) {
+        double ad = angle(nodepos[n2].x-nodepos[n1].x,nodepos[n2].y-nodepos[n1].
+          y);
+        double pd, nd;
+        pd = angle(nodepos[target(p_u)].x-nodepos[n1].x,nodepos[target(p_u)].y-
+          nodepos[n1].y);
+        nd = pd;
+        if (ad < pd) ad += 2*PI;
+        Arc n_u = p_u;
+        while (n_u != INVALID && ad > nd) {
+          p_u = n_u;
+          pd = nd;
+          nextOut(n_u);
+          if (n_u != INVALID)
+            nd = angle(nodepos[target(n_u)].x-nodepos[n1].x,nodepos[target(n_u)]
+              .y-nodepos[n1].y);
         }
-    public:
+      }
 
-        PlaneDigraph() : nodepos(*this) {}
-
-        Node addNode() = delete;
-        Arc addArc(Node, Node, Arc, Arc) = delete;
-
-        Node addNode(const dim2::Point<double> &point) {
-            Node n = PlanarDigraph::addNode();
-            nodepos[n] = point;
-            return n;
+      firstOut(p_v,n2);
+      if (p_v != INVALID) {
+        double ad = angle(nodepos[n1].x-nodepos[n2].x,nodepos[n1].y-nodepos[n2].
+          y);
+        double pd, nd;
+        pd = angle(nodepos[target(p_v)].x-nodepos[n2].x,nodepos[target(p_v)].y-
+          nodepos[n2].y);
+        nd = pd;
+        if (ad < pd) ad += 2*PI;
+        Arc n_v = p_v;
+        while (n_v != INVALID && ad > nd) {
+          p_v = n_v;
+          pd = nd;
+          nextOut(n_v);
+          if (n_v != INVALID)
+            nd = angle(nodepos[target(n_v)].x-nodepos[n2].x,nodepos[target(n_v)]
+            .y-nodepos[n2].y);
         }
-
-        Arc addArc(Node n1, Node n2) {
-            Arc p_u, p_v;
-            firstCcw(p_u,n1);
-            if (p_u != INVALID) {
-                double ad = angle(nodepos[n2].x-nodepos[n1].x,nodepos[n2].y-nodepos[n1].y);
-                double pd, nd;
-                pd = angle(nodepos[oppositeNode(n1,p_u)].x-nodepos[n1].x,nodepos[oppositeNode(n1,p_u)].y-nodepos[n1].y);
-                nd = pd;
-                if (ad < pd) ad += 2*PI;
-                Arc n_u = p_u;
-                while (n_u != INVALID && ad > nd) {
-                    p_u = n_u;
-                    pd = nd;
-                    nextCcw(n_u,n1);
-                    if (n_u != INVALID)
-                        nd = angle(nodepos[oppositeNode(n1,n_u)].x-nodepos[n1].x,nodepos[oppositeNode(n1,n_u)].y-nodepos[n1].y);
-                }
-            }
-
-            firstCcw(p_v,n2);
-            if (p_v != INVALID) {
-                double ad = angle(nodepos[n1].x-nodepos[n2].x,nodepos[n1].y-nodepos[n2].y);
-                double pd, nd;
-                pd = angle(nodepos[oppositeNode(n2,p_v)].x-nodepos[n2].x,nodepos[oppositeNode(n2,p_v)].y-nodepos[n2].y);
-                nd = pd;
-                if (ad < pd) ad += 2*PI;
-                Arc n_v = p_v;
-                while (n_v != INVALID && ad > nd) {
-                    p_v = n_v;
-                    pd = nd;
-                    nextCcw(n_v,n2);
-                    if (n_v != INVALID)
-                        nd = angle(nodepos[oppositeNode(n2,n_v)].x-nodepos[n2].x,nodepos[oppositeNode(n2,n_v)].y-nodepos[n2].y);
-                }
-            }
+      }
 
 #ifdef REMOVE_BEFORE_RELEASE
-//            std::cout << "AddArc: " << id(n1) << "->" << id(n2) << ", p_u: " << id(p_u) << ", p_v: " << id(p_v) << std::endl;
+//            std::cout << "AddArc: " << id(n1) << "->" << id(n2) << ", p_u: "
+//             << id(p_u) << ", p_v: " << id(p_v) << std::endl;
 #endif
-            return PlanarDigraph::addArc(n1,n2,p_u,p_v);
+      return PlanarGraph::addEdge(n1,n2,p_u,p_v);
 
-        }
+    }
 
-    };
-
-    class PlaneGraph : public PlanarGraph {
-
-    private:
-        NodeMap<dim2::Point<double> > nodepos;
-
-        static double angle(double x, double y) {
-            return atan2(y,x);
-        }
-    public:
-
-        PlaneGraph() : nodepos(*this) {}
-
-        Node addNode() = delete;
-        Arc addArc(Node, Node, Arc, Arc) = delete;
-
-        Node addNode(const dim2::Point<double> &point) {
-            Node n = PlanarGraph::addNode();
-            nodepos[n] = point;
-            return n;
-        }
-
-        Edge addEdge(Node n1, Node n2) {
-            Arc p_u, p_v;
-            firstOut(p_u,n1);
-            if (p_u != INVALID) {
-                double ad = angle(nodepos[n2].x-nodepos[n1].x,nodepos[n2].y-nodepos[n1].y);
-                double pd, nd;
-                pd = angle(nodepos[target(p_u)].x-nodepos[n1].x,nodepos[target(p_u)].y-nodepos[n1].y);
-                nd = pd;
-                if (ad < pd) ad += 2*PI;
-                Arc n_u = p_u;
-                while (n_u != INVALID && ad > nd) {
-                    p_u = n_u;
-                    pd = nd;
-                    nextOut(n_u);
-                    if (n_u != INVALID)
-                        nd = angle(nodepos[target(n_u)].x-nodepos[n1].x,nodepos[target(n_u)].y-nodepos[n1].y);
-                }
-            }
-
-            firstOut(p_v,n2);
-            if (p_v != INVALID) {
-                double ad = angle(nodepos[n1].x-nodepos[n2].x,nodepos[n1].y-nodepos[n2].y);
-                double pd, nd;
-                pd = angle(nodepos[target(p_v)].x-nodepos[n2].x,nodepos[target(p_v)].y-nodepos[n2].y);
-                nd = pd;
-                if (ad < pd) ad += 2*PI;
-                Arc n_v = p_v;
-                while (n_v != INVALID && ad > nd) {
-                    p_v = n_v;
-                    pd = nd;
-                    nextOut(n_v);
-                    if (n_v != INVALID)
-                        nd = angle(nodepos[target(n_v)].x-nodepos[n2].x,nodepos[target(n_v)].y-nodepos[n2].y);
-                }
-            }
-
-#ifdef REMOVE_BEFORE_RELEASE
-//            std::cout << "AddArc: " << id(n1) << "->" << id(n2) << ", p_u: " << id(p_u) << ", p_v: " << id(p_v) << std::endl;
-#endif
-            return PlanarGraph::addEdge(n1,n2,p_u,p_v);
-
-        }
-
-    };
+  };
+  /// @}
 
 } //namespace lemon