COIN-OR::LEMON - Graph Library

Ticket #363: lemon-sikgraf_ea56aa8243b1.patch

File lemon-sikgraf_ea56aa8243b1.patch, 128.0 KB (added by gyorokp, 9 years ago)
  • lemon/Makefile.am

    # HG changeset patch
    # User Peter Gyorok
    # Date 1275937071 -7200
    # Node ID ea56aa8243b16ac39e1834d1c83a72724d6553ee
    # Parent  0bca98cbebbbc65ed7baecaf7a4f68ab881e28e9
    Add PlanarGraph and PlaneGraph (#363)
    
    diff -r 0bca98cbebbb -r ea56aa8243b1 lemon/Makefile.am
    a b  
    110110        lemon/network_simplex.h \
    111111        lemon/pairing_heap.h \
    112112        lemon/path.h \
     113        lemon/planar_graph.h \
    113114        lemon/planarity.h \
     115        lemon/plane_graph.h \
    114116        lemon/preflow.h \
    115117        lemon/quad_heap.h \
    116118        lemon/radix_heap.h \
  • lemon/bits/bezier.h

    diff -r 0bca98cbebbb -r ea56aa8243b1 lemon/bits/bezier.h
    a b  
    1919#ifndef LEMON_BEZIER_H
    2020#define LEMON_BEZIER_H
    2121
     22#include <cmath>
     23
    2224//\ingroup misc
    2325//\file
    2426//\brief Classes to compute with Bezier curves.
     
    6668  Point norm() const { return rot90(p2-p1); }
    6769  Point grad(double) const { return grad(); }
    6870  Point norm(double t) const { return rot90(grad(t)); }
     71  Box<double> boundingBox() const {
     72    double minx, miny;
     73    double maxx, maxy;
     74    if (p1.x < p2.x) {
     75      minx = p1.x; maxx = p2.x;
     76    } else {
     77      minx = p2.x; maxx = p1.x;
     78    }
     79    if (p1.y < p2.y) {
     80      miny = p1.y; maxy = p2.y;
     81    } else {
     82      miny = p2.y; maxy = p1.y;
     83    }
     84    return Box<double>(minx, miny, maxx, maxy);
     85  }
    6986};
    7087
    7188class Bezier2 : public BezierBase
     
    100117  Bezier1 norm() const { return Bezier1(2.0*rot90(p2-p1),2.0*rot90(p3-p2)); }
    101118  Point grad(double t) const { return grad()(t); }
    102119  Point norm(double t) const { return rot90(grad(t)); }
     120  Box<double> boundingBox() const {
     121    double minx, miny;
     122    double maxx, maxy;
     123    if (p1.x < p3.x) {
     124      minx = p1.x; maxx = p3.x;
     125    } else {
     126      minx = p3.x; maxx = p1.x;
     127    }
     128    if (p1.y < p3.y) {
     129      miny = p1.y; maxy = p3.y;
     130    } else {
     131      miny = p3.y; maxy = p1.y;
     132    }
     133    if (p1.x + p3.x - 2 * p2.x != 0.0) {
     134      double t = (p1.x - p2.x) / (p1.x + p3.x - 2 * p2.x);
     135      if (t >= 0.0 && t <= 1.0) {
     136        double x = ((1-t)*(1-t))*p1.x+(2*(1-t)*t)*p2.x+(t*t)*p3.x;
     137        if (x < minx) minx = x;
     138        if (x > maxx) maxx = x;
     139      }
     140    }
     141    if (p1.y + p3.y - 2 * p2.y != 0.0) {
     142      double t = (p1.y - p2.y) / (p1.y + p3.y - 2 * p2.y);
     143      if (t >= 0.0 && t <= 1.0) {
     144        double y = ((1-t)*(1-t))*p1.y+(2*(1-t)*t)*p2.y+(t*t)*p3.y;
     145        if (y < miny) miny = y;
     146        if (y > maxy) maxy = y;
     147      }
     148    }
     149    return Box<double>(minx, miny, maxx, maxy);
     150  }
    103151};
    104152
    105153class Bezier3 : public BezierBase
     
    150198                                  3.0*rot90(p4-p3)); }
    151199  Point grad(double t) const { return grad()(t); }
    152200  Point norm(double t) const { return rot90(grad(t)); }
     201  Box<double> boundingBox() const {
     202    double minx, miny;
     203    double maxx, maxy;
     204    if (p1.x < p4.x) {
     205      minx = p1.x; maxx = p4.x;
     206    } else {
     207      minx = p4.x; maxx = p1.x;
     208    }
     209    if (p1.y < p4.y) {
     210      miny = p1.y; maxy = p4.y;
     211    } else {
     212      miny = p4.y; maxy = p1.y;
     213    }
     214    if (p1.x - 3 * p2.x + 3 * p3.x - p4.x != 0.0) {
     215      double disc = p2.x * p2.x + p3.x * p3.x - p1.x * p3.x
     216        + p1.x * p4.x - p2.x * p3.x - p2.x * p4.x;
     217      if (disc >= 0.0) {
     218        double t = (p1.x + p3.x - 2 * p2.x + std::sqrt(disc)) /
     219          (p1.x - 3 * p2.x + 3 * p3.x - p4.x);
     220        if (t >= 0.0 && t <= 1.0) {
     221          double x = ((1-t)*(1-t)*(1-t))*p1.x+(3*t*(1-t)*(1-t))*p2.x+
     222            (3*t*t*(1-t))*p3.x+(t*t*t)*p4.x;
     223          if (x < minx) minx = x;
     224          if (x > maxx) maxx = x;
     225        }
     226        t = (p1.x + p3.x - 2 * p2.x - std::sqrt(disc)) /
     227          (p1.x - 3 * p2.x + 3 * p3.x - p4.x);
     228        if (t >= 0.0 && t <= 1.0) {
     229          double x = ((1-t)*(1-t)*(1-t))*p1.x+(3*t*(1-t)*(1-t))*p2.x+
     230            (3*t*t*(1-t))*p3.x+(t*t*t)*p4.x;
     231          if (x < minx) minx = x;
     232          if (x > maxx) maxx = x;
     233        }
     234      }
     235    } else if (p1.x + p3.x - 2 * p2.x != 0.0) {
     236      double t = (p1.x - p2.x) / (2 * p1.x + 2 * p3.x - 4 * p2.x);
     237      if (t >= 0.0 && t <= 1.0) {
     238        double x = ((1-t)*(1-t)*(1-t))*p1.x+(3*t*(1-t)*(1-t))*p2.x+
     239          (3*t*t*(1-t))*p3.x+(t*t*t)*p4.x;
     240        if (x < minx) minx = x;
     241        if (x > maxx) maxx = x;
     242      }     
     243    }
     244    if (p1.y - 3 * p2.y + 3 * p3.y - p4.y != 0.0) {
     245      double disc = p2.y * p2.y + p3.y * p3.y - p1.y * p3.y
     246        + p1.y * p4.y - p2.y * p3.y - p2.y * p4.y;
     247      if (disc >= 0.0) {
     248        double t = (p1.y + p3.y - 2 * p2.y + std::sqrt(disc)) /
     249          (p1.y - 3 * p2.y + 3 * p3.y - p4.y);
     250        if (t >= 0.0 && t <= 1.0) {
     251          double y = ((1-t)*(1-t)*(1-t))*p1.y+(3*t*(1-t)*(1-t))*p2.y+
     252            (3*t*t*(1-t))*p3.y+(t*t*t)*p4.y;
     253          if (y < miny) miny = y;
     254          if (y > maxy) maxy = y;
     255        }
     256        t = (p1.y + p3.y - 2 * p2.y - std::sqrt(disc)) /
     257          (p1.y - 3 * p2.y + 3 * p3.y - p4.y);
     258        if (t >= 0.0 && t <= 1.0) {
     259          double y = ((1-t)*(1-t)*(1-t))*p1.y+(3*t*(1-t)*(1-t))*p2.y+
     260            (3*t*t*(1-t))*p3.y+(t*t*t)*p4.y;
     261          if (y < miny) miny = y;
     262          if (y > maxy) maxy = y;
     263        }
     264      }
     265    } else if (p1.y + p3.y - 2 * p2.y != 0.0) {
     266      double t = (p1.y - p2.y) / (2 * p1.y + 2 * p3.y - 4 * p2.y);
     267      if (t >= 0.0 && t <= 1.0) {
     268        double y = ((1-t)*(1-t)*(1-t))*p1.y+(3*t*(1-t)*(1-t))*p2.y+
     269          (3*t*t*(1-t))*p3.y+(t*t*t)*p4.y;
     270        if (y < miny) miny = y;
     271        if (y > maxy) maxy = y;
     272      }     
     273    }
     274    return Box<double>(minx, miny, maxx, maxy);
     275  }
    153276
    154277  template<class R,class F,class S,class D>
    155278  R recSplit(F &_f,const S &_s,D _d) const
  • lemon/graph_to_eps.h

    diff -r 0bca98cbebbb -r ea56aa8243b1 lemon/graph_to_eps.h
    a b  
    103103  double _arrowLength, _arrowWidth;
    104104
    105105  bool _showNodes, _showArcs;
     106  bool _customEdgeShapes;
    106107
    107108  bool _enableParallel;
    108109  double _parArcDist;
     
    154155    _nodeScale(.01), _xBorder(10), _yBorder(10), _scale(1.0),
    155156    _nodeBorderQuotient(.1),
    156157    _drawArrows(false), _arrowLength(1), _arrowWidth(0.3),
    157     _showNodes(true), _showArcs(true),
     158    _showNodes(true), _showArcs(true), _customEdgeShapes(false),
    158159    _enableParallel(false), _parArcDist(1),
    159160    _showNodeText(false), _nodeTexts(false), _nodeTextSize(1),
    160161    _showNodePsText(false), _nodePsTexts(false), _nodePsTextsPreamble(0),
     
    201202
    202203  using T::_showNodes;
    203204  using T::_showArcs;
     205  using T::_customEdgeShapes;
    204206
    205207  using T::_enableParallel;
    206208  using T::_parArcDist;
     
    579581  GraphToEps<T> &hideArcs(bool b=true) {_showArcs=!b;return *this;}
    580582  ///Hides the nodes
    581583  GraphToEps<T> &hideNodes(bool b=true) {_showNodes=!b;return *this;}
     584  ///Uses default edge shapes
     585  GraphToEps<T> &customEdgeShapes(bool b=false) {
     586    _customEdgeShapes=b;return *this;
     587  }
     588
    582589
    583590  ///Sets the size of the node texts
    584591  GraphToEps<T> &nodeTextSize(double d) {_nodeTextSize=d;return *this;}
     
    655662public:
    656663  ~GraphToEps() { }
    657664
    658   ///Draws the graph.
     665  template<typename GR, typename D>
     666  void edgeToEps(const GR &g, std::ostream &os, const typename GR::Arc e,
     667    const Color &col, double width, D t) {
     668  }
    659669
    660   ///Like other functions using
    661   ///\ref named-templ-func-param "named template parameters",
    662   ///this function calls the algorithm itself, i.e. in this case
    663   ///it draws the graph.
    664   void run() {
     670  template<typename GR>
     671  void edgeToEps(const GR &g, std::ostream &os, const typename GR::Arc e,
     672    const Color &col, double width, bool b) {
     673    g.edgeToEps(os,e,_arcColors[e],_arcWidths[e]*_arcWidthScale);
     674  }
     675
     676  template<typename D>
     677  void inner_run(D t) {
    665678    const double EPSILON=1e-9;
    666679    if(dontPrint) return;
    667680
     
    868881
    869882    if(_showArcs) {
    870883      os << "%Arcs:\ngsave\n";
    871       if(_enableParallel) {
     884      if (_customEdgeShapes) {
     885        for(ArcIt e(g);e!=INVALID;++e) {
     886          if((!_undirected||g.source(e) == g.target(e)||g.source(e)<g.target(e))
     887              &&_arcWidths[e]>0) {
     888            edgeToEps(g,os,e,_arcColors[e],_arcWidths[e]*_arcWidthScale,t);
     889          }
     890        }
     891      } else if (_enableParallel) {
    872892        std::vector<Arc> el;
    873893        for(ArcIt e(g);e!=INVALID;++e)
    874894          if((!_undirected||g.source(e)<g.target(e))&&_arcWidths[e]>0
     
    10561076    if(_pleaseRemoveOsStream) {delete &os;}
    10571077  }
    10581078
     1079  ///Draws the graph.
     1080
     1081  ///Like other functions using
     1082  ///\ref named-templ-func-param "named template parameters",
     1083  ///this function calls the algorithm itself, i.e. in this case
     1084  ///it draws the graph.
     1085  void run(int i = 0) {
     1086    this->inner_run(i); //for normal graphs
     1087  }
     1088
     1089  void run(bool b) {
     1090    this->inner_run(b); //for plane graphs with custom edge shape support
     1091  }
     1092
    10591093  ///\name Aliases
    10601094  ///These are just some aliases to other parameter setting functions.
    10611095
     
    11501184///\sa graphToEps(GR &g, std::ostream& os)
    11511185template<class GR>
    11521186GraphToEps<DefaultGraphToEpsTraits<GR> >
    1153 graphToEps(GR &g,const char *file_name)
     1187graphToEps(const GR &g,const char *file_name)
    11541188{
    11551189  std::ostream* os = new std::ofstream(file_name);
    11561190  if (!(*os)) {
  • new file lemon/planar_graph.h

    diff -r 0bca98cbebbb -r ea56aa8243b1 lemon/planar_graph.h
    - +  
     1/* -*- mode: C++; indent-tabs-mode: nil; -*-
     2 *
     3 * This file is a part of LEMON, a generic C++ optimization library.
     4 *
     5 * Copyright (C) 2003-2010
     6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8 *
     9 * Permission to use, modify and distribute this software is granted
     10 * provided that this copyright notice appears in all copies. For
     11 * precise terms see the accompanying LICENSE file.
     12 *
     13 * This software is provided "AS IS" with no warranty of any kind,
     14 * express or implied, and with no claim as to its suitability for any
     15 * purpose.
     16 *
     17 */
     18
     19#ifndef LEMON_PLANAR_GRAPH_H
     20#define LEMON_PLANAR_GRAPH_H
     21
     22///\ingroup graphs
     23///\file
     24///\brief PlanarGraph classes. UNDER CONSTRUCTION.
     25
     26#include <lemon/core.h>
     27#include <lemon/error.h>
     28#include <lemon/bits/graph_extender.h>
     29#include <lemon/planarity.h>
     30
     31#include <vector>
     32#include <list>
     33#include <set>
     34
     35#ifdef REMOVE_BEFORE_RELEASE
     36#include <iostream>
     37#endif
     38
     39namespace lemon {
     40
     41  class PlanarGraphBase {
     42
     43  protected:
     44
     45    struct NodeT {
     46      int first_out, last_out;
     47      int prev, next;
     48      int outer_face;
     49      int component;
     50    };
     51
     52    struct ArcT {
     53      int target;
     54      int left_face;
     55      int prev_out, next_out;
     56    };
     57
     58    struct FaceT {
     59      int first_arc;
     60      int prev, next;
     61    };
     62
     63    struct ComponentT {
     64      int num_nodes;
     65      int prev, next;
     66    };
     67
     68    std::vector<NodeT> nodes;
     69    int first_node;
     70    int first_free_node;
     71
     72    std::vector<ArcT> arcs;
     73    int first_free_arc;
     74
     75    std::vector<FaceT> faces;
     76    int first_face;
     77    int first_free_face;
     78
     79    std::vector<ComponentT> components;
     80    int first_component;
     81    int first_free_component;
     82
     83  public:
     84
     85    typedef PlanarGraphBase Graph;
     86
     87    class Node {
     88      friend class PlanarGraphBase;
     89    protected:
     90
     91      int id;
     92      explicit Node(int pid) {
     93        id = pid;
     94      }
     95
     96    public:
     97      Node() {}
     98      Node (Invalid) {
     99        id = -1;
     100      }
     101      bool operator==(const Node& node) const {
     102        return id == node.id;
     103      }
     104      bool operator!=(const Node& node) const {
     105        return id != node.id;
     106      }
     107      bool operator<(const Node& node) const {
     108        return id < node.id;
     109      }
     110    };
     111
     112    class Edge {
     113      friend class PlanarGraphBase;
     114    protected:
     115
     116      int id;
     117      explicit Edge(int pid) {
     118        id = pid;
     119      }
     120
     121    public:
     122      Edge() {}
     123      Edge (Invalid) {
     124        id = -1;
     125      }
     126      bool operator==(const Edge& edge) const {
     127        return id == edge.id;
     128      }
     129      bool operator!=(const Edge& edge) const {
     130        return id != edge.id;
     131      }
     132      bool operator<(const Edge& edge) const {
     133        return id < edge.id;
     134      }
     135    };
     136
     137    class Arc {
     138      friend class PlanarGraphBase;
     139    protected:
     140
     141      int id;
     142      explicit Arc(int pid) {
     143        id = pid;
     144      }
     145
     146    public:
     147      operator Edge() const {
     148        return id != -1 ? edgeFromId(id / 2) : INVALID;
     149      }
     150
     151      Arc() {}
     152      Arc (Invalid) {
     153        id = -1;
     154      }
     155      bool operator==(const Arc& arc) const {
     156        return id == arc.id;
     157      }
     158      bool operator!=(const Arc& arc) const {
     159        return id != arc.id;
     160      }
     161      bool operator<(const Arc& arc) const {
     162        return id < arc.id;
     163      }
     164    };
     165
     166    class Face {
     167      friend class PlanarGraphBase;
     168    protected:
     169
     170      int id;
     171      explicit Face(int pid) {
     172        id = pid;
     173      }
     174
     175    public:
     176      Face() {}
     177      Face (Invalid) {
     178        id = -1;
     179      }
     180      bool operator==(const Face& face) const {
     181        return id == face.id;
     182      }
     183      bool operator!=(const Face& face) const {
     184        return id != face.id;
     185      }
     186      bool operator<(const Face& face) const {
     187        return id < face.id;
     188      }
     189    };
     190
     191    PlanarGraphBase()
     192        : nodes(), first_node(-1), first_free_node(-1),
     193        arcs(), first_free_arc(-1),
     194        faces(), first_face(-1), first_free_face(-1),
     195        components(), first_component(-1), first_free_component(-1) {
     196    }
     197
     198
     199    int maxNodeId() const {
     200      return nodes.size()-1;
     201    }
     202    int maxEdgeId() const {
     203      return arcs.size() / 2 - 1;
     204    }
     205    int maxArcId() const {
     206      return arcs.size()-1;
     207    }
     208    int maxFaceId() const {
     209      return faces.size()-1;
     210    }
     211
     212    int numComponents() const {
     213      int result = 0;
     214      int c = first_component;
     215      while (c != -1) {
     216        ++result;
     217        c = components[c].next;
     218      }
     219      return result;
     220    }
     221
     222    Node source(Arc e) const {
     223      return Node(arcs[e.id ^ 1].target);
     224    }
     225    Node target(Arc e) const {
     226      return Node(arcs[e.id].target);
     227    }
     228    Face leftFace(Arc e) const {
     229      return Face(arcs[e.id].left_face);
     230    }
     231    Face rightFace(Arc e) const {
     232      return Face(arcs[e.id ^ 1].left_face);
     233    }
     234    Face outerFace(Node n) const {
     235      return Face(nodes[n.id].outer_face);
     236    }
     237
     238    Node u(Edge e) const {
     239      return Node(arcs[2 * e.id].target);
     240    }
     241    Node v(Edge e) const {
     242      return Node(arcs[2 * e.id + 1].target);
     243    }
     244    Face w1(Edge e) const {
     245      return Face(arcs[2 * e.id].left_face);
     246    }
     247    Face w2(Edge e) const {
     248      return Face(arcs[2 * e.id + 1].left_face);
     249    }
     250
     251    static bool direction(Arc e) {
     252      return (e.id & 1) == 1;
     253    }
     254
     255    static Arc direct(Edge e, bool d) {
     256      return Arc(e.id * 2 + (d ? 1 : 0));
     257    }
     258
     259    //Primitives to use in iterators
     260    void first(Node& node) const {
     261      node.id = first_node;
     262    }
     263
     264    void next(Node& node) const {
     265      node.id = nodes[node.id].next;
     266    }
     267
     268    void first(Arc& e) const {
     269      int n = first_node;
     270      while (n != -1 && nodes[n].first_out == -1) {
     271        n = nodes[n].next;
     272      }
     273      e.id = (n == -1) ? -1 : nodes[n].first_out;
     274    }
     275
     276    void next(Arc& e) const {
     277      if (arcs[e.id].next_out != -1) {
     278        e.id = arcs[e.id].next_out;
     279      } else {
     280        int n = nodes[arcs[e.id ^ 1].target].next;
     281        while (n != -1 && nodes[n].first_out == -1) {
     282          n = nodes[n].next;
     283        }
     284        e.id = (n == -1) ? -1 : nodes[n].first_out;
     285      }
     286    }
     287
     288    void first(Edge& e) const {
     289      int n = first_node;
     290      while (n != -1) {
     291        e.id = nodes[n].first_out;
     292        while ((e.id & 1) != 1) {
     293          e.id = arcs[e.id].next_out;
     294        }
     295        if (e.id != -1) {
     296          e.id /= 2;
     297          return;
     298        }
     299        n = nodes[n].next;
     300      }
     301      e.id = -1;
     302    }
     303
     304    void next(Edge& e) const {
     305      int n = arcs[e.id * 2].target;
     306      e.id = arcs[(e.id * 2) | 1].next_out;
     307      while ((e.id & 1) != 1) {
     308        e.id = arcs[e.id].next_out;
     309      }
     310      if (e.id != -1) {
     311        e.id /= 2;
     312        return;
     313      }
     314      n = nodes[n].next;
     315      while (n != -1) {
     316        e.id = nodes[n].first_out;
     317        while ((e.id & 1) != 1) {
     318          e.id = arcs[e.id].next_out;
     319        }
     320        if (e.id != -1) {
     321          e.id /= 2;
     322          return;
     323        }
     324        n = nodes[n].next;
     325      }
     326      e.id = -1;
     327    }
     328
     329    void firstOut(Arc &e, const Node& v) const {
     330      e.id = nodes[v.id].first_out;
     331    }
     332    void nextOut(Arc &e) const {
     333      e.id = arcs[e.id].next_out;
     334    }
     335
     336    void lastOut(Arc &e, const Node& v) const {
     337      e.id = nodes[v.id].last_out;
     338    }
     339    void prevOut(Arc &e) const {
     340      e.id = arcs[e.id].prev_out;
     341    }
     342
     343    void firstIn(Arc &e, const Node& v) const {
     344      e.id = ((nodes[v.id].first_out) ^ 1);
     345      if (e.id == -2) e.id = -1;
     346    }
     347    void nextIn(Arc &e) const {
     348      e.id = ((arcs[e.id ^ 1].next_out) ^ 1);
     349      if (e.id == -2) e.id = -1;
     350    }
     351    void lastIn(Arc &e, const Node& v) const {
     352      e.id = ((nodes[v.id].last_out) ^ 1);
     353      if (e.id == -2) e.id = -1;
     354    }
     355    void prevIn(Arc &e) const {
     356      e.id = ((arcs[e.id ^ 1].prev_out) ^ 1);
     357      if (e.id == -2) e.id = -1;
     358    }
     359
     360    void firstCcwF(Arc &e, const Face &f) const {
     361      e.id = faces[f.id].first_arc;
     362      turnRight(e);
     363      setToOpposite(e);
     364    }
     365    void nextCcwF(Arc &e) const {
     366      if (e.id == faces[arcs[e.id].left_face].first_arc) e = INVALID;
     367      else {
     368        turnRight(e);
     369        setToOpposite(e);
     370      }
     371    }
     372    void firstCwF(Arc &e, const Face &f) const {
     373      e.id = faces[f.id].first_arc;
     374    }
     375    void nextCwF(Arc &e) const {
     376      turnLeft(e);
     377      setToOpposite(e);
     378      if (e.id == faces[arcs[e.id].left_face].first_arc) e = INVALID;
     379    }
     380
     381    void firstInF(Arc &e, const Face &f) const {
     382      e.id = faces[f.id].first_arc;
     383    }
     384    void nextInF(Arc &e) const {
     385      setToOpposite(e);
     386      turnRight(e);
     387      if (e.id == faces[arcs[e.id].left_face].first_arc) e = INVALID;
     388    }
     389    void lastInF(Arc &e, const Face &f) const {
     390      e.id = faces[f.id].first_arc;
     391      setToOpposite(e);
     392      turnRight(e);
     393    }
     394    void prevInF(Arc &e) const {
     395      if (e.id == faces[arcs[e.id].left_face].first_arc) {
     396        e = INVALID;
     397        return;
     398      }
     399      setToOpposite(e);
     400      turnRight(e);
     401    }
     402
     403    void firstOutF(Arc &e, const Face &f) const {
     404      e.id = faces[e.id].first_arc ^ 1;
     405    }
     406    void nextOutF(Arc &e) const {
     407      turnRight(e);
     408      setToOpposite(e);
     409      if (e.id == faces[arcs[e.id ^ 1].left_face].first_arc) e = INVALID;
     410    }
     411
     412    void first(Arc &arc, const Face& face) const {
     413      arc.id = faces[face.id].first_arc;
     414    }
     415
     416    void turnLeftF(Arc &e) const {
     417      setToOpposite(e);
     418      turnRight(e);
     419    }
     420
     421    void turnRightF(Arc &e) const {
     422      turnLeft(e);
     423      setToOpposite(e);
     424    }
     425
     426    void turnLeft(Arc &e) const {
     427      if (arcs[e.id].next_out > -1) {
     428        e.id = arcs[e.id].next_out;
     429      } else {
     430        e.id = nodes[arcs[e.id ^ 1].target].first_out;
     431      }
     432    }
     433    void turnRight(Arc &e) const {
     434      if (arcs[e.id].prev_out > -1) {
     435        e.id = arcs[e.id].prev_out;
     436      } else {
     437        e.id = nodes[arcs[e.id ^ 1].target].last_out;
     438      }
     439    }
     440    void setToOpposite(Arc &a) const {
     441      if (a.id != -1) a.id ^= 1;
     442    }
     443
     444    void firstInc(Edge &e, bool& d, const Node& v) const {
     445      int a = nodes[v.id].first_out;
     446      if (a != -1 ) {
     447        e.id = a / 2;
     448        d = ((a & 1) == 1);
     449      } else {
     450        e.id = -1;
     451        d = true;
     452      }
     453    }
     454    void nextInc(Edge &e, bool& d) const {
     455      int a = (arcs[(e.id * 2) | (d ? 1 : 0)].next_out);
     456      if (a != -1 ) {
     457        e.id = a / 2;
     458        d = ((a & 1) == 1);
     459      } else {
     460        e.id = -1;
     461        d = true;
     462      }
     463    }
     464
     465    void first(Face& face) const {
     466      face.id = first_face;
     467    }
     468
     469    void next(Face& face) const {
     470      face.id = faces[face.id].next;
     471    }
     472
     473    Arc arcAt(Node n, Edge e) {
     474      if (e.id == -1) return INVALID;
     475      return Arc((e.id*2) | (arcs[e.id*2].target == n.id?1:0));
     476    }
     477
     478    static int id(Node v) {
     479      return v.id;
     480    }
     481    static int id(Arc e) {
     482      return e.id;
     483    }
     484    static int id(Edge e) {
     485      return e.id;
     486    }
     487    static int id(Face f) {
     488      return f.id;
     489    }
     490
     491    static Node nodeFromId(int id) {
     492      return Node(id);
     493    }
     494    static Arc arcFromId(int id) {
     495      return Arc(id);
     496    }
     497    static Edge edgeFromId(int id) {
     498      return Edge(id);
     499    }
     500    static Face faceFromId(int id) {
     501      return Face(id);
     502    }
     503
     504    bool valid(Node n) const {
     505      return n.id >= 0 && n.id < static_cast<int>(nodes.size()) &&
     506             nodes[n.id].prev != -2;
     507    }
     508
     509    bool valid(Arc a) const {
     510      return a.id >= 0 && a.id < static_cast<int>(arcs.size()) &&
     511             arcs[a.id].prev_out != -2;
     512    }
     513
     514    bool valid(Edge e) const {
     515      return e.id >= 0 && 2 * e.id < static_cast<int>(arcs.size()) &&
     516             arcs[2 * e.id].prev_out != -2;
     517    }
     518
     519    bool valid(Face f) const {
     520      return f.id >= 0 && f.id < static_cast<int>(faces.size()) &&
     521             faces[f.id].prev != -2;
     522    }
     523
     524    Node addNode() {
     525      int n = link_node();
     526
     527      nodes[n].component = addComponent();
     528      components[nodes[n].component].num_nodes = 1;
     529      nodes[n].outer_face = addFace().id;
     530
     531      return Node(n);
     532    }
     533
     534    Edge addEdge(Node u, Node v, Edge e_u, Edge e_v) {
     535
     536      Arc p_u = arcAt(u,e_u);
     537      Arc p_v = arcAt(v,e_v);
     538
     539      if (p_u.id > -1 && p_v.id > -1 && arcs[p_u.id].left_face != arcs[p_v.id].
     540          left_face && nodes[u.id].component == nodes[v.id].component) return
     541          INVALID;
     542
     543      int n = link_edge();
     544
     545      arcs[n].target = u.id;
     546      arcs[n | 1].target = v.id;
     547
     548      arcs[n].prev_out = p_v.id;
     549      if (p_v.id > -1) {
     550        arcs[n].next_out = arcs[p_v.id].next_out;
     551        arcs[p_v.id].next_out = n;
     552      } else {
     553        arcs[n].next_out = nodes[v.id].first_out;
     554        nodes[v.id].first_out = n;
     555      }
     556      if (arcs[n].next_out > -1) {
     557        arcs[arcs[n].next_out].prev_out = n;
     558      } else {
     559        nodes[v.id].last_out = n;
     560      }
     561
     562      arcs[n | 1].prev_out = p_u.id;
     563      if (p_u.id > -1) {
     564        arcs[n | 1].next_out = arcs[p_u.id].next_out;
     565        arcs[p_u.id].next_out = n | 1;
     566      } else {
     567        arcs[n | 1].next_out = nodes[u.id].first_out;
     568        nodes[u.id].first_out = n | 1;
     569      }
     570      if (arcs[n | 1].next_out > -1) {
     571        arcs[arcs[n | 1].next_out].prev_out = n | 1;
     572      } else {
     573        nodes[u.id].last_out = n | 1;
     574      }
     575
     576      int ca = nodes[u.id].component;
     577      int cb = nodes[v.id].component;
     578
     579      //Add the extra face, if needed
     580      if (p_u.id > -1 & p_v.id > -1) {
     581        int oldf = arcs[p_u.id].left_face;
     582        int oldfb = arcs[p_v.id].left_face;
     583        arcs[n].left_face = arcs[n | 1].left_face = oldf;
     584        Face f = addFace();
     585        faces[f.id].first_arc = n | 1;
     586        faces[oldf].first_arc = n;
     587        Arc arc(n | 1);
     588        wall_paint(arc,f.id,arc);
     589        if (nodes[v.id].component != nodes[u.id].component) {
     590          erase(Face(oldf));
     591          erase(Face(oldfb));
     592          if (components[ca].num_nodes > components[cb].num_nodes) {
     593            component_relabel(v,ca);
     594            eraseComponent(cb);
     595          } else {
     596            component_relabel(u,cb);
     597            eraseComponent(ca);
     598          }
     599        }
     600      } else if (p_u.id > -1) {
     601        arcs[n].left_face = arcs[n | 1].left_face = arcs[p_u.id].left_face;
     602        faces[arcs[n].left_face].first_arc = n | 1;
     603        nodes[v.id].component = nodes[u.id].component;
     604        ++components[ca].num_nodes;
     605        eraseComponent(cb);
     606        erase(Face(nodes[v.id].outer_face));
     607      } else if (p_v.id > -1) {
     608        arcs[n].left_face = arcs[n | 1].left_face = arcs[p_v.id].left_face;
     609        faces[arcs[n].left_face].first_arc = n | 1;
     610        ++components[cb].num_nodes;
     611        eraseComponent(ca);
     612        nodes[u.id].component = nodes[v.id].component;
     613        erase(Face(nodes[u.id].outer_face));
     614      } else if (u.id != v.id) {    //both prevs are INVALID
     615        Face f(nodes[u.id].outer_face);
     616        arcs[n].left_face = arcs[n | 1].left_face = f.id;
     617        faces[f.id].first_arc = n | 1;
     618        ++components[ca].num_nodes;
     619        eraseComponent(cb);
     620        nodes[v.id].component = nodes[u.id].component;
     621        erase(Face(nodes[v.id].outer_face));
     622      } else {  //insert a loop to a new node
     623        Face f(nodes[u.id].outer_face);
     624        Face f2 = addFace();
     625        arcs[n].left_face = f2.id;
     626        arcs[n | 1].left_face = f.id;
     627        faces[f.id].first_arc = n | 1;
     628        faces[f2.id].first_arc = n;
     629      }
     630
     631      nodes[u.id].outer_face = nodes[v.id].outer_face = -1;
     632      return Edge(n / 2);
     633    }
     634
     635    void erase(const Node& node) {
     636      int n = node.id;
     637
     638      eraseComponent(nodes[n].component);
     639      erase(Face(nodes[n].outer_face));
     640
     641      unlink_node(n);
     642    }
     643
     644    void erase(const Edge& edge) {
     645      int n = edge.id * 2;
     646
     647      //"retreat" the incident faces' first arcs
     648      int fl = arcs[n].left_face;
     649      if ((faces[fl].first_arc | 1) == (n | 1)) {
     650        Arc e(faces[fl].first_arc);
     651        turnRightF(e);
     652        if ((e.id | 1) == (n | 1)) turnRightF(e);
     653        faces[fl].first_arc = e.id;
     654      }
     655
     656      int fr = arcs[n | 1].left_face;
     657
     658      bool comp_split = false;
     659      if (fr != fl) {
     660        Arc arc(faces[fr].first_arc);
     661        wall_paint(arc,fl,arc);
     662        if ((faces[fl].first_arc | 1) == (n | 1))
     663          faces[fl].first_arc = faces[fr].first_arc;
     664        erase(Face(fr));
     665      } else {
     666        comp_split = true;
     667        if ((arcs[n].next_out > -1 || arcs[n].prev_out > -1) &&
     668                 (arcs[n | 1].next_out > -1 || arcs[n | 1].prev_out > -1)) {
     669          Arc arc(n);
     670          Arc ed = arc;
     671          ed.id ^= 1;
     672          turnRightF(arc);
     673          Face f = addFace();
     674          wall_paint(arc,f.id,ed);
     675          faces[f.id].first_arc = arc.id;
     676          turnRightF(ed);
     677          faces[fr].first_arc = ed.id;
     678        }
     679      }
     680
     681      if (arcs[n].next_out != -1) {
     682        arcs[arcs[n].next_out].prev_out = arcs[n].prev_out;
     683      } else {
     684        nodes[arcs[n].target].last_out = arcs[n].prev_out;
     685      }
     686
     687      if (arcs[n].prev_out != -1) {
     688        arcs[arcs[n].prev_out].next_out = arcs[n].next_out;
     689      } else {
     690        nodes[arcs[n | 1].target].first_out = arcs[n].next_out;
     691      }
     692
     693      if (arcs[n | 1].next_out != -1) {
     694        arcs[arcs[n | 1].next_out].prev_out = arcs[n | 1].prev_out;
     695      } else {
     696        nodes[arcs[n | 1].target].last_out = arcs[n | 1].prev_out;
     697      }
     698
     699      if (arcs[n | 1].prev_out != -1) {
     700        arcs[arcs[n | 1].prev_out].next_out = arcs[n | 1].next_out;
     701      } else {
     702        nodes[arcs[n].target].first_out = arcs[n | 1].next_out;
     703      }
     704
     705      unlink_edge(edge);
     706
     707      if (comp_split) component_relabel(Node(arcs[n | 1].target),
     708        addComponent());
     709      if (nodes[arcs[n].target].first_out == -1 && nodes[arcs[n | 1].target].
     710              first_out == -1) {
     711        int t1 = arcs[n].target;
     712        int t2 = arcs[n | 1].target;
     713        if (t1 != t2) {
     714          nodes[t1].outer_face = addFace().id;
     715          nodes[t2].outer_face = fr;
     716        } else {
     717          faces[fl].first_arc = -1;
     718          nodes[t1].outer_face = fl;
     719        }
     720      } else if (nodes[arcs[n].target].first_out == -1)
     721        nodes[arcs[n].target].outer_face = addFace().id;
     722      else if (nodes[arcs[n | 1].target].first_out == -1)
     723        nodes[arcs[n | 1].target].outer_face = addFace().id;
     724
     725    }
     726
     727    void clear() {
     728      arcs.clear();
     729      nodes.clear();
     730      faces.clear();
     731      components.clear();
     732      first_node = first_free_node = first_free_arc = first_face =
     733        first_free_face = first_component = first_free_component = -1;
     734    }
     735
     736    Node split(Edge e) {
     737      Node v;
     738      Edge e2;
     739      inner_split1(v,e2);
     740      return inner_split2(e,v,e2);
     741    }
     742
     743    Node split(Node n1, Edge e1, Edge e2, bool b) {
     744      Node n2;
     745      inner_split1(n1,n2);
     746      Edge a3;
     747      inner_split2(n1,n2,e1,e2,b,a3);
     748      if (!b && a3 != INVALID) erase(a3);
     749      return n2;
     750    }
     751
     752    void contract(Edge e) {
     753      Node n1 = u(e);
     754      Node n2 = v(e);
     755      Node nd;
     756      bool simple;
     757      inner_contract1(n1,n2,nd,simple);
     758      inner_contract2(e,n1,n2,nd,simple);
     759    }
     760
     761  protected:
     762
     763    //Adds a node to the linked list and sets default properties. Used by
     764    //addNode and other functions that take care of the planar properties by
     765    //themselves.
     766    int link_node() {
     767      int n;
     768      if (first_free_node==-1) {
     769        n = nodes.size();
     770        nodes.push_back(NodeT());
     771      } else {
     772        n = first_free_node;
     773        first_free_node = nodes[n].next;
     774      }
     775
     776      nodes[n].next = first_node;
     777      if (first_node != -1) nodes[first_node].prev = n;
     778      first_node = n;
     779      nodes[n].prev = -1;
     780      nodes[n].first_out = nodes[n].last_out = -1;
     781      return n;
     782    }
     783
     784    //Removes a node from the linked list. Used by
     785    //erase(Node) and other functions that take care of the planar properties by
     786    //themselves.
     787    void unlink_node(int n) {
     788      if (nodes[n].next != -1) {
     789        nodes[nodes[n].next].prev = nodes[n].prev;
     790      }
     791
     792      if (nodes[n].prev != -1) {
     793        nodes[nodes[n].prev].next = nodes[n].next;
     794      } else {
     795        first_node = nodes[n].next;
     796      }
     797
     798      nodes[n].next = first_free_node;
     799      first_free_node = n;
     800      nodes[n].prev = -2;
     801    }
     802
     803    //Adds an edge to the linked list. Used by
     804    //addEdge and other functions that take care of the planar properties by
     805    //themselves.
     806    int link_edge() {
     807      int n;
     808      if (first_free_arc == -1) {
     809        n = arcs.size();
     810        arcs.push_back(ArcT());
     811        arcs.push_back(ArcT());
     812      } else {
     813        n = first_free_arc;
     814        first_free_arc = arcs[n].next_out;
     815      }
     816      return n;
     817    }
     818
     819    //Removes an edge from the linked list. Used by
     820    //erase(Edge) and other functions that take care of the planar properties by
     821    //themselves.
     822    void unlink_edge(Edge e) {
     823      int n = e.id*2;
     824      arcs[n].next_out = first_free_arc;
     825      first_free_arc = n;
     826      arcs[n].prev_out = -2;
     827      arcs[n | 1].prev_out = -2;
     828    }
     829
     830    //Primitives of split(Edge)
     831    void inner_split1(Node &v, Edge &e2) {
     832      v = Node(link_node());
     833      e2 = Arc(link_edge());
     834    }
     835
     836    Node inner_split2(Edge e, Node v, Edge e2) {
     837      Arc a(e.id*2);
     838      int b = e2.id*2;
     839      nodes[v.id].component = nodes[arcs[a.id].target].component;
     840      ++components[nodes[v.id].component].num_nodes;
     841      nodes[v.id].first_out = a.id;
     842      nodes[v.id].last_out = b | 1;
     843
     844      arcs[b] = arcs[a.id];
     845      arcs[b].target = v.id;
     846      if (arcs[a.id].next_out > -1)
     847        arcs[arcs[a.id].next_out].prev_out = b;
     848      else
     849        nodes[arcs[a.id | 1].target].last_out = b;
     850      if (arcs[a.id].prev_out > -1)
     851        arcs[arcs[a.id].prev_out].next_out = b;
     852      else
     853        nodes[arcs[a.id | 1].target].first_out = b;
     854
     855      arcs[b | 1] = arcs[a.id | 1];
     856      arcs[b | 1].next_out = -1;
     857      arcs[b | 1].prev_out = a.id;
     858
     859      arcs[a.id].next_out = b | 1;
     860      arcs[a.id].prev_out = -1;
     861      arcs[a.id | 1].target = v.id;
     862
     863      return v;
     864    }
     865
     866    //Primitives of contract(Edge)
     867    void inner_contract1(Node n1, Node n2, Node &to_delete, bool &simple) {
     868      if (nodes[n1.id].first_out == nodes[n1.id].last_out) {
     869        simple = true;
     870        to_delete = n1;
     871      } else if (nodes[n2.id].first_out == nodes[n2.id].last_out) {
     872        simple = true;
     873        to_delete = n2;
     874      } else {
     875        simple = false;
     876        to_delete = n2;
     877      }
     878    }
     879
     880    void inner_contract2(Edge e, Node n1, Node n2, Node &to_delete, bool
     881      &simple) {
     882      if (simple) {
     883        erase(e);
     884        erase(to_delete);
     885        return;
     886      }
     887      Arc a(e.id*2);
     888      int fl = arcs[a.id].left_face;
     889      if ((faces[fl].first_arc | 1) == a.id | 1) {
     890        Arc b = a;
     891        turnRightF(b);
     892        faces[fl].first_arc = b.id;
     893      }
     894      int fr = arcs[a.id | 1].left_face;
     895      if ((faces[fr].first_arc | 1) == a.id | 1) {
     896        Arc b(a.id | 1);
     897        turnRightF(b);
     898        faces[fr].first_arc = b.id;
     899      }
     900      Arc b = a;
     901      turnLeft(b);
     902      while (b != a) {
     903        arcs[b.id ^ 1].target = n1.id;
     904        if (arcs[b.id].next_out > -1)
     905          b.id = arcs[b.id].next_out;
     906        else
     907          b.id = nodes[n2.id].first_out;
     908      }
     909      arcs[nodes[n2.id].last_out].next_out = nodes[n2.id].first_out;
     910      arcs[nodes[n2.id].first_out].prev_out = nodes[n2.id].last_out;
     911      if (arcs[a.id | 1].prev_out > -1) {
     912        arcs[arcs[a.id | 1].prev_out].next_out = arcs[a.id].next_out;
     913      } else {
     914        nodes[n1.id].first_out = arcs[a.id].next_out;
     915      }
     916      arcs[arcs[a.id].next_out].prev_out = arcs[a.id | 1].prev_out;
     917      if (arcs[a.id | 1].next_out > -1) {
     918        arcs[arcs[a.id | 1].next_out].prev_out = arcs[a.id].prev_out;
     919      } else {
     920        nodes[n1.id].last_out = arcs[a.id].prev_out;
     921      }
     922      arcs[arcs[a.id].prev_out].next_out = arcs[a.id | 1].next_out;
     923
     924      unlink_edge(e);
     925      --components[nodes[n2.id].component].num_nodes;
     926      unlink_node(n2.id);
     927    }
     928
     929    //Primitives of split(Node)
     930    void inner_split1(Node n1, Node &n2) {
     931      n2 = Node(link_node());
     932    }
     933
     934    void inner_split2(Node n1, Node n2, Edge e1, Edge e2, bool b, Edge &e3) {
     935      e3 = INVALID;
     936      if (nodes[n1.id].first_out == -1) {
     937        if (b) e3 = addEdge(n1,n2,INVALID,INVALID);
     938        return;
     939      }
     940      arcs[nodes[n1.id].first_out].prev_out = nodes[n1.id].last_out;
     941      arcs[nodes[n1.id].last_out].next_out = nodes[n1.id].first_out;
     942      nodes[n2.id].component = nodes[n1.id].component;
     943      ++components[nodes[n2.id].component].num_nodes;
     944      Arc a1(e1.id*2);
     945      if (arcs[a1.id].target != n1.id) a1.id |= 1;
     946      Arc a2(e2.id*2);
     947      if (arcs[a2.id].target != n1.id) a2.id |= 1;
     948      Arc a3(link_edge());
     949      e3 = a3;
     950      arcs[a3.id].target = n1.id;
     951      arcs[a3.id | 1].target = n2.id;
     952      arcs[a3.id].left_face = arcs[a1.id].left_face;
     953      arcs[a3.id | 1].left_face = arcs[a2.id].left_face;
     954
     955      arcs[a3.id].prev_out = arcs[a2.id ^ 1].prev_out;
     956      arcs[arcs[a2.id ^ 1].prev_out].next_out = a3.id;
     957      arcs[a3.id | 1].prev_out = arcs[a1.id ^ 1].prev_out;
     958      arcs[arcs[a1.id ^ 1].prev_out].next_out = a3.id | 1;
     959      nodes[n1.id].first_out = a2.id ^ 1;
     960      arcs[a1.id ^ 1].prev_out = -1;
     961      nodes[n1.id].last_out = a3.id | 1;
     962      arcs[a3.id | 1].next_out = -1;
     963      nodes[n2.id].first_out = a1.id ^ 1;
     964      arcs[a2.id ^ 1].prev_out = -1;
     965      nodes[n2.id].last_out = a3.id;
     966      arcs[a3.id].next_out = -1;
     967
     968      Arc a4(a1.id);
     969      while (a4.id != (a3.id | 1)) {
     970        arcs[a4.id].target = n2.id;
     971        a4.id = arcs[a4.id ^ 1].next_out ^ 1;
     972      }
     973    }
     974
     975    //Relabels the "wall" of a face with a new face ID.
     976    void wall_paint(Arc arc, int f_id, Arc ed) {
     977      do {
     978        arcs[arc.id].left_face = f_id;
     979        turnRightF(arc);
     980      } while (arc != ed);
     981    }
     982
     983    //Relabels a component with a new component ID.
     984    void component_relabel(Node node, int comp_id) {
     985      std::vector<int> ns(nodes.size());
     986      std::list<int> q;
     987      q.push_back(node.id);
     988      ns[node.id] = 1;
     989      int cnt = 0;
     990      while (!q.empty()) {
     991        int n = q.front();
     992        ns[n] = 2;
     993        --components[nodes[n].component].num_nodes;
     994        nodes[n].component = comp_id;
     995        ++cnt;
     996        q.pop_front();
     997        Arc arc;
     998        firstOut(arc,Node(n));
     999        while (arc.id > -1) {
     1000          int m = arcs[arc.id].target;
     1001          if (ns[m] == 0 && nodes[m].component != comp_id) {
     1002            ns[m] = 1;
     1003            q.push_back(m);
     1004          }
     1005          nextOut(arc);
     1006        }
     1007      }
     1008      components[comp_id].num_nodes += cnt;
     1009    }
     1010
     1011    Face addFace() {
     1012      int n;
     1013
     1014      if (first_free_face==-1) {
     1015        n = faces.size();
     1016        faces.push_back(FaceT());
     1017      } else {
     1018        n = first_free_face;
     1019        first_free_face = faces[n].next;
     1020      }
     1021
     1022      faces[n].next = first_face;
     1023      if (first_face != -1) faces[first_face].prev = n;
     1024      first_face = n;
     1025      faces[n].prev = -1;
     1026
     1027      faces[n].first_arc = -1;
     1028
     1029      return Face(n);
     1030    }
     1031
     1032    void erase(const Face& face) {
     1033      int n = face.id;
     1034
     1035      if (faces[n].next != -1) {
     1036        faces[faces[n].next].prev = faces[n].prev;
     1037      }
     1038
     1039      if (faces[n].prev != -1) {
     1040        faces[faces[n].prev].next = faces[n].next;
     1041      } else {
     1042        first_face = faces[n].next;
     1043      }
     1044
     1045      faces[n].next = first_free_face;
     1046      first_free_face = n;
     1047      faces[n].prev = -2;
     1048    }
     1049
     1050    int addComponent() {
     1051      int n;
     1052
     1053      if (first_free_component==-1) {
     1054        n = components.size();
     1055        components.push_back(ComponentT());
     1056      } else {
     1057        n = first_free_component;
     1058        first_free_component = components[n].next;
     1059      }
     1060
     1061      components[n].next = first_component;
     1062      if (first_component != -1) components[first_component].prev = n;
     1063      first_component = n;
     1064      components[n].prev = -1;
     1065
     1066      components[n].num_nodes = 0;
     1067
     1068      return n;
     1069    }
     1070
     1071    void eraseComponent(const int n) {
     1072
     1073      if (components[n].next != -1) {
     1074        components[components[n].next].prev = components[n].prev;
     1075      }
     1076
     1077      if (components[n].prev != -1) {
     1078        components[components[n].prev].next = components[n].next;
     1079      } else {
     1080        first_component = components[n].next;
     1081      }
     1082
     1083      components[n].next = first_free_component;
     1084      first_free_component = n;
     1085      components[n].prev = -2;
     1086    }
     1087
     1088    int idComponent(Node n)
     1089    {
     1090      return nodes[n.id].component;
     1091    }
     1092
     1093#ifdef REMOVE_BEFORE_RELEASE
     1094  public:
     1095    void print() {
     1096      std::cout << "Nodes: " << std::endl;
     1097      for (int i=0; i<nodes.size(); ++i)
     1098        std::cout << i << ":"
     1099        << " fo=" << nodes[i].first_out
     1100        << " lo=" << nodes[i].last_out
     1101        << " pr=" << nodes[i].prev
     1102        << " nx=" << nodes[i].next
     1103        << " co=" << nodes[i].component
     1104        << " of=" << nodes[i].outer_face
     1105        <<std::endl;
     1106      std::cout << "Arcs: " << std::endl;
     1107      for (int i=0; i<arcs.size(); ++i) {
     1108        if (arcs[i].next_out > -2) {
     1109          std::cout << i << ":"
     1110          << " tg=" << arcs[i].target
     1111          << " po=" << arcs[i].prev_out
     1112          << " no=" << arcs[i].next_out
     1113          << " lf=" << arcs[i].left_face
     1114          <<std::endl;
     1115        } else std::cout << i << ": (deleted)" << std::endl;
     1116      }
     1117      std::cout << "Faces: " << std::endl;
     1118      for (int i=0; i<faces.size(); ++i) {
     1119        if (faces[i].prev > -2) {
     1120          std::cout << i
     1121            << " pr=" << faces[i].prev
     1122            << " nx=" << faces[i].next
     1123            << " fa=" << faces[i].first_arc
     1124            <<std::endl;
     1125        } else std::cout << i << ": (deleted)" << std::endl;
     1126      }
     1127      std::cout << "Components: " << std::endl;
     1128      for (int i=0; i<components.size(); ++i) {
     1129        std::cout << i << ':';
     1130        if (components[i].prev > -2) {
     1131          std::cout
     1132            << " pr=" << components[i].prev
     1133            << " nx=" << components[i].next
     1134            << " nn=" << components[i].num_nodes
     1135            <<std::endl;
     1136        } else std::cout << " (deleted)" << std::endl;
     1137      }
     1138    }
     1139#endif
     1140
     1141  };
     1142
     1143  // Face counting:
     1144
     1145  namespace _core_bits {
     1146
     1147    template <typename Graph, typename Enable = void>
     1148    struct CountFacesSelector {
     1149      static int count(const Graph &g) {
     1150        return countItems<Graph, typename Graph::Face>(g);
     1151      }
     1152    };
     1153
     1154    template <typename Graph>
     1155    struct CountFacesSelector<
     1156      Graph, typename
     1157      enable_if<typename Graph::FaceNumTag, void>::type>
     1158    {
     1159      static int count(const Graph &g) {
     1160        return g.faceNum();
     1161      }
     1162    };
     1163  }
     1164
     1165  /// \brief Function to count the faces in the graph.
     1166  ///
     1167  /// This function counts the faces in the graph.
     1168  /// The complexity of the function is <em>O</em>(<em>n</em>), but for some
     1169  /// graph structures it is specialized to run in <em>O</em>(1).
     1170  ///
     1171  /// \note If the graph contains a \c faceNum() member function and a
     1172  /// \c FaceNumTag tag then this function calls directly the member
     1173  /// function to query the cardinality of the face set.
     1174  template <typename Graph>
     1175  inline int countFaces(const Graph& g) {
     1176    return _core_bits::CountFacesSelector<Graph>::count(g);
     1177  }
     1178
     1179  template <typename Graph, typename DegIt>
     1180  inline int countFaceDegree(const Graph& _g, const typename Graph::Face& _n) {
     1181    int num = 0;
     1182    for (DegIt it(_g, _n); it != INVALID; ++it) {
     1183      ++num;
     1184    }
     1185    return num;
     1186  }
     1187  /// \brief Function to count the number of the boundary arcs of face \c f.
     1188  ///
     1189  /// This function counts the number of the boundary arcs of face \c f
     1190  /// in the undirected graph \c g.
     1191  template <typename Graph>
     1192  inline int countBoundaryArcs(const Graph& g,  const typename Graph::Face& f) {
     1193    return countFaceDegree<Graph, typename Graph::CwBoundaryArcIt>(g, f);
     1194  }
     1195
     1196
     1197  template <typename GR, typename Enable = void>
     1198  struct FaceNotifierIndicator {
     1199    typedef InvalidType Type;
     1200  };
     1201  template <typename GR>
     1202  struct FaceNotifierIndicator<
     1203    GR,
     1204    typename enable_if<typename GR::FaceNotifier::Notifier, void>::type
     1205  > {
     1206    typedef typename GR::FaceNotifier Type;
     1207  };
     1208
     1209  template <typename GR>
     1210  class ItemSetTraits<GR, typename GR::Face> {
     1211  public:
     1212
     1213    typedef GR Graph;
     1214    typedef GR Digraph;
     1215
     1216    typedef typename GR::Face Item;
     1217    typedef typename GR::FaceIt ItemIt;
     1218
     1219    typedef typename FaceNotifierIndicator<GR>::Type ItemNotifier;
     1220
     1221    template <typename V>
     1222    class Map : public GR::template FaceMap<V> {
     1223      typedef typename GR::template FaceMap<V> Parent;
     1224
     1225    public:
     1226      typedef typename GR::template FaceMap<V> Type;
     1227      typedef typename Parent::Value Value;
     1228
     1229      Map(const GR& _digraph) : Parent(_digraph) {}
     1230      Map(const GR& _digraph, const Value& _value)
     1231        : Parent(_digraph, _value) {}
     1232
     1233     };
     1234
     1235  };
     1236
     1237  template<typename Base>
     1238  class PlanarGraphExtender : public Base{
     1239
     1240    typedef Base Parent;
     1241
     1242  public:
     1243    typedef PlanarGraphExtender Graph;
     1244
     1245    typedef typename Parent::Node Node;
     1246    typedef typename Parent::Arc Arc;
     1247    typedef typename Parent::Edge Edge;
     1248    typedef typename Parent::Face Face;
     1249
     1250    typedef typename Parent::NodeNotifier NodeNotifier;
     1251    typedef typename Parent::EdgeNotifier EdgeNotifier;
     1252    typedef typename Parent::ArcNotifier ArcNotifier;
     1253    typedef AlterationNotifier<PlanarGraphExtender, Face> FaceNotifier;
     1254
     1255  protected:
     1256    mutable FaceNotifier face_notifier;
     1257
     1258  public:
     1259
     1260    int maxId(Face) const {
     1261      return Parent::maxFaceId();
     1262    }
     1263
     1264    NodeNotifier& notifier(Node) const {
     1265      return Parent::node_notifier;
     1266    }
     1267
     1268    ArcNotifier& notifier(Arc) const {
     1269      return Parent::arc_notifier;
     1270    }
     1271
     1272    EdgeNotifier& notifier(Edge) const {
     1273      return Parent::edge_notifier;
     1274    }
     1275
     1276    FaceNotifier& notifier(Face) const {
     1277      return face_notifier;
     1278    }
     1279
     1280    template <typename _Value>
     1281    class FaceMap
     1282    : public MapExtender<DefaultMap<Graph, Face, _Value> > {
     1283      typedef MapExtender<DefaultMap<Graph, Face, _Value> > Parent;
     1284
     1285    public:
     1286      explicit FaceMap(const Graph& graph)
     1287      : Parent(graph) {}
     1288      FaceMap(const Graph& graph, const _Value& value)
     1289      : Parent(graph, value) {}
     1290
     1291    private:
     1292      FaceMap& operator=(const FaceMap& cmap) {
     1293        return operator=<FaceMap>(cmap);
     1294      }
     1295
     1296      template <typename CMap>
     1297      FaceMap& operator=(const CMap& cmap) {
     1298        Parent::operator=(cmap);
     1299        return *this;
     1300      }
     1301
     1302    };
     1303
     1304
     1305  /// Iterator class for the faces.
     1306
     1307  /// This iterator goes through the faces of a planar graph.
     1308  class FaceIt : public Face {
     1309      const Graph* _graph;
     1310    public:
     1311
     1312      FaceIt() {}
     1313
     1314      FaceIt(Invalid i) : Face(i) { }
     1315
     1316      explicit FaceIt(const Graph& graph) : _graph(&graph) {
     1317        _graph->first(static_cast<Face&>(*this));
     1318      }
     1319
     1320      FaceIt(const Graph& graph, const Face& face)
     1321          : Face(face), _graph(&graph) {}
     1322
     1323      FaceIt& operator++() {
     1324        _graph->next(*this);
     1325        return *this;
     1326      }
     1327
     1328    };
     1329
     1330
     1331  /// Iterator class for the common faces of two nodes
     1332
     1333  /// This iterator goes through the common faces of two nodes
     1334  class CommonFacesIt : public Face {
     1335      const Graph* _graph;
     1336      const Node _n1, _n2;
     1337      std::set<Face> _f1, _f2;
     1338      typename std::set<Face>::const_iterator _fi2;
     1339    public:
     1340
     1341      CommonFacesIt() {}
     1342
     1343      CommonFacesIt(Invalid i) : Face(i) { }
     1344
     1345      explicit CommonFacesIt(const Graph& graph, Node n1, Node n2) : _graph(
     1346        &graph), _n1(n1), _n2(n2), _f1(), _f2(){
     1347        Arc _i1;
     1348        _graph->firstOut(_i1,_n1);
     1349        while (_i1 != INVALID) {
     1350          _f1.insert(_graph->leftFace(_i1));
     1351          _graph->nextOut(_i1);
     1352        }
     1353        Arc _i2;
     1354        _graph->firstOut(_i2,_n2);
     1355        while (_i2 != INVALID) {
     1356          _f2.insert(_graph->leftFace(_i2));
     1357          _graph->nextOut(_i2);
     1358        }
     1359
     1360        _fi2 = _f2.begin();
     1361        while (_fi2 != _f2.end() && _f1.count(*_fi2) == 0)
     1362          ++_fi2;
     1363        static_cast<Face&>(*this) = (_fi2 != _f2.end())?*_fi2:INVALID;
     1364      }
     1365
     1366      CommonFacesIt& operator++() {
     1367        ++_fi2;
     1368        while (_fi2 != _f2.end() && _f1.count(*_fi2) == 0)
     1369          ++_fi2;
     1370        static_cast<Face&>(*this) = (_fi2 != _f2.end())?*_fi2:INVALID;
     1371        return *this;
     1372      }
     1373
     1374    };
     1375
     1376  /// Iterator class for the common nodes of two facess
     1377
     1378  /// This iterator goes through the common nodes of two faces
     1379  class CommonNodesIt : public Node {
     1380      const Graph* _graph;
     1381      const Face _f1, _f2;
     1382      std::set<Node> _ns1,_ns2;
     1383      typename std::set<Node>::const_iterator _ni2;
     1384    public:
     1385
     1386      CommonNodesIt() {}
     1387
     1388      CommonNodesIt(Invalid i) : Face(i) { }
     1389
     1390      explicit CommonNodesIt(const Graph& graph, Face f1, Face f2) : _graph(
     1391        &graph), _f1(f1), _f2(f2), _ns1(), _ns2(){
     1392        Arc _i1;
     1393        _graph->firstCwF(_i1,_f1);
     1394        while (_i1 != INVALID) {
     1395          _ns1.insert(_graph->source(_i1));
     1396          _graph->nextCwF(_i1);
     1397        }
     1398        Arc _i2;
     1399        _graph->firstCwF(_i2,_f2);
     1400        while (_i2 != INVALID) {
     1401          _ns2.insert(_graph->source(_i2));
     1402          _graph->nextCwF(_i2);
     1403        }
     1404
     1405        _ni2 = _ns2.begin();
     1406        while (_ni2 != _ns2.end() && _ns1.count(*_ni2) == 0)
     1407          ++_ni2;
     1408        static_cast<Node&>(*this) = (_ni2 != _ns2.end())?*_ni2:INVALID;
     1409      }
     1410
     1411      CommonNodesIt& operator++() {
     1412        ++_ni2;
     1413        while (_ni2 != _ns2.end() && _ns1.count(*_ni2) == 0)
     1414          ++_ni2;
     1415        static_cast<Node&>(*this) = (_ni2 != _ns2.end())?*_ni2:INVALID;
     1416        return *this;
     1417      }
     1418
     1419    };
     1420
     1421  /// Iterator class for the arcs on the boundary of a face.
     1422
     1423  /// This iterator goes through the arcs on the boundary of a face clockwise.
     1424  class CwBoundaryArcIt : public Arc {
     1425      const Graph* _graph;
     1426      Face _face;
     1427      Arc f_arc;
     1428    public:
     1429
     1430      CwBoundaryArcIt() { }
     1431
     1432      CwBoundaryArcIt(Invalid i) : Arc(i) { }
     1433
     1434      CwBoundaryArcIt(const Graph& graph, const Face& face)
     1435          : _graph(&graph), _face(face) {
     1436        _graph->firstCwF(static_cast<Arc&>(*this),face);
     1437        f_arc = *this;
     1438      }
     1439
     1440      CwBoundaryArcIt(const Graph& graph, const Arc& arc)
     1441          : Arc(arc), _graph(&graph) {}
     1442
     1443      CwBoundaryArcIt& operator++() {
     1444        _graph->nextCwF(*this);
     1445        return *this;
     1446      }
     1447
     1448    };
     1449
     1450  /// Iterator class for the arcs around a node.
     1451
     1452  /// This iterator goes through the arcs around a node anticlockwise.
     1453  class CcwArcIt : public Arc {
     1454      const Graph* _graph;
     1455      const Node _node;
     1456    public:
     1457
     1458      CcwArcIt() { }
     1459
     1460      CcwArcIt(Invalid i) : Arc(i) { }
     1461
     1462      CcwArcIt(const Graph& graph, const Node& node)
     1463          : _graph(&graph), _node(node) {
     1464        _graph->firstCcw(*this, node);
     1465      }
     1466
     1467      CcwArcIt(const Graph& graph, const Arc& arc)
     1468          : Arc(arc), _graph(&graph) {}
     1469
     1470      CcwArcIt& operator++() {
     1471        _graph->nextCcw(*this, _node);
     1472        return *this;
     1473      }
     1474
     1475    };
     1476
     1477    void clear() {
     1478      notifier(Face()).clear();
     1479      Parent::clear();
     1480    }
     1481
     1482    template <typename Graph, typename NodeRefMap, typename EdgeRefMap>
     1483    void build(const Graph& graph, NodeRefMap& nodeRef,
     1484            EdgeRefMap& edgeRef) {
     1485      Parent::build(graph, nodeRef, edgeRef);
     1486      notifier(Face()).build();
     1487    }
     1488
     1489    PlanarGraphExtender() {
     1490      face_notifier.setContainer(*this);
     1491    }
     1492
     1493    ~PlanarGraphExtender() {
     1494      face_notifier.clear();
     1495    }
     1496
     1497  };
     1498
     1499  typedef PlanarGraphExtender<GraphExtender<PlanarGraphBase> >
     1500    ExtendedPlanarGraphBase;
     1501
     1502
     1503  /// \addtogroup graphs
     1504  /// @{
     1505
     1506  ///An undirected planar graph structure.
     1507
     1508  ///\ref PlanarGraph is a versatile and fast undirected planar graph
     1509  ///implementation based on linked lists that are stored in
     1510  ///\c std::vector structures. It maintains a topology of nodes, edges
     1511  ///and faces.
     1512  ///
     1513  ///This type fully conforms to the \ref concepts::Graph "Graph concept"
     1514  ///and it also provides several useful additional functionalities, including
     1515  ///those specific to planar graphs.
     1516  ///Most of its member functions and nested classes are documented
     1517  ///only in the concept class.
     1518  ///
     1519  ///This class provides only linear time counting for nodes, edges, arcs and
     1520  ///faces.
     1521  ///
     1522  ///A disconnected planar graph has have an outer face for each of its
     1523  ///components, effectively turning them into separate graphs. Each component
     1524  ///has a corresponding component in the dual.
     1525  ///\sa concepts::Graph
     1526  ///\sa PlaneGraph
     1527  class PlanarGraph : public ExtendedPlanarGraphBase {
     1528    typedef ExtendedPlanarGraphBase Parent;
     1529
     1530  public:
     1531    /// Graphs are \e not copy constructible. Use GraphCopy instead.
     1532    PlanarGraph(const PlanarGraph &) = delete;
     1533    /// \brief Assignment of a graph to another one is \e not allowed.
     1534    /// Use GraphCopy instead.
     1535    void operator=(const PlanarGraph &) = delete;
     1536
     1537    /// \brief Constructor
     1538
     1539    /// Constructor.
     1540    ///
     1541    PlanarGraph() {}
     1542
     1543
     1544    ///\brief Constructor that copies a planar embedding
     1545
     1546    ///Constructor that copies a planar embedding.
     1547    template<typename Graph>
     1548    PlanarGraph(const Graph &g, const lemon::PlanarEmbedding<Graph> &em) {
     1549      typename Graph::template NodeMap<Node> nm(g);
     1550      typename Graph::template ArcMap<Arc> am(g);
     1551      this->copyEmbedding(g,em,nm,am);
     1552    }
     1553
     1554    template<typename Graph>
     1555    PlanarGraph(const Graph &g, const lemon::PlanarEmbedding<Graph> &em,
     1556        typename Graph::template NodeMap<Node> &nm,
     1557        typename Graph::template ArcMap<Arc> &am) {
     1558      this->copyEmbedding(g,em,nm,am);
     1559    }
     1560
     1561    typedef Parent::OutArcIt IncEdgeIt;
     1562
     1563  protected:
     1564
     1565    template<typename Graph>
     1566    void copyEmbedding(const Graph &g, const lemon::PlanarEmbedding<Graph> &em,
     1567        typename Graph::template NodeMap<Node> &nm,
     1568        typename Graph::template ArcMap<Arc> &am) {
     1569      for (typename Graph::NodeIt it(g); it != INVALID; ++it) {
     1570        nm[it] = addNode();
     1571      }
     1572      for (typename Graph::ArcIt it(g); it != INVALID; ++it) {
     1573        am[it] = INVALID;
     1574      }
     1575      for (typename Graph::NodeIt it(g); it != INVALID; ++it) {
     1576        for (typename Graph::OutArcIt it2(g,it); it2 != INVALID; ++it2) {
     1577          if (am[it2] == INVALID) {
     1578            typename Graph::Arc p_u = em.next(it2);
     1579            while (am[p_u] == INVALID && it2 != p_u) {
     1580              p_u = em.next(p_u);
     1581            }
     1582            typename Graph::Arc ra = g.oppositeArc(it2);
     1583            typename Graph::Arc p_v = em.next(ra);
     1584            while (am[p_v] == INVALID && ra != p_v) {
     1585              p_v = em.next(p_v);
     1586            }
     1587            am[it2] = direct(addEdge(nm[g.source(it2)],nm[g.target(it2)],
     1588                    am[p_u],am[p_v]),nm[g.source(it2)]);
     1589            am[g.oppositeArc(it2)] = oppositeArc(am[it2]);
     1590          }
     1591        }
     1592      }
     1593    }
     1594
     1595    void edge_add_notify(Edge edge) {
     1596      notifier(Edge()).add(edge);
     1597      std::vector<Arc> ev;
     1598      ev.push_back(Parent::direct(edge, true));
     1599      ev.push_back(Parent::direct(edge, false));
     1600      notifier(Arc()).add(ev);
     1601    }
     1602
     1603    void edge_erase_notify(Edge edge) {
     1604      std::vector<Arc> av;
     1605      av.push_back(Parent::direct(edge, true));
     1606      av.push_back(Parent::direct(edge, false));
     1607      notifier(Arc()).erase(av);
     1608      notifier(Edge()).erase(edge);
     1609    }
     1610
     1611  public:
     1612    /// \brief Add a new node to the graph.
     1613    ///
     1614    /// This function adds a new node to the graph. A new outer face is created
     1615    /// for the node.
     1616    /// \return The new node.
     1617    Node addNode() {
     1618      Node n = PlanarGraphBase::addNode();
     1619      notifier(Node()).add(n);
     1620      notifier(Face()).add(outerFace(n));
     1621      return n;
     1622    }
     1623
     1624    /// \brief Add a new edge to the graph.
     1625    ///
     1626    /// This function adds a new edge to the graph between nodes
     1627    /// \c u and \c v with inherent orientation from node \c u to
     1628    /// node \c v. \c p_u and \c p_v are the edges that directly precede the new
     1629    /// edge in anticlockwise order at the nodes \c u and \c v, respectively.
     1630    /// INVALID should be passed as \c p_u or \c p_v if and only if the
     1631    /// respective node is isolated.
     1632    ///
     1633    /// If \c u and \c v are in the same component, \c p_u and \c p_v must share
     1634    /// the same left face (when looking from \c u or \c v). This face will be
     1635    /// split in two. If \c u and \c v are in different components, one of the
     1636    /// outer faces will be deleted.
     1637    /// \note It can be costly to join components unless one of them is an
     1638    /// isolated node.
     1639    /// \note Even if a face is not deleted in the process, BoundaryArcIt
     1640    /// iterators might run an incomplete lap or revisit previously passed edges
     1641    /// if edges are added to or removed from the boundary of the face.
     1642    /// \return The new edge, or INVALID if the parameters don't meet the
     1643    /// preconditions.
     1644    Edge addEdge(Node u, Node v, Edge p_u, Edge p_v) {
     1645      Face f_u = INVALID, f_v = INVALID;
     1646      if (p_u != INVALID) {
     1647        Arc a_u = direct(p_u,u);
     1648        f_u = leftFace(a_u);
     1649      }
     1650      if (p_v != INVALID) {
     1651        Arc a_v = direct(p_v,v);
     1652        f_v = leftFace(a_v);
     1653      }
     1654      Face o_u = outerFace(u);
     1655      Face o_v = outerFace(v);
     1656      Edge edge = PlanarGraphBase::addEdge(u, v, p_u, p_v);
     1657      if (edge == INVALID) return edge;
     1658      if (!valid(f_u)) notifier(Face()).erase(f_u);
     1659      if (!valid(f_v) && f_u != f_v) notifier(Face()).erase(f_v);
     1660      if (!valid(o_u)) notifier(Face()).erase(o_u);
     1661      if (!valid(o_v)) notifier(Face()).erase(o_v);
     1662      edge_add_notify(edge);
     1663      Face e_l = leftFace(direct(edge,u));
     1664      Face e_r = rightFace(direct(edge,u));
     1665      if (e_l != f_u && e_l != f_v && e_l != o_u && e_l != o_v)
     1666        notifier(Face()).add(e_l);
     1667      if (e_r != f_u && e_r != f_v && e_r != o_u && e_r != o_v && e_r != e_l)
     1668        notifier(Face()).add(e_r);
     1669      return edge;
     1670    }
     1671
     1672    ///\brief Erase a node from the graph.
     1673    ///
     1674    /// This function erases the given node along with its incident arcs
     1675    /// from the graph.
     1676    ///
     1677    /// \note All iterators referencing the removed node or the incident
     1678    /// edges are invalidated, of course.
     1679    void erase(Node n) {
     1680      notifier(Face()).erase(outerFace(n));
     1681      Parent::erase(n);
     1682    }
     1683
     1684    ///\brief Erase an edge from the graph.
     1685    ///
     1686    /// This function erases the given edge from the graph. The faces on the two
     1687    /// sides are merged into one, unless the edge holds two components
     1688    /// together. In the latter case, a new outer face is added and the
     1689    /// component is split in two.
     1690    ///
     1691    /// \note It can be costly to split a component, unless one of the resulting
     1692    /// components is an isolated node.
     1693    /// \note All iterators referencing the removed edge are invalidated,
     1694    /// of course.
     1695    /// \note Even if a face is not deleted in the process, BoundaryArcIt
     1696    /// iterators might run an incomplete lap or revisit previously passed edges
     1697    /// if edges are added to or removed from the boundary of the face.
     1698    void erase(Edge e) {
     1699      Node n1 = u(e);
     1700      Node n2 = v(e);
     1701      Face f_l = leftFace(direct(e,n1));
     1702      Face f_r = rightFace(direct(e,n1));
     1703      Parent::erase(e);
     1704      Face o_u = outerFace(n1);
     1705      Face o_v = outerFace(n2);
     1706      if (!valid(f_l)) notifier(Face()).erase(f_l);
     1707      if (!valid(f_r)) notifier(Face()).erase(f_r);
     1708      if (valid(o_u) && o_u != f_l && o_u != f_r)
     1709        notifier(Face()).add(o_u);
     1710      if (valid(o_v) && o_v != f_l && o_v != f_r)
     1711        notifier(Face()).add(o_v);
     1712    }
     1713    /// Node validity check
     1714
     1715    /// This function gives back \c true if the given node is valid,
     1716    /// i.e. it is a real node of the graph.
     1717    ///
     1718    /// \warning A removed node could become valid again if new nodes are
     1719    /// added to the graph.
     1720    bool valid(Node n) const {
     1721      return Parent::valid(n);
     1722    }
     1723    /// Edge validity check
     1724
     1725    /// This function gives back \c true if the given edge is valid,
     1726    /// i.e. it is a real edge of the graph.
     1727    ///
     1728    /// \warning A removed edge could become valid again if new edges are
     1729    /// added to the graph.
     1730    bool valid(Edge e) const {
     1731      return Parent::valid(e);
     1732    }
     1733    /// Arc validity check
     1734
     1735    /// This function gives back \c true if the given arc is valid,
     1736    /// i.e. it is a real arc of the graph.
     1737    ///
     1738    /// \warning A removed arc could become valid again if new edges are
     1739    /// added to the graph.
     1740    bool valid(Arc a) const {
     1741      return Parent::valid(a);
     1742    }
     1743
     1744    /// Face validity check
     1745
     1746    /// This function gives back \c true if the given face is valid,
     1747    /// i.e. it is a real face of the graph.
     1748    ///
     1749    /// \warning A removed face could become valid again if new faces are
     1750    /// added to the graph.
     1751    bool valid(Face f) const {
     1752      return Parent::valid(f);
     1753    }
     1754
     1755    /// \brief Change the first node of an edge.
     1756    ///
     1757    /// Planar graphs don't support changing the endpoints of edges.
     1758    void changeU(Edge e, Node n) = delete;
     1759    /// \brief Change the second node of an edge.
     1760    ///
     1761    /// Planar graphs don't support changing the endpoints of edges.
     1762    void changeV(Edge e, Node n) = delete;
     1763
     1764    /// \brief Contract two nodes.
     1765    ///
     1766    /// This function contracts the two ends of the given edge.
     1767    /// One of the nodes on \c e is removed, but instead of deleting
     1768    /// its incident edges, they are joined to the other node.
     1769    /// \c e will be deleted.
     1770    ///
     1771    /// The deleted node is v(e) unless the degree of v(e) is higher than 1 and
     1772    /// the degree of u(e) is 1, in which case u(e) is deleted.
     1773    ///
     1774    /// \note All edge and arc iterators whose base node is
     1775    /// the deleted one are invalidated.
     1776    /// Moreover all iterators referencing the removed node or
     1777    /// edge are also invalidated. Other iterators remain valid. However, they
     1778    /// might run an incomplete lap or revisit previously passed edges if
     1779    /// incremented.
     1780    ///
     1781    ///\warning This functionality cannot be used together with the
     1782    ///Snapshot feature.
     1783    void contract(Edge e) {
     1784      Node n1 = u(e);
     1785      Node n2 = v(e);
     1786      if (n1 == n2) return;
     1787      Node nd;
     1788      bool simple;
     1789      inner_contract1(n1,n2,nd,simple);
     1790      notifier(Node()).erase(nd);
     1791      edge_erase_notify(e);
     1792      inner_contract2(e,n1,n2,nd,simple);
     1793    }
     1794
     1795    /// \brief Split an edge.
     1796    ///
     1797    ///This function splits the given edge. First, a new node \c v is
     1798    ///added to the graph, then the target node of the original edge
     1799    ///is set to \c v. Finally, an edge from \c v to the original target
     1800    ///is added.
     1801    ///\return The newly created node.
     1802    ///
     1803    ///\note Iterators referencing the original edge are
     1804    ///invalidated. Other iterators remain valid.
     1805    ///
     1806    ///\warning This functionality cannot be used together with the
     1807    ///Snapshot feature.
     1808    Node split(Edge e) {
     1809      Node v;
     1810      Edge e2;
     1811      inner_split1(v,e2);
     1812      notifier(Node()).add(v);
     1813      edge_add_notify(e2);
     1814      return inner_split2(e,v,e2);
     1815    }
     1816
     1817    ///Split a node.
     1818
     1819    ///This function splits the given node. First, a new node is added
     1820    ///to the graph, then all edges in anticlockwise order from \c e1 until but
     1821    ///not including \c e2
     1822    ///are re-anchored from \c n to the new node.
     1823    ///If the second parameter \c connect is \c true (this is the default
     1824    ///value), then a new edge from node \c n to the newly created node
     1825    ///is also added. \c e1 and \c e2 must be INVALID if and only if \c n is
     1826    ///isolated. If \c connect is false, the component of \c n might be split in
     1827    ///two. It can be costly to split components unless one of the resulting
     1828    ///components is an isolated node.
     1829    ///\return The newly created node.
     1830    ///
     1831    ///\note All iterators remain valid but some might run an incomplete lap or
     1832    ///revisit previously passed nodes when incremented.
     1833    ///
     1834    ///\warning This functionality cannot be used together with the
     1835    ///Snapshot feature.
     1836    Node split(Node n1, Edge e1, Edge e2, bool connect) {
     1837      Node n2;
     1838      inner_split1(n1,n2);
     1839      notifier(Node()).add(n2);
     1840      Edge a3;
     1841      inner_split2(n1,n2,e1,e2,connect,a3);
     1842      if (!connect) {
     1843        if (a3 != INVALID)
     1844          erase(a3);
     1845      } else {
     1846        edge_add_notify(a3);
     1847      }
     1848      return n2;
     1849    }
     1850
     1851
     1852    ///Clear the graph.
     1853
     1854    ///This function erases all nodes and arcs from the graph.
     1855    ///
     1856    ///\note All iterators of the graph are invalidated, of course.
     1857    void clear() {
     1858      Parent::clear();
     1859    }
     1860
     1861    /// Reserve memory for nodes.
     1862
     1863    /// Using this function, it is possible to avoid superfluous memory
     1864    /// allocation: if you know that the graph you want to build will
     1865    /// be large (e.g. it will contain millions of nodes and/or edges),
     1866    /// then it is worth reserving space for this amount before starting
     1867    /// to build the graph.
     1868    /// \sa reserveEdge()
     1869    void reserveNode(int n) {
     1870      nodes.reserve(n);
     1871    };
     1872
     1873    /// Reserve memory for edges.
     1874
     1875    /// Using this function, it is possible to avoid superfluous memory
     1876    /// allocation: if you know that the graph you want to build will
     1877    /// be large (e.g. it will contain millions of nodes and/or edges),
     1878    /// then it is worth reserving space for this amount before starting
     1879    /// to build the graph.
     1880    /// \sa reserveNode()
     1881    void reserveEdge(int m) {
     1882      arcs.reserve(2 * m);
     1883    };
     1884
     1885    /// Reserve memory for faces.
     1886
     1887    /// Using this function, it is possible to avoid superfluous memory
     1888    /// allocation: if you know that the graph you want to build will
     1889    /// be large (e.g. it will contain millions of nodes and/or edges),
     1890    /// then it is worth reserving space for this amount before starting
     1891    /// to build the graph.
     1892    /// \sa reserveFace()
     1893    void reserveFace(int n) {
     1894      faces.reserve(n);
     1895    };
     1896
     1897    class DualBase {
     1898      const Graph *_graph;
     1899    protected:
     1900      void initialize(const Graph *graph) {
     1901        _graph = graph;
     1902      }
     1903    public:
     1904
     1905      typedef PlanarGraph::Face Node;
     1906      typedef PlanarGraph::Arc Arc;
     1907      typedef PlanarGraph::Edge Edge;
     1908      typedef PlanarGraph::Node Face;
     1909
     1910      int maxNodeId() const {
     1911        return _graph->maxFaceId();
     1912      }
     1913      int maxArcId() const {
     1914        return _graph->maxArcId();
     1915      }
     1916      int maxFaceId() const {
     1917        return _graph->maxNodeId();
     1918      }
     1919
     1920      Node source(Arc e) const {
     1921        return _graph->leftFace(e);
     1922      }
     1923      Node target(Arc e) const {
     1924        return _graph->rightFace(e);
     1925      }
     1926      Face leftFace(Arc e) const {
     1927        return _graph->target(e);
     1928      }
     1929      Face rightFace(Arc e) const {
     1930        return _graph->source(e);
     1931      }
     1932      Arc direct(const Edge &edge, const Node &node) const {
     1933        return _graph->direct(edge, _graph->w1(edge) == node);
     1934      }
     1935
     1936      void first(Node &i) const {
     1937        _graph->first(i);
     1938      }
     1939      void next(Node &i) const {
     1940        _graph->next(i);
     1941      }
     1942      void first(Arc &i) const {
     1943        _graph->first(i);
     1944      }
     1945      void next(Arc &i) const {
     1946        _graph->next(i);
     1947      }
     1948      void firstCcw(Arc& i, const Node& n) const {
     1949        _graph->lastInF(i, n);
     1950      }
     1951      void nextCcw(Arc& i, const Node &n) const {
     1952        _graph->prevInF(i);
     1953      }
     1954      void firstIn(Arc& i, const Node& n) const {
     1955        _graph->firstInF(i, n);
     1956      }
     1957      void nextIn(Arc& i) const {
     1958        _graph->nextInF(i);
     1959      }
     1960      void firstCwF(Arc& i, const Face& n) const {
     1961        _graph->lastIn(i, n);
     1962      }
     1963      void nextCwF(Arc& i) const {
     1964        _graph->prevIn(i);
     1965      }
     1966      void firstOut(Arc& i, const Node& n ) const {
     1967        _graph->firstOutF(i, n);
     1968      }
     1969      void nextOut(Arc& i) const {
     1970        _graph->nextOutF(i);
     1971      }
     1972      void first(Face &i) const {
     1973        _graph->first(i);
     1974      }
     1975      void next(Face &i) const {
     1976        _graph->next(i);
     1977      }
     1978
     1979      static int id(Node v) {
     1980        return PlanarGraph::id(v);
     1981      }
     1982      static int id(Arc e) {
     1983        return PlanarGraph::id(e);
     1984      }
     1985      static int id(Face f) {
     1986        return PlanarGraph::id(f);
     1987      }
     1988      static Node nodeFromId(int id) {
     1989        return PlanarGraph::faceFromId(id);
     1990      }
     1991      static Arc arcFromId(int id) {
     1992        return PlanarGraph::arcFromId(id);
     1993      }
     1994      static Face faceFromId(int id) {
     1995        return PlanarGraph::nodeFromId(id);
     1996      }
     1997
     1998      bool valid(Node n) const {
     1999        return _graph->valid(n);
     2000      }
     2001      bool valid(Arc n) const {
     2002        return _graph->valid(n);
     2003      }
     2004      bool valid(Face n) const {
     2005        return _graph->valid(n);
     2006      }
     2007
     2008    };
     2009
     2010    typedef PlanarGraphExtender<GraphExtender<DualBase> > ExtendedDualBase;
     2011
     2012  /// Adaptor class for the dual of a planar graph.
     2013
     2014  /// This is an adaptor class for the dual of a planar graph.
     2015  class Dual : public ExtendedDualBase {
     2016    public:
     2017      Dual(const PlanarGraph &graph) {
     2018        initialize(&graph);
     2019      }
     2020
     2021    };
     2022    /// \brief Class to make a snapshot of the graph and restore
     2023    /// it later.
     2024    ///
     2025    /// Class to make a snapshot of the graph and restore it later.
     2026    ///
     2027    /// The newly added nodes and edges can be removed
     2028    /// using the restore() function.
     2029    ///
     2030    /// \note After a state is restored, you cannot restore a later state,
     2031    /// i.e. you cannot add the removed nodes and edges again using
     2032    /// another Snapshot instance.
     2033    ///
     2034    /// \warning Node and edge deletions and other modifications
     2035    /// (e.g. changing the end-nodes of edges or contracting nodes)
     2036    /// cannot be restored. These events invalidate the snapshot.
     2037    /// However, the edges and nodes that were added to the graph after
     2038    /// making the current snapshot can be removed without invalidating it.
     2039    class Snapshot {
     2040    protected:
     2041
     2042      typedef Parent::NodeNotifier NodeNotifier;
     2043
     2044    class NodeObserverProxy : public NodeNotifier::ObserverBase {
     2045      public:
     2046
     2047        NodeObserverProxy(Snapshot& _snapshot)
     2048            : snapshot(_snapshot) {}
     2049
     2050        using NodeNotifier::ObserverBase::attach;
     2051        using NodeNotifier::ObserverBase::detach;
     2052        using NodeNotifier::ObserverBase::attached;
     2053
     2054      protected:
     2055
     2056        virtual void add(const Node& node) {
     2057          snapshot.addNode(node);
     2058        }
     2059        virtual void add(const std::vector<Node>& nodes) {
     2060          for (int i = nodes.size() - 1; i >= 0; ++i) {
     2061            snapshot.addNode(nodes[i]);
     2062          }
     2063        }
     2064        virtual void erase(const Node& node) {
     2065          snapshot.eraseNode(node);
     2066        }
     2067        virtual void erase(const std::vector<Node>& nodes) {
     2068          for (int i = 0; i < int(nodes.size()); ++i) {
     2069            snapshot.eraseNode(nodes[i]);
     2070          }
     2071        }
     2072        virtual void build() {
     2073          Node node;
     2074          std::vector<Node> nodes;
     2075          for (notifier()->first(node); node != INVALID;
     2076               notifier()->next(node)) {
     2077            nodes.push_back(node);
     2078          }
     2079          for (int i = nodes.size() - 1; i >= 0; --i) {
     2080            snapshot.addNode(nodes[i]);
     2081          }
     2082        }
     2083        virtual void clear() {
     2084          Node node;
     2085          for (notifier()->first(node); node != INVALID;
     2086               notifier()->next(node)) {
     2087            snapshot.eraseNode(node);
     2088          }
     2089        }
     2090
     2091        Snapshot& snapshot;
     2092      };
     2093
     2094    class EdgeObserverProxy : public EdgeNotifier::ObserverBase {
     2095      public:
     2096
     2097        EdgeObserverProxy(Snapshot& _snapshot)
     2098            : snapshot(_snapshot) {}
     2099
     2100        using EdgeNotifier::ObserverBase::attach;
     2101        using EdgeNotifier::ObserverBase::detach;
     2102        using EdgeNotifier::ObserverBase::attached;
     2103
     2104      protected:
     2105
     2106        virtual void add(const Edge& edge) {
     2107          snapshot.addEdge(edge);
     2108        }
     2109        virtual void add(const std::vector<Edge>& edges) {
     2110          for (int i = edges.size() - 1; i >= 0; ++i) {
     2111            snapshot.addEdge(edges[i]);
     2112          }
     2113        }
     2114        virtual void erase(const Edge& edge) {
     2115          snapshot.eraseEdge(edge);
     2116        }
     2117        virtual void erase(const std::vector<Edge>& edges) {
     2118          for (int i = 0; i < int(edges.size()); ++i) {
     2119            snapshot.eraseEdge(edges[i]);
     2120          }
     2121        }
     2122        virtual void build() {
     2123          Edge edge;
     2124          std::vector<Edge> edges;
     2125          for (notifier()->first(edge); edge != INVALID;
     2126               notifier()->next(edge)) {
     2127            edges.push_back(edge);
     2128          }
     2129          for (int i = edges.size() - 1; i >= 0; --i) {
     2130            snapshot.addEdge(edges[i]);
     2131          }
     2132        }
     2133        virtual void clear() {
     2134          Edge edge;
     2135          for (notifier()->first(edge); edge != INVALID;
     2136               notifier()->next(edge)) {
     2137            snapshot.eraseEdge(edge);
     2138          }
     2139        }
     2140
     2141        Snapshot& snapshot;
     2142      };
     2143
     2144      PlanarGraph *graph;
     2145
     2146      NodeObserverProxy node_observer_proxy;
     2147      EdgeObserverProxy edge_observer_proxy;
     2148
     2149      std::list<Node> added_nodes;
     2150      std::list<Edge> added_edges;
     2151
     2152
     2153      void addNode(const Node& node) {
     2154        added_nodes.push_front(node);
     2155      }
     2156      void eraseNode(const Node& node) {
     2157        std::list<Node>::iterator it =
     2158          std::find(added_nodes.begin(), added_nodes.end(), node);
     2159        if (it == added_nodes.end()) {
     2160          clear();
     2161          edge_observer_proxy.detach();
     2162          throw NodeNotifier::ImmediateDetach();
     2163        } else {
     2164          added_nodes.erase(it);
     2165        }
     2166      }
     2167
     2168      void addEdge(const Edge& edge) {
     2169        added_edges.push_front(edge);
     2170      }
     2171      void eraseEdge(const Edge& edge) {
     2172        std::list<Edge>::iterator it =
     2173          std::find(added_edges.begin(), added_edges.end(), edge);
     2174        if (it == added_edges.end()) {
     2175          clear();
     2176          node_observer_proxy.detach();
     2177          throw EdgeNotifier::ImmediateDetach();
     2178        } else {
     2179          added_edges.erase(it);
     2180        }
     2181      }
     2182
     2183      void attach(PlanarGraph &_graph) {
     2184        graph = &_graph;
     2185        node_observer_proxy.attach(graph->notifier(Node()));
     2186        edge_observer_proxy.attach(graph->notifier(Edge()));
     2187      }
     2188
     2189      void detach() {
     2190        node_observer_proxy.detach();
     2191        edge_observer_proxy.detach();
     2192      }
     2193
     2194      bool attached() const {
     2195        return node_observer_proxy.attached();
     2196      }
     2197
     2198      void clear() {
     2199        added_nodes.clear();
     2200        added_edges.clear();
     2201      }
     2202
     2203    public:
     2204
     2205      /// \brief Default constructor.
     2206      ///
     2207      /// Default constructor.
     2208      /// You have to call save() to actually make a snapshot.
     2209      Snapshot()
     2210          : graph(0), node_observer_proxy(*this),
     2211          edge_observer_proxy(*this) {}
     2212
     2213      /// \brief Constructor that immediately makes a snapshot.
     2214      ///
     2215      /// This constructor immediately makes a snapshot of the given graph.
     2216      Snapshot(PlanarGraph &gr)
     2217          : node_observer_proxy(*this),
     2218          edge_observer_proxy(*this) {
     2219        attach(gr);
     2220      }
     2221
     2222      /// \brief Make a snapshot.
     2223      ///
     2224      /// This function makes a snapshot of the given graph.
     2225      /// It can be called more than once. In case of a repeated
     2226      /// call, the previous snapshot gets lost.
     2227      void save(PlanarGraph &gr) {
     2228        if (attached()) {
     2229          detach();
     2230          clear();
     2231        }
     2232        attach(gr);
     2233      }
     2234
     2235      /// \brief Undo the changes until the last snapshot.
     2236      ///
     2237      /// This function undos the changes until the last snapshot
     2238      /// created by save() or Snapshot(PlanarGraph&).
     2239      ///
     2240      /// \warning This method invalidates the snapshot, i.e. repeated
     2241      /// restoring is not supported unless you call save() again.
     2242      void restore() {
     2243        detach();
     2244        for (std::list<Edge>::iterator it = added_edges.begin();
     2245             it != added_edges.end(); ++it) {
     2246          graph->erase(*it);
     2247        }
     2248        for (std::list<Node>::iterator it = added_nodes.begin();
     2249             it != added_nodes.end(); ++it) {
     2250          graph->erase(*it);
     2251        }
     2252        clear();
     2253      }
     2254
     2255      /// \brief Returns \c true if the snapshot is valid.
     2256      ///
     2257      /// This function returns \c true if the snapshot is valid.
     2258      bool valid() const {
     2259        return attached();
     2260      }
     2261    };
     2262  };
     2263
     2264  /// @}
     2265} //namespace lemon
     2266
     2267
     2268#endif
  • new file lemon/plane_graph.h

    diff -r 0bca98cbebbb -r ea56aa8243b1 lemon/plane_graph.h
    - +  
     1/* -*- mode: C++; indent-tabs-mode: nil; -*-
     2 *
     3 * This file is a part of LEMON, a generic C++ optimization library.
     4 *
     5 * Copyright (C) 2003-2010
     6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8 *
     9 * Permission to use, modify and distribute this software is granted
     10 * provided that this copyright notice appears in all copies. For
     11 * precise terms see the accompanying LICENSE file.
     12 *
     13 * This software is provided "AS IS" with no warranty of any kind,
     14 * express or implied, and with no claim as to its suitability for any
     15 * purpose.
     16 *
     17 */
     18
     19#ifndef LEMON_PLANE_GRAPH_H
     20#define LEMON_PLANE_GRAPH_H
     21
     22///\ingroup graphs
     23///\file
     24///\brief PlaneGraph classes. UNDER CONSTRUCTION.
     25
     26#include <lemon/planarity.h>
     27#include <lemon/planar_graph.h>
     28#include <lemon/dim2.h>
     29#include <cmath>
     30#include <lemon/math.h>
     31#include <map>
     32
     33#include "bits\graph_extender.h"
     34#include "bits\bezier.h"
     35
     36
     37namespace lemon {
     38
     39  template<typename Point>
     40  struct LexicographicPointOrdering {
     41    bool operator()(const Point &p1, const Point &p2) {
     42      return p1.x < p2.x || (p1.x == p2.x && p1.y < p2.y);
     43    }
     44  };
     45
     46  ///Traits class for straight edges in a plane graph.
     47
     48  ///This class represents straight edges. It can be used as a template
     49  ///parameter to \ref PlaneGraph.
     50  class PlaneGraphStraightEdge {
     51
     52  private:
     53    struct Data {
     54      dim2::Point<double> p1, p2;
     55      double angle;
     56    };
     57
     58    static const double epsilon = 0.00000000000001;
     59
     60    static void getEquation(const Data &d, double &a, double &b, double &c) {
     61      double vx = d.p2.x-d.p1.x, vy = d.p2.y-d.p1.y;
     62      a = vy;
     63      b = -vx;
     64      c = vy*d.p1.x-vx*d.p1.y;
     65    }
     66
     67    static bool hasPoint(const Data &d, dim2::Point<double> p, bool os, bool oe
     68      ) {
     69      if (!os && d.p1.x-epsilon < p.x && p.x < d.p1.x+epsilon &&
     70        d.p1.y-epsilon < p.y && p.y < d.p1.y+epsilon) return true;
     71      if (!oe && d.p2.x-epsilon < p.x && p.x < d.p2.x+epsilon &&
     72        d.p2.y-epsilon < p.y && p.y < d.p2.y+epsilon) return true;
     73      double x1,x2,y1,y2;
     74      d.p1.x <= d.p2.x ? (x1 = d.p1.x, x2 = d.p2.x) : (x2 = d.p1.x, x1 =
     75        d.p2.x);
     76      d.p1.y <= d.p2.y ? (y1 = d.p1.y, y2 = d.p2.y) : (y2 = d.p1.y, y1 =
     77        d.p2.y);
     78      return ((x1+epsilon < p.x && p.x < x2-epsilon) ||
     79             (x1-epsilon <= p.x && p.x <= x1+epsilon &&
     80              x2-epsilon <= p.x && p.x <= x2+epsilon)) &&
     81             ((y1+epsilon < p.y && p.y < y2-epsilon) ||
     82             (y1-epsilon <= p.y && p.y <= y1+epsilon &&
     83              y2-epsilon <= p.y && p.y <= y2+epsilon));
     84    }
     85
     86  public:
     87
     88    typedef Data Shape;
     89    typedef dim2::Point<double> Point;
     90
     91    static Data create() { return Data(); }
     92    static Data create(dim2::Point<double> p1_, dim2::Point<double> p2_) {
     93      Data d = {p1_, p2_, atan2(p2_.y-p1_.y,p2_.x-p1_.x)};
     94      return d;
     95    }
     96
     97    static Point source(const Data &d) {
     98      return d.p1;
     99    }
     100
     101    static Point target(const Data &d) {
     102      return d.p2;
     103    }
     104
     105    //Determine if c comes between a and b or not, in positive order
     106    static bool between(const Data &c, const Data &a, const Data &b)
     107    {
     108      double aa = a.angle, ba = b.angle, ca = c.angle;
     109      if (ba < aa) ba += 2*PI;
     110      if (ca < aa) ca += 2*PI;
     111      return ba > ca;
     112    }
     113
     114    static Data reverse(const Data &d) {
     115      return create(d.p2, d.p1);
     116    }
     117
     118    static bool match(const Data &a, const Data &b) {
     119      return (a.p1 == b.p1 && a.p2 == b.p2) || (a.p1 == b.p2 && a.p2 == b.p1);
     120    }
     121
     122    static bool intersect(const Data &q, const Data &e, bool qos = true,
     123      bool qoe = true, bool eos = true, bool eoe = true) {
     124      if (match(q,e)) return true;
     125      double a1,b1,c1,a2,b2,c2;
     126      getEquation(q,a1,b1,c1);
     127      getEquation(e,a2,b2,c2);
     128      double d = b1*a2-b2*a1;
     129      if (d == 0) {
     130        if (a2*q.p1.x+b2*q.p1.y == c2) {
     131          return hasPoint(q,e.p1,qos,qoe) || hasPoint(q,e.p2,qos,qoe) ||
     132            hasPoint(e,q.p1,eos,eoe) || hasPoint(e,q.p2,eos,eoe);
     133        } else {
     134          return false;
     135        }
     136      }
     137      double x = (c2*b1-c1*b2)/(d);
     138      double y = (c1*a2-c2*a1)/(d);
     139      dim2::Point<double> pt(x,y);
     140      return hasPoint(q,pt,qos,qoe) && hasPoint(e,pt,eos,eoe);
     141    }
     142
     143    static bool degenerate(const Data &d) {
     144      return d.p1 == d.p2;
     145    }
     146
     147    static void toEps(std::ostream &os, const Data &d, const Color &col,
     148      double width) {
     149      os << d.p1.x << ' '
     150         << d.p1.y << ' '
     151         << d.p2.x << ' '
     152         << d.p2.y << ' '
     153         << col.red() << ' '
     154         << col.green() << ' '
     155         << col.blue() << ' '
     156         << width << " l\n";
     157    }
     158  };
     159
     160  ///Class for a polyline edge in a plane graph.
     161
     162  ///This class represents a polyline edge. It can be used as a template
     163  ///parameter to \ref PlaneGraph.
     164  class PlaneGraphLineStripEdge {
     165  private:
     166    struct Data {
     167      std::vector<dim2::Point<double> > ps;
     168      double angle;
     169    };
     170
     171  public:
     172
     173    typedef Data Shape;
     174    typedef dim2::Point<double> Point;
     175
     176    static Data create() { return Data(); }
     177    static Data create(const std::vector<Point> &ps_) {
     178      Data d;
     179      d.ps = ps_;
     180      if (ps_.size() > 1) {
     181        const Point &p1 = ps_[0], &p2 = ps_[1];
     182        d.angle = atan2(p2.y-p1.y,p2.x-p1.x);
     183      }
     184      return d;
     185    }
     186
     187    static Point source(const Data &d) {
     188      return d.ps.front();
     189    }
     190
     191    static Point target(const Data &d) {
     192      return d.ps.back();
     193    }
     194
     195    //Determine if c comes between a and b or not, in positive order
     196    static bool between(const Data &c, const Data &a, const Data &b)
     197    {
     198      double aa = a.angle, ba = b.angle, ca = c.angle;
     199      if (ba < aa) ba += 2*PI;
     200      if (ca < aa) ca += 2*PI;
     201      return ba > ca;
     202    }
     203
     204    static Data reverse(const Data &d) {
     205      std::vector<Point> e(d.ps.size());
     206      reverse_copy(d.ps.begin(), d.ps.end(), e.begin());
     207      return create(e);
     208    }
     209
     210    static bool intersect(const Data &q, const Data &e) {
     211      for (int i=0; i<q.ps.size()-1; ++i) {
     212        PlaneGraphStraightEdge::Shape s1 { q.ps[i], q.ps[i+1] };
     213        for (int j=0; j<e.ps.size()-1; ++j) {
     214          PlaneGraphStraightEdge::Shape s2 { e.ps[j], e.ps[j+1] };
     215          if (PlaneGraphStraightEdge::intersect(s1,s2,i==0,i==q.ps.size()-2,
     216            j==0,j==e.ps.size()-2)) return true;
     217        }
     218      }
     219      return false;
     220    }
     221
     222    static bool degenerate(const Data &d) {
     223      if (d.ps.size() < 2) return true;
     224      for (int i=0; i<d.ps.size(); ++i)
     225        for (int j=i+1; j<d.ps.size(); ++j) {
     226          if (d.ps[i] == d.ps[j] && (i > 0 || j<d.ps.size()-1)) return true;
     227        }
     228
     229      for (int i=0; i<d.ps.size()-1; ++i) {
     230        PlaneGraphStraightEdge::Shape s1 { d.ps[i], d.ps[i+1] };
     231        for (int j=i+1; j<d.ps.size()-1; ++j) {
     232          PlaneGraphStraightEdge::Shape s2 { d.ps[j], d.ps[j+1] };
     233          if (PlaneGraphStraightEdge::intersect(s1,s2,true,true,true,true))
     234            return true;
     235        }
     236      }
     237
     238      return false;
     239    }
     240
     241    static void toEps(std::ostream &os, const Data &d, const Color &col,
     242      double width) {
     243      for (int i=0; i<d.ps.size()-1; ++i) {
     244        os << d.ps[i].x << ' '
     245           << d.ps[i].y << ' '
     246           << d.ps[i+1].x << ' '
     247           << d.ps[i+1].y << ' '
     248           << col.red() << ' '
     249           << col.green() << ' '
     250           << col.blue() << ' '
     251           << width << " l\n";
     252      }
     253  };
     254
     255  };
     256
     257  ///Class for a B&eacute;zier curve edge in a plane graph.
     258
     259  ///This class represents a B&eacute;zier curve edge of degree 2. It can be
     260  ///used as a
     261  ///template parameter to \ref PlaneGraph.
     262  class PlaneGraphBezier2Edge {
     263
     264  private:
     265    struct Data {
     266      dim2::Bezier2 ps;
     267      double angle2, angle1;
     268    };
     269
     270  public:
     271
     272    typedef Data Shape;
     273    typedef dim2::Point<double> Point;
     274
     275    static Data create() { return Data(); }
     276    static Data create(const Point &p1, const Point &p2, const Point &p3) {
     277      Data d;
     278      d.ps = dim2::Bezier2(p1,p2,p3);
     279      d.angle2 = atan2(p2.y-p1.y,p2.x-p1.x);
     280      dim2::Bezier1 pps = d.ps.grad();
     281      d.angle1 = atan2(pps.p2.y-pps.p1.y,pps.p2.x-pps.p1.x);
     282      return d;
     283    }
     284
     285    static Point source(const Data &d) {
     286      return d.ps.p1;
     287    }
     288
     289    static Point target(const Data &d) {
     290      return d.ps.p3;
     291    }
     292
     293    //Determine if c comes between a and b or not, in positive order
     294    static bool between(const Data &c, const Data &a, const Data &b)
     295    {
     296      double aa = a.angle2, ba = b.angle2, ca = c.angle2;
     297      if (ba < aa) ba += 2*PI;
     298      if (ca < aa) ca += 2*PI;
     299      if (ba > ca) return true;
     300      else if (ba < ca) return false;
     301      else {
     302        double aa = a.angle1, ba = b.angle1, ca = c.angle1;
     303        if (ba < aa) ba += 2*PI;
     304        if (ca < aa) ca += 2*PI;
     305        if (ba > ca) return true;
     306        return ba > ca;
     307      }
     308    }
     309
     310    static Data reverse(const Data &d) {
     311      return create(d.ps.p3,d.ps.p2,d.ps.p1);
     312    }
     313
     314    static bool intersect(const Data &q, const Data &e) {
     315      dim2::Box<double> box1 = q.ps.boundingBox();
     316      dim2::Box<double> box2 = q.ps.boundingBox();
     317      if ((box1 & box2).empty()) return false;
     318      std::vector<Point> v1, v2;
     319      for (int i=0; i<=100; ++i) {
     320        v1.push_back(q.ps(1.0*i/100));
     321        v2.push_back(e.ps(1.0*i/100));
     322      }
     323      return PlaneGraphLineStripEdge::intersect(
     324              PlaneGraphLineStripEdge::create(v1),
     325              PlaneGraphLineStripEdge::create(v2));
     326    }
     327
     328    static bool degenerate(const Data &d) {
     329      return d.ps.p1 == d.ps.p2 || d.ps.p2 == d.ps.p3 || d.ps.p1 == d.ps.p3;
     330    }
     331
     332    static void toEps(std::ostream &os, const Data &d, const Color &col,
     333      double width) {
     334      dim2::Bezier3 b = d.ps;
     335      os << width << " setlinewidth "
     336         << col.red() << ' '
     337         << col.green() << ' '
     338         << col.blue() << " setrgbcolor newpath\n"
     339         << b.p1.x << ' ' << b.p1.y << " moveto\n"
     340         << b.p2.x << ' ' << b.p2.y << ' '
     341         << b.p3.x << ' ' << b.p3.y << ' '
     342         << b.p4.x << ' ' << b.p4.y << " curveto stroke\n";
     343    }
     344  };
     345
     346  ///Class for a B&eacute;zier curve edge in a plane graph.
     347
     348  ///This class represents a B&eacute;zier curve edge of degree 3. It can be
     349  ///used as a
     350  ///template parameter to \ref PlaneGraph.
     351  class PlaneGraphBezier3Edge {
     352
     353  private:
     354    struct Data {
     355      dim2::Bezier3 ps;
     356      double angle3, angle2, angle1;
     357    };
     358
     359  public:
     360
     361    typedef Data Shape;
     362    typedef dim2::Point<double> Point;
     363
     364    static Data create() { return Data(); }
     365    static Data create(const Point &p1, const Point &p2, const Point &p3,
     366        const Point &p4) {
     367      Data d;
     368      d.ps = dim2::Bezier3(p1,p2,p3,p4);
     369      d.angle3 = atan2(p2.y-p1.y,p2.x-p1.x);
     370      dim2::Bezier2 pps = d.ps.grad();
     371      d.angle2 = atan2(pps.p2.y-pps.p1.y,pps.p2.x-pps.p1.x);
     372      dim2::Bezier1 ppps = pps.grad();
     373      d.angle3 = atan2(ppps.p2.y-ppps.p1.y,ppps.p2.x-ppps.p1.x);
     374      return d;
     375    }
     376
     377    static Point source(const Data &d) {
     378      return d.ps.p1;
     379    }
     380
     381    static Point target(const Data &d) {
     382      return d.ps.p4;
     383    }
     384
     385    //Determine if c comes between a and b or not, in positive order
     386    static bool between(const Data &c, const Data &a, const Data &b)
     387    {
     388      double aa = a.angle3, ba = b.angle3, ca = c.angle3;
     389      if (ba < aa) ba += 2*PI;
     390      if (ca < aa) ca += 2*PI;
     391      if (ba > ca) return true;
     392      else if (ba < ca) return false;
     393      else {
     394        double aa = a.angle2, ba = b.angle2, ca = c.angle2;
     395        if (ba < aa) ba += 2*PI;
     396        if (ca < aa) ca += 2*PI;
     397        if (ba > ca) return true;
     398        else if (ba < ca) return false;
     399        else {
     400          double aa = a.angle1, ba = b.angle1, ca = c.angle1;
     401          if (ba < aa) ba += 2*PI;
     402          if (ca < aa) ca += 2*PI;
     403          return (ba > ca);
     404        }
     405      }
     406    }
     407
     408    static Data reverse(const Data &d) {
     409      return create(d.ps.p4,d.ps.p3,d.ps.p2,d.ps.p1);
     410    }
     411
     412    static bool intersect(const Data &q, const Data &e) {
     413      dim2::Box<double> box1 = q.ps.boundingBox();
     414      dim2::Box<double> box2 = q.ps.boundingBox();
     415      if ((box1 & box2).empty()) return false;
     416      std::vector<Point> v1, v2;
     417      for (int i=0; i<=100; ++i) {
     418        v1.push_back(q.ps(1.0*i/100));
     419        v2.push_back(e.ps(1.0*i/100));
     420      }
     421      return PlaneGraphLineStripEdge::intersect(
     422              PlaneGraphLineStripEdge::create(v1),
     423              PlaneGraphLineStripEdge::create(v2));
     424    }
     425
     426    static bool degenerate(const Data &d) {
     427      if (d.ps.p1 == d.ps.p2 || d.ps.p3 == d.ps.p4)
     428        return true;
     429      return PlaneGraphStraightEdge::intersect(
     430              PlaneGraphStraightEdge::create(d.ps.p1,d.ps.p2),
     431              PlaneGraphStraightEdge::create(d.ps.p3,d.ps.p4),
     432              true,true,true,true);
     433    }
     434
     435    static void toEps(std::ostream &os, const Data &d, const Color &col,
     436      double width) {
     437      Point p1 = d.ps.p1;
     438/*      if (d.ps.p1 == d.ps.p2) {
     439        p1 = d.ps(0.5);
     440      }*/
     441      os << width << " setlinewidth "
     442         << col.red() << ' '
     443         << col.green() << ' '
     444         << col.blue() << " setrgbcolor newpath\n"
     445         << p1.x << ' ' << p1.y << " moveto\n"
     446         << d.ps.p2.x << ' ' << d.ps.p2.y << ' '
     447         << d.ps.p3.x << ' ' << d.ps.p3.y << ' '
     448         << d.ps.p4.x << ' ' << d.ps.p4.y << " curveto stroke\n";
     449    }
     450  };
     451
     452  /// \addtogroup graphs
     453  /// @{
     454
     455  ///A data structure for a graph embedded in the plane.
     456
     457  ///\ref PlaneGraph adds geometry features to PlanarGraph. While PlanarGraph
     458  ///must be built from the topology, PlaneGraph requires coordinates for nodes
     459  ///and determines the topology that fits the coordinates.
     460  ///\sa PlanarGraph
     461  template <typename EdgeType = PlaneGraphStraightEdge>
     462  class PlaneGraph : public PlanarGraph {
     463
     464  private:
     465
     466    typedef typename EdgeType::Shape EdgeShape;
     467    typedef typename EdgeType::Point EdgePoint;
     468
     469    NodeMap<EdgePoint> nodepos;
     470    ArcMap<EdgeShape> edgeshapes;
     471   
     472    typedef std::map<EdgePoint, Node, LexicographicPointOrdering<EdgePoint> >
     473      PNMap;
     474    PNMap posnode;
     475
     476    static double angle(double x, double y) {
     477      return atan2(y,x);
     478    }
     479
     480    bool checkIntersections(const EdgeShape &es, Arc b) {
     481      Arc a;
     482      firstCwF(a, leftFace(b));
     483      while (a != INVALID) {
     484        if (EdgeType::intersect(es, edgeshapes[a])) return true;
     485        nextCwF(a);
     486      }
     487      return false;
     488    }
     489
     490  protected:
     491
     492    //Intersection between two faces. O(n*m)
     493    bool polygonIntersect(Face f1, Face f2) {
     494      Arc e1, e2;
     495      for(firstCwF(e1,f1); e1 != INVALID; nextCwF(e1)) {
     496        for(firstCwF(e2,f2); e2 != INVALID; nextCwF(e2)) {
     497          if (EdgeType::intersect(edgeshapes[e1], edgeshapes[e2])) return true;
     498        }
     499      }
     500      return false;
     501    }
     502
     503  public:
     504
     505    ///\brief Constructor.
     506
     507    ///Constructor.
     508    PlaneGraph() : nodepos(*this), edgeshapes(*this) {}
     509
     510    ///\brief Constructor that copies a planar drawing
     511
     512    ///Constructor that copies a planar drawing.
     513    template<typename Graph>
     514    PlaneGraph(const Graph &graph, const lemon::PlanarDrawing<Graph> &drawing) :
     515      nodepos(*this), edgeshapes(*this) {
     516        typename Graph::template NodeMap<Node> nm(graph);
     517        for (typename Graph::NodeIt it(graph); it != INVALID; ++it) {
     518          nm[it] = addNode(drawing[it]);
     519        }
     520        for (typename Graph::EdgeIt it(graph); it != INVALID; ++it) {
     521          addEdge(nm[graph.u(it)], nm[graph.v(it)]);
     522        }
     523    }
     524
     525    template<typename Graph>
     526    PlaneGraph(const Graph &graph, const lemon::PlanarDrawing<Graph> &drawing,
     527        typename Graph::template NodeMap<Node> &nm,
     528        typename Graph::template ArcMap<Arc> &em) :
     529      nodepos(*this), edgeshapes(*this) {
     530        for (typename Graph::NodeIt it(graph); it != INVALID; ++it) {
     531          nm[it] = addNode(drawing[it]);
     532        }
     533        for (typename Graph::EdgeIt it(graph); it != INVALID; ++it) {
     534          em[graph.direct(it,graph.u(it))] = direct(addEdge(nm[graph.u(it)
     535            ], nm[graph.v(it)]),nm[graph.u(it)]);
     536          em[graph.direct(it,graph.v(it))] = oppositeArc(em[graph.direct(
     537            it,graph.u(it))]);
     538        }
     539    }
     540
     541    ///Add a new node to the graph at the given point.
     542    ///\return The new node, or INVALID if the given point already contains a
     543    ///node.
     544    Node addNode(const EdgePoint &point) {
     545      typename PNMap::iterator it = posnode.find(point);
     546      if (it != posnode.end()) return INVALID;
     547      Node n = PlanarGraph::addNode();
     548      nodepos[n] = point;
     549      posnode[point] = n;
     550      return n;
     551    }
     552
     553    ///Add a new node to the graph at the given point.
     554    ///\return The new node, or INVALID if the given point already contains a
     555    ///node.
     556    Node addNode(double x, double y) {
     557      return addNode(EdgePoint(x,y));
     558    }
     559
     560    ///Add a new edge to the graph.
     561
     562    ///Add a new edge to the graph between the nodes \c n1 and \c n2.
     563    ///This function should only be used if the edge type doesn't
     564    ///need parameters other than the locations of the two nodes.
     565    ///\return The new edge, or INVALID if the topology doesn't allow the
     566    ///addition.
     567    Edge addEdge(Node n1, Node n2) {
     568      EdgeShape es = EdgeType::create(nodepos[n1],nodepos[n2]);
     569      return addEdge(n1, n2, es);
     570    }
     571
     572    ///Add a new edge to the graph.
     573
     574    ///Add a new edge to the graph with the given edge shape \c es.
     575    ///\return The new edge, or INVALID if \c es doesn't connect the locations
     576    ///of two existing nodes or the topology doesn't allow the
     577    ///addition.
     578    Edge addEdge(const EdgeShape &es) {
     579      typename PNMap::iterator it = posnode.find(EdgeType::source(es));
     580      if (it == posnode.end()) return INVALID;
     581      Node n1 = it->second;
     582      it = posnode.find(EdgeType::target(es));
     583      if (it == posnode.end()) return INVALID;
     584      Node n2 = it->second;
     585      return addEdge(n1, n2, es);
     586    }
     587
     588    ///Add a new edge to the graph.
     589
     590    ///Add a new edge to the graph between the nodes \c n1 and \c n2, using the
     591    ///edge shape \c es. The coordinates of \c n1, \c n2 and \c es will be used
     592    ///to determine the correct topology.
     593    ///\note For performance reasons, this function doesn't check if the source
     594    ///and target coordinates of \c es match those of \c n1 and \c n2.
     595    ///\return The new edge, or INVALID if the topology doesn't allow the
     596    ///addition.
     597    Edge addEdge(Node n1, Node n2, const EdgeShape &es) {
     598      if (EdgeType::degenerate(es)) return INVALID;
     599      EdgeShape esr = EdgeType::reverse(es);
     600      Arc p_u;
     601      firstOut(p_u,n1);
     602      if (p_u != INVALID) {
     603        Arc p_u0 = p_u;
     604        nextOut(p_u0);
     605        while (p_u0 != INVALID && !EdgeType::between(es, edgeshapes[p_u],
     606            edgeshapes[p_u0])) {
     607          p_u = p_u0;
     608          nextOut(p_u0);
     609        }
     610      }
     611      Arc p_v;
     612      firstOut(p_v,n2);
     613      if (p_v != INVALID) {
     614        Arc p_v0 = p_v;
     615        nextOut(p_v0);
     616        while (p_v0 != INVALID && !EdgeType::between(esr, edgeshapes[p_v],
     617            edgeshapes[p_v0])) {
     618          p_v = p_v0;
     619          nextOut(p_v0);
     620        }
     621      }
     622      if (p_u != INVALID && p_v != INVALID) {
     623        if (idComponent(n1) == idComponent(n2)) {
     624          if (leftFace(p_u) != leftFace(p_v))
     625            return INVALID;
     626        } else {
     627          if (polygonIntersect(leftFace(p_u),leftFace(p_v)))
     628            return INVALID;
     629        }
     630      }
     631
     632      if (p_u != INVALID && checkIntersections(es,p_u)) return INVALID;
     633      if (p_v != INVALID && (p_u == INVALID || leftFace(p_u) != leftFace(p_v))
     634        && checkIntersections(es,p_v)) return INVALID;
     635
     636#ifdef REMOVE_BEFORE_RELEASE
     637            std::cout << "AddArc: " << id(n1) << "->" << id(n2) << ", p_u: "
     638             << id(p_u) << ", p_v: " << id(p_v) << std::endl;
     639#endif
     640      Edge e = PlanarGraph::addEdge(n1,n2,p_u,p_v);
     641      if (e != INVALID) {
     642        if (n1 != n2 || EdgeType::between(esr,edgeshapes[p_u],es)) {
     643          edgeshapes[direct(e,true)] = es;
     644          edgeshapes[direct(e,false)] = esr;
     645        } else {
     646          edgeshapes[direct(e,true)] = esr;
     647          edgeshapes[direct(e,false)] = es;
     648        }
     649      }
     650      return e;
     651
     652    }
     653
     654    void erase(Node n) {
     655      EdgePoint p = nodepos[n];
     656      posnode.remove(p);
     657      PlanarGraph::erase(n);
     658    }
     659
     660    void erase(Edge e) {
     661      PlanarGraph::erase(e);
     662    }
     663
     664    void clear() {
     665      posnode.clear();
     666      PlanarGraph::clear();
     667    }
     668
     669    ///Get a map of the node locations.
     670    ///Get a map of the node locations.
     671    const NodeMap<EdgePoint> &nodePosMap() const {
     672      return nodepos;
     673    }
     674
     675    ///Tell the location of a node.
     676    ///Tell the location of a node.
     677    const EdgePoint &locate(const Node &n) const {
     678      return nodepos[n];
     679    }
     680
     681    ///Tell the location of an arc.
     682    ///Tell the location of an arc.
     683    const EdgeShape &locate(const Arc &a) const {
     684      return edgeshapes[a];
     685    }
     686
     687
     688    void edgeToEps(std::ostream &os, const Arc &e, const Color &col,
     689      double width) const {
     690      EdgeType::toEps(os,edgeshapes[e],col,width);
     691    }
     692
     693
     694    /// Iterator class for the edge shapes around a node.
     695
     696    /// This iterator goes through the arcs around a node anticlockwise.
     697    class CcwIncIt : public Arc {
     698      const PlaneGraph* _graph;
     699      const Node _node;
     700    public:
     701
     702      CcwIncIt() { }
     703
     704      CcwIncIt(Invalid i) : Arc(i) { }
     705
     706      CcwIncIt(const PlaneGraph& graph, const Node& node)
     707          : _graph(&graph), _node(node) {
     708        _graph->firstOut(*this, node);
     709      }
     710
     711      CcwIncIt(const Graph& graph, const Arc& arc)
     712          : Arc(arc), _graph(&graph) {}
     713
     714      CcwIncIt& operator++() {
     715        _graph->nextOut(*this, _node);
     716        return *this;
     717      }
     718
     719      const EdgeShape &operator*() {
     720        return _graph->edgeshapes[*this];
     721      }
     722
     723    };
     724
     725    /// Iterator class for the edge shapes around a node.
     726
     727    /// This iterator goes through the arcs around a node clockwise.
     728    class CwIncIt : public Arc {
     729      const PlaneGraph* _graph;
     730      const Node _node;
     731    public:
     732
     733      CwIncIt() { }
     734
     735      CwIncIt(Invalid i) : Arc(i) { }
     736
     737      CwIncIt(const PlaneGraph& graph, const Node& node)
     738          : _graph(&graph), _node(node) {
     739        _graph->lastOut(*this, node);
     740      }
     741
     742      CwIncIt(const Graph& graph, const Arc& arc)
     743          : Arc(arc), _graph(&graph) {}
     744
     745      CwIncIt& operator++() {
     746        _graph->prevOut(*this, _node);
     747        return *this;
     748      }
     749
     750      const EdgeShape &operator*() {
     751        return _graph->edgeshapes[*this];
     752      }
     753
     754    };
     755
     756    /// Iterator class for the edge shapes on the boundary of a face.
     757
     758    /// This iterator goes through the arcs on the boundary of a face clockwise.
     759    class CwBoundaryIt {
     760      const PlaneGraph* _graph;
     761      Face _face;
     762      Arc f_arc;
     763    public:
     764
     765      CwBoundaryIt() { }
     766
     767      CwBoundaryIt(const PlaneGraph& graph, const Face& face)
     768          : _graph(&graph), _face(face) {
     769        _graph->firstCwF(f_arc,face);
     770      }
     771
     772      CwBoundaryIt(const Graph& graph, const Arc& arc)
     773          : f_arc(arc), _graph(&graph) {}
     774
     775      CwBoundaryIt& operator++() {
     776        _graph->nextCwF(f_arc);
     777        return *this;
     778      }
     779
     780      const EdgeShape &operator*() {
     781        return _graph->edgeshapes[f_arc];
     782      }
     783
     784    };
     785
     786    /// Iterator class for the edge shapes on the boundary of a face.
     787
     788    /// This iterator goes through the arcs on the boundary of a face
     789    /// anticlockwise.
     790    class CcwBoundaryIt {
     791      const PlaneGraph* _graph;
     792      Face _face;
     793      Arc f_arc;
     794    public:
     795
     796      CcwBoundaryIt() { }
     797
     798      CcwBoundaryIt(const PlaneGraph& graph, const Face& face)
     799          : _graph(&graph), _face(face) {
     800        _graph->firstCcwF(f_arc,face);
     801      }
     802
     803      CcwBoundaryIt(const Graph& graph, const Arc& arc)
     804          : f_arc(arc), _graph(&graph) {}
     805
     806      CcwBoundaryIt& operator++() {
     807        _graph->nextCcwF(f_arc);
     808        return *this;
     809      }
     810
     811      const EdgeShape &operator*() {
     812        return _graph->edgeshapes[f_arc];
     813      }
     814
     815    };
     816
     817  };
     818  /// @}
     819
     820} //namespace lemon
     821
     822#endif
  • test/CMakeLists.txt

    diff -r 0bca98cbebbb -r ea56aa8243b1 test/CMakeLists.txt
    a b  
    3535  min_cost_flow_test
    3636  min_mean_cycle_test
    3737  path_test
     38  planar_graph_test
    3839  planarity_test
    3940  preflow_test
    4041  radix_sort_test
  • test/Makefile.am

    diff -r 0bca98cbebbb -r ea56aa8243b1 test/Makefile.am
    a b  
    3737        test/min_cost_flow_test \
    3838        test/min_mean_cycle_test \
    3939        test/path_test \
     40        test/planar_graph_test \
    4041        test/planarity_test \
    4142        test/preflow_test \
    4243        test/radix_sort_test \
     
    8889test_min_cost_flow_test_SOURCES = test/min_cost_flow_test.cc
    8990test_min_mean_cycle_test_SOURCES = test/min_mean_cycle_test.cc
    9091test_path_test_SOURCES = test/path_test.cc
     92test_planar_graph_test_SOURCES = test/planar_graph_test.cc
    9193test_planarity_test_SOURCES = test/planarity_test.cc
    9294test_preflow_test_SOURCES = test/preflow_test.cc
    9395test_radix_sort_test_SOURCES = test/radix_sort_test.cc
  • new file test/planar_graph_test.cc

    diff -r 0bca98cbebbb -r ea56aa8243b1 test/planar_graph_test.cc
    - +  
     1/* -*- mode: C++; indent-tabs-mode: nil; -*-
     2 *
     3 * This file is a part of LEMON, a generic C++ optimization library.
     4 *
     5 * Copyright (C) 2003-2010
     6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8 *
     9 * Permission to use, modify and distribute this software is granted
     10 * provided that this copyright notice appears in all copies. For
     11 * precise terms see the accompanying LICENSE file.
     12 *
     13 * This software is provided "AS IS" with no warranty of any kind,
     14 * express or implied, and with no claim as to its suitability for any
     15 * purpose.
     16 *
     17 */
     18
     19#include <lemon/concepts/graph.h>
     20#include <lemon/planar_graph.h>
     21#include <lemon/plane_graph.h>
     22
     23#include "test_tools.h"
     24#include "planar_graph_test.h"
     25
     26using namespace lemon;
     27using namespace lemon::concepts;
     28
     29template <class Graph>
     30void checkGraphBuild() {
     31  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     32
     33  Graph G;
     34  checkGraphNodeList(G, 0);
     35  checkGraphEdgeList(G, 0);
     36  checkGraphArcList(G, 0);
     37  checkGraphFaceList(G, 0);
     38
     39  G.reserveNode(3);
     40  G.reserveEdge(3);
     41
     42  Node
     43    n1 = G.addNode(),
     44    n2 = G.addNode(),
     45    n3 = G.addNode();
     46  checkGraphNodeList(G, 3);
     47  checkGraphEdgeList(G, 0);
     48  checkGraphArcList(G, 0);
     49  checkGraphFaceList(G, 3);
     50
     51  Edge e1 = G.addEdge(n1, n2, INVALID, INVALID);
     52  check((G.u(e1) == n1 && G.v(e1) == n2) || (G.u(e1) == n2 && G.v(e1) == n1),
     53        "Wrong edge");
     54
     55  checkGraphNodeList(G, 3);
     56  checkGraphEdgeList(G, 1);
     57  checkGraphArcList(G, 2);
     58  checkGraphFaceList(G, 2);
     59
     60  checkGraphIncEdgeArcLists(G, n1, 1);
     61  checkGraphIncEdgeArcLists(G, n2, 1);
     62  checkGraphIncEdgeArcLists(G, n3, 0);
     63
     64  checkGraphConEdgeList(G, 1);
     65  checkGraphConArcList(G, 2);
     66
     67  Edge e2 = G.addEdge(n2, n1, e1, e1),
     68       e3 = G.addEdge(n2, n3, e1, INVALID);
     69
     70  checkGraphNodeList(G, 3);
     71  checkGraphEdgeList(G, 3);
     72  checkGraphArcList(G, 6);
     73  checkGraphFaceList(G, 2);
     74
     75  checkGraphIncEdgeArcLists(G, n1, 2);
     76  checkGraphIncEdgeArcLists(G, n2, 3);
     77  checkGraphIncEdgeArcLists(G, n3, 1);
     78
     79  checkGraphConEdgeList(G, 3);
     80  checkGraphConArcList(G, 6);
     81
     82  checkArcDirections(G);
     83
     84  checkNodeIds(G);
     85  checkArcIds(G);
     86  checkEdgeIds(G);
     87  checkFaceIds(G);
     88  checkGraphNodeMap(G);
     89  checkGraphArcMap(G);
     90  checkGraphFaceMap(G);
     91}
     92
     93template <class Graph>
     94void checkPlaneGraphBuild() {
     95  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     96
     97  Graph G;
     98  checkGraphNodeList(G, 0);
     99  checkGraphEdgeList(G, 0);
     100  checkGraphArcList(G, 0);
     101  checkGraphFaceList(G, 0);
     102
     103  G.reserveNode(3);
     104  G.reserveEdge(3);
     105
     106  Node
     107    n1 = G.addNode(dim2::Point<double>(0.0,0.0)),
     108    n2 = G.addNode(dim2::Point<double>(1.0,0.0)),
     109    n3 = G.addNode(dim2::Point<double>(0.5,1.0));
     110  checkGraphNodeList(G, 3);
     111  checkGraphEdgeList(G, 0);
     112  checkGraphArcList(G, 0);
     113  checkGraphFaceList(G, 3);
     114
     115  Edge e1 = G.addEdge(n1, n2);
     116  check((G.u(e1) == n1 && G.v(e1) == n2) || (G.u(e1) == n2 && G.v(e1) == n1),
     117        "Wrong edge");
     118
     119  checkGraphNodeList(G, 3);
     120  checkGraphEdgeList(G, 1);
     121  checkGraphArcList(G, 2);
     122  checkGraphFaceList(G, 2);
     123
     124  checkGraphIncEdgeArcLists(G, n1, 1);
     125  checkGraphIncEdgeArcLists(G, n2, 1);
     126  checkGraphIncEdgeArcLists(G, n3, 0);
     127
     128  checkGraphConEdgeList(G, 1);
     129  checkGraphConArcList(G, 2);
     130
     131  Edge e2 = G.addEdge(n2, n3),
     132       e3 = G.addEdge(n3, n1);
     133
     134  checkGraphNodeList(G, 3);
     135  checkGraphEdgeList(G, 3);
     136  checkGraphArcList(G, 6);
     137  checkGraphFaceList(G, 2);
     138
     139  checkGraphIncEdgeArcLists(G, n1, 2);
     140  checkGraphIncEdgeArcLists(G, n2, 2);
     141  checkGraphIncEdgeArcLists(G, n3, 2);
     142
     143  checkGraphConEdgeList(G, 3);
     144  checkGraphConArcList(G, 6);
     145
     146  checkArcDirections(G);
     147
     148  checkNodeIds(G);
     149  checkArcIds(G);
     150  checkEdgeIds(G);
     151  checkFaceIds(G);
     152  checkGraphNodeMap(G);
     153  checkGraphArcMap(G);
     154  checkGraphFaceMap(G);
     155}
     156
     157template <class Graph>
     158void checkGraphArcSplit() {
     159  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     160
     161  Graph G;
     162  Node n0 = G.addNode(),
     163       n1 = G.addNode(),
     164       n2 = G.addNode(),
     165       n3 = G.addNode();
     166  Edge a0 = G.addEdge(n0,n1,INVALID,INVALID),
     167       a1 = G.addEdge(n0,n3,a0,INVALID),
     168       a2 = G.addEdge(n0,n2,a0,INVALID),
     169       a3 = G.addEdge(n1,n3,a0,a1),
     170       a4 = G.addEdge(n3,n2,a1,a2);
     171
     172  checkGraphNodeList(G, 4);
     173  checkGraphEdgeList(G, 5);
     174  checkGraphArcList(G, 10);
     175  checkGraphFaceList(G, 3);
     176
     177  G.split(a1);
     178
     179  checkGraphNodeList(G, 5);
     180  checkGraphEdgeList(G, 6);
     181  checkGraphArcList(G, 12);
     182  checkGraphFaceList(G, 3);
     183
     184  checkGraphIncEdgeArcLists(G, n0, 3);
     185  checkGraphIncEdgeArcLists(G, n1, 2);
     186  checkGraphIncEdgeArcLists(G, n2, 2);
     187  checkGraphIncEdgeArcLists(G, n3, 3);
     188
     189  checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,true)), 4);
     190  checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,false)), 4);
     191  checkGraphBoundaryArcList(G, G.leftFace(G.direct(a1,false)), 4);
     192
     193  checkGraphConEdgeList(G, 6);
     194  checkGraphConArcList(G, 12);
     195
     196}
     197
     198template <class Graph>
     199void checkGraphNodeSplit() {
     200  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     201
     202    Graph G;
     203    Node n0 = G.addNode(),
     204         n1 = G.addNode(),
     205         n2 = G.addNode(),
     206         n3 = G.addNode(),
     207         n4 = G.addNode();
     208    Edge a0 = G.addEdge(n4,n3,INVALID,INVALID),
     209         a1 = G.addEdge(n4,n0,a0,INVALID),
     210         a2 = G.addEdge(n1,n4,INVALID,a1),
     211         a3 = G.addEdge(n4,n2,a2,INVALID),
     212         a4 = G.addEdge(n2,n3,a3,a0);
     213
     214    checkGraphNodeList(G, 5);
     215    checkGraphEdgeList(G, 5);
     216    checkGraphArcList(G, 10);
     217    checkGraphFaceList(G, 2);
     218
     219    checkGraphIncEdgeArcLists(G, n0, 1);
     220    checkGraphIncEdgeArcLists(G, n1, 1);
     221    checkGraphIncEdgeArcLists(G, n2, 2);
     222    checkGraphIncEdgeArcLists(G, n3, 2);
     223    checkGraphIncEdgeArcLists(G, n4, 4);
     224
     225    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,false)), 3);
     226    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,true)), 7);
     227
     228    Node n5 = G.split(n4,a3,a1,true);
     229
     230    checkGraphNodeList(G, 6);
     231    checkGraphEdgeList(G, 6);
     232    checkGraphArcList(G, 12);
     233    checkGraphFaceList(G, 2);
     234
     235    checkGraphIncEdgeArcLists(G, n0, 1);
     236    checkGraphIncEdgeArcLists(G, n1, 1);
     237    checkGraphIncEdgeArcLists(G, n2, 2);
     238    checkGraphIncEdgeArcLists(G, n3, 2);
     239    checkGraphIncEdgeArcLists(G, n4, 3);
     240    checkGraphIncEdgeArcLists(G, n5, 3);
     241
     242    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,false)), 3);
     243    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,true)), 9);
     244}
     245
     246template <class Graph>
     247void checkGraphContract() {
     248  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     249
     250    Graph G;
     251    Node n0 = G.addNode(),
     252         n1 = G.addNode(),
     253         n2 = G.addNode(),
     254         n3 = G.addNode();
     255    Edge a0 = G.addEdge(n0,n1,INVALID,INVALID),
     256         a1 = G.addEdge(n0,n1,a0,a0),
     257         a2 = G.addEdge(n0,n1,a1,a0),
     258         a3 = G.addEdge(n0,n2,a1,INVALID),
     259         a4 = G.addEdge(n0,n3,a3,INVALID),
     260         a5 = G.addEdge(n2,n1,a3,a2),
     261         a6 = G.addEdge(n3,n1,a4,a2);
     262
     263    checkGraphNodeList(G, 4);
     264    checkGraphEdgeList(G, 7);
     265    checkGraphArcList(G, 14);
     266    checkGraphFaceList(G, 5);
     267
     268    checkGraphIncEdgeArcLists(G, n0, 5);
     269    checkGraphIncEdgeArcLists(G, n1, 5);
     270    checkGraphIncEdgeArcLists(G, n2, 2);
     271    checkGraphIncEdgeArcLists(G, n3, 2);
     272
     273    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,false)), 2);
     274    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a0,true)), 2);
     275    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a1,true)), 3);
     276    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a2,false)), 3);
     277    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a3,true)), 4);
     278
     279    G.contract(a0);
     280
     281    checkGraphNodeList(G, 3);
     282    checkGraphEdgeList(G, 6);
     283    checkGraphArcList(G, 12);
     284    checkGraphFaceList(G, 5);
     285
     286    checkGraphIncEdgeArcLists(G, n0, 8);
     287    checkGraphIncEdgeArcLists(G, n2, 2);
     288    checkGraphIncEdgeArcLists(G, n3, 2);
     289
     290    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a1,false)), 1);
     291    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a2,true)), 1);
     292    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a1,true)), 3);
     293    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a2,false)), 3);
     294    checkGraphBoundaryArcList(G, G.leftFace(G.direct(a3,true)), 4);
     295
     296    checkGraphCommonFaceList(G, n0, n2, 2);
     297    checkGraphCommonFaceList(G, n0, n3, 2);
     298    checkGraphCommonFaceList(G, n2, n3, 1);
     299
     300    checkGraphCommonNodeList(G, G.leftFace(G.direct(a3,true)), G.leftFace(G.
     301            direct(a3,false)), 2);
     302    checkGraphCommonNodeList(G, G.leftFace(G.direct(a1,false)), G.leftFace(G.
     303            direct(a2,true)), 1);
     304
     305}
     306
     307template <class Graph>
     308void checkGraphErase() {
     309  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     310
     311  Graph G;
     312  Node n1 = G.addNode(),
     313       n2 = G.addNode(),
     314       n3 = G.addNode(),
     315       n4 = G.addNode();
     316  Edge e1 = G.addEdge(n1, n2, INVALID, INVALID),
     317       e2 = G.addEdge(n2, n1, e1, e1),
     318       e3 = G.addEdge(n2, n3, e1, INVALID),
     319       e4 = G.addEdge(n1, n4, e2, INVALID),
     320       e5 = G.addEdge(n4, n3, e4, e3);
     321
     322  // Check edge deletion
     323  G.erase(e2);
     324
     325  checkGraphNodeList(G, 4);
     326  checkGraphEdgeList(G, 4);
     327  checkGraphArcList(G, 8);
     328  checkGraphFaceList(G, 2);
     329
     330  checkGraphIncEdgeArcLists(G, n1, 2);
     331  checkGraphIncEdgeArcLists(G, n2, 2);
     332  checkGraphIncEdgeArcLists(G, n3, 2);
     333  checkGraphIncEdgeArcLists(G, n4, 2);
     334
     335  checkGraphBoundaryArcList(G, G.leftFace(G.direct(e1,true)), 4);
     336  checkGraphBoundaryArcList(G, G.leftFace(G.direct(e1,false)), 4);
     337
     338  checkGraphConEdgeList(G, 4);
     339  checkGraphConArcList(G, 8);
     340
     341  // Check node deletion
     342  G.erase(n3);
     343
     344  checkGraphNodeList(G, 3);
     345  checkGraphEdgeList(G, 2);
     346  checkGraphArcList(G, 4);
     347  checkGraphFaceList(G, 1);
     348
     349  checkGraphIncEdgeArcLists(G, n1, 2);
     350  checkGraphIncEdgeArcLists(G, n2, 1);
     351  checkGraphIncEdgeArcLists(G, n4, 1);
     352
     353  checkGraphBoundaryArcList(G, G.leftFace(G.direct(e1,true)), 4);
     354
     355  checkGraphConEdgeList(G, 2);
     356  checkGraphConArcList(G, 4);
     357}
     358
     359
     360template <class Graph>
     361void checkGraphSnapshot() {
     362  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     363
     364  Graph G;
     365  Node n1 = G.addNode(),
     366       n2 = G.addNode(),
     367       n3 = G.addNode();
     368  Edge e1 = G.addEdge(n1, n2, INVALID, INVALID),
     369       e2 = G.addEdge(n2, n1, e1, e1),
     370       e3 = G.addEdge(n2, n3, e1, INVALID);
     371
     372  checkGraphNodeList(G, 3);
     373  checkGraphEdgeList(G, 3);
     374  checkGraphArcList(G, 6);
     375  checkGraphFaceList(G, 2);
     376
     377  typename Graph::Snapshot snapshot(G);
     378
     379  Node n = G.addNode();
     380  Edge ea4 = G.addEdge(n3, n, e3, INVALID),
     381       ea5 = G.addEdge(n, n3, ea4, ea4),
     382       ea6 = G.addEdge(n3, n2, ea5, e3);
     383
     384  checkGraphNodeList(G, 4);
     385  checkGraphEdgeList(G, 6);
     386  checkGraphArcList(G, 12);
     387  checkGraphFaceList(G, 4);
     388
     389  snapshot.restore();
     390
     391  checkGraphNodeList(G, 3);
     392  checkGraphEdgeList(G, 3);
     393  checkGraphArcList(G, 6);
     394  checkGraphFaceList(G, 2);
     395
     396  checkGraphIncEdgeArcLists(G, n1, 2);
     397  checkGraphIncEdgeArcLists(G, n2, 3);
     398  checkGraphIncEdgeArcLists(G, n3, 1);
     399  checkGraphBoundaryArcList(G, G.leftFace(G.direct(e1,true)), 2);
     400  checkGraphBoundaryArcList(G, G.leftFace(G.direct(e1,false)), 4);
     401
     402  checkGraphConEdgeList(G, 3);
     403  checkGraphConArcList(G, 6);
     404
     405  checkNodeIds(G);
     406  checkEdgeIds(G);
     407  checkArcIds(G);
     408  checkFaceIds(G);
     409  checkGraphNodeMap(G);
     410  checkGraphEdgeMap(G);
     411  checkGraphArcMap(G);
     412  checkGraphFaceMap(G);
     413
     414  G.addNode();
     415  snapshot.save(G);
     416
     417  G.addEdge(G.addNode(), G.addNode(), INVALID, INVALID);
     418
     419  snapshot.restore();
     420  snapshot.save(G);
     421
     422  checkGraphNodeList(G, 4);
     423  checkGraphEdgeList(G, 3);
     424  checkGraphArcList(G, 6);
     425  checkGraphFaceList(G, 3);
     426
     427  G.addEdge(G.addNode(), G.addNode(), INVALID, INVALID);
     428
     429  snapshot.restore();
     430
     431  checkGraphNodeList(G, 4);
     432  checkGraphEdgeList(G, 3);
     433  checkGraphArcList(G, 6);
     434  checkGraphFaceList(G, 3);
     435}
     436
     437template <typename Graph>
     438void checkGraphValidity() {
     439  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     440  Graph g;
     441
     442  Node
     443    n1 = g.addNode(),
     444    n2 = g.addNode(),
     445    n3 = g.addNode();
     446
     447  Edge
     448    e1 = g.addEdge(n1, n2, INVALID, INVALID),
     449    e2 = g.addEdge(n2, n3, e1, INVALID);
     450
     451  check(g.valid(n1), "Wrong validity check");
     452  check(g.valid(e1), "Wrong validity check");
     453  check(g.valid(g.direct(e1, true)), "Wrong validity check");
     454
     455  check(!g.valid(g.nodeFromId(-1)), "Wrong validity check");
     456  check(!g.valid(g.edgeFromId(-1)), "Wrong validity check");
     457  check(!g.valid(g.arcFromId(-1)), "Wrong validity check");
     458  check(!g.valid(g.faceFromId(-1)), "Wrong validity check");
     459}
     460
     461template <typename Graph>
     462void checkGraphValidityErase() {
     463  TEMPLATE_GRAPH_TYPEDEFS(Graph);
     464  Graph g;
     465
     466  Node
     467    n1 = g.addNode(),
     468    n2 = g.addNode(),
     469    n3 = g.addNode();
     470
     471  Edge
     472    e1 = g.addEdge(n1, n2, INVALID, INVALID),
     473    e2 = g.addEdge(n2, n3, e1, INVALID);
     474
     475  check(g.valid(n1), "Wrong validity check");
     476  check(g.valid(e1), "Wrong validity check");
     477  check(g.valid(g.direct(e1, true)), "Wrong validity check");
     478
     479  g.erase(n1);
     480
     481  check(!g.valid(n1), "Wrong validity check");
     482  check(g.valid(n2), "Wrong validity check");
     483  check(g.valid(n3), "Wrong validity check");
     484  check(!g.valid(e1), "Wrong validity check");
     485  check(g.valid(e2), "Wrong validity check");
     486
     487  check(!g.valid(g.nodeFromId(-1)), "Wrong validity check");
     488  check(!g.valid(g.edgeFromId(-1)), "Wrong validity check");
     489  check(!g.valid(g.arcFromId(-1)), "Wrong validity check");
     490  check(!g.valid(g.faceFromId(-1)), "Wrong validity check");
     491}
     492
     493void checkGraphs() {
     494  { // Checking PlanarGraph
     495    checkGraphBuild<PlanarGraph>();
     496    checkPlaneGraphBuild<PlaneGraph<> >();
     497    checkGraphArcSplit<PlanarGraph>();
     498    checkGraphNodeSplit<PlanarGraph>();
     499    checkGraphContract<PlanarGraph>();
     500    checkGraphErase<PlanarGraph>();
     501    checkGraphSnapshot<PlanarGraph>();
     502    checkGraphValidityErase<PlanarGraph>();
     503  }
     504}
     505
     506int main() {
     507  checkGraphs();
     508  return 0;
     509}
  • new file test/planar_graph_test.h

    diff -r 0bca98cbebbb -r ea56aa8243b1 test/planar_graph_test.h
    - +  
     1/* -*- mode: C++; indent-tabs-mode: nil; -*-
     2 *
     3 * This file is a part of LEMON, a generic C++ optimization library.
     4 *
     5 * Copyright (C) 2003-2009
     6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8 *
     9 * Permission to use, modify and distribute this software is granted
     10 * provided that this copyright notice appears in all copies. For
     11 * precise terms see the accompanying LICENSE file.
     12 *
     13 * This software is provided "AS IS" with no warranty of any kind,
     14 * express or implied, and with no claim as to its suitability for any
     15 * purpose.
     16 *
     17 */
     18
     19#ifndef LEMON_TEST_GRAPH_TEST_H
     20#define LEMON_TEST_GRAPH_TEST_H
     21
     22#include <set>
     23
     24#include <lemon/core.h>
     25#include <lemon/maps.h>
     26
     27#include "test_tools.h"
     28
     29namespace lemon {
     30
     31  template<class Graph>
     32  void checkGraphNodeList(const Graph &G, int cnt)
     33  {
     34    typename Graph::NodeIt n(G);
     35    for(int i=0;i<cnt;i++) {
     36      check(n!=INVALID,"Wrong Node list linking.");
     37      ++n;
     38    }
     39    check(n==INVALID,"Wrong Node list linking.");
     40    check(countNodes(G)==cnt,"Wrong Node number.");
     41  }
     42
     43  template<class Graph>
     44  void checkGraphFaceList(const Graph &G, int cnt)
     45  {
     46    typename Graph::FaceIt n(G);
     47    for(int i=0;i<cnt;i++) {
     48      check(n!=INVALID,"Wrong Face list linking.");
     49      ++n;
     50    }
     51    check(n==INVALID,"Wrong Face list linking.");
     52    check(countFaces(G)==cnt,"Wrong Face number.");
     53  }
     54
     55  template<class Graph>
     56  void checkGraphArcList(const Graph &G, int cnt)
     57  {
     58    typename Graph::ArcIt e(G);
     59    for(int i=0;i<cnt;i++) {
     60      check(e!=INVALID,"Wrong Arc list linking.");
     61      check(G.oppositeNode(G.source(e), e) == G.target(e),
     62            "Wrong opposite node");
     63      check(G.oppositeNode(G.target(e), e) == G.source(e),
     64            "Wrong opposite node");
     65      ++e;
     66    }
     67    check(e==INVALID,"Wrong Arc list linking.");
     68    check(countArcs(G)==cnt,"Wrong Arc number.");
     69  }
     70
     71  template<class Graph>
     72  void checkGraphOutArcList(const Graph &G, typename Graph::Node n, int cnt)
     73  {
     74    typename Graph::OutArcIt e(G,n);
     75    for(int i=0;i<cnt;i++) {
     76      check(e!=INVALID,"Wrong OutArc list linking.");
     77      check(n==G.source(e),"Wrong OutArc list linking.");
     78      check(n==G.baseNode(e),"Wrong OutArc list linking.");
     79      check(G.target(e)==G.runningNode(e),"Wrong OutArc list linking.");
     80      ++e;
     81    }
     82    check(e==INVALID,"Wrong OutArc list linking.");
     83    check(countOutArcs(G,n)==cnt,"Wrong OutArc number.");
     84  }
     85
     86  template<class Graph>
     87  void checkGraphInArcList(const Graph &G, typename Graph::Node n, int cnt)
     88  {
     89    typename Graph::InArcIt e(G,n);
     90    for(int i=0;i<cnt;i++) {
     91      check(e!=INVALID,"Wrong InArc list linking.");
     92      check(n==G.target(e),"Wrong InArc list linking.");
     93      check(n==G.baseNode(e),"Wrong OutArc list linking.");
     94      check(G.source(e)==G.runningNode(e),"Wrong OutArc list linking.");
     95      ++e;
     96    }
     97    check(e==INVALID,"Wrong InArc list linking.");
     98    check(countInArcs(G,n)==cnt,"Wrong InArc number.");
     99  }
     100
     101  template<class Graph>
     102  void checkGraphEdgeList(const Graph &G, int cnt)
     103  {
     104    typename Graph::EdgeIt e(G);
     105    for(int i=0;i<cnt;i++) {
     106      check(e!=INVALID,"Wrong Edge list linking.");
     107      check(G.oppositeNode(G.u(e), e) == G.v(e), "Wrong opposite node");
     108      check(G.oppositeNode(G.v(e), e) == G.u(e), "Wrong opposite node");
     109      ++e;
     110    }
     111    check(e==INVALID,"Wrong Edge list linking.");
     112    check(countEdges(G)==cnt,"Wrong Edge number.");
     113  }
     114
     115  template<class Graph>
     116  void checkGraphIncEdgeList(const Graph &G, typename Graph::Node n, int cnt)
     117  {
     118    typename Graph::IncEdgeIt e(G,n);
     119    for(int i=0;i<cnt;i++) {
     120      check(e!=INVALID,"Wrong IncEdge list linking.");
     121      check(n==G.u(e) || n==G.v(e),"Wrong IncEdge list linking.");
     122      check(n==G.baseNode(e),"Wrong OutArc list linking.");
     123      check(G.u(e)==G.runningNode(e) || G.v(e)==G.runningNode(e),
     124            "Wrong OutArc list linking.");
     125      ++e;
     126    }
     127    check(e==INVALID,"Wrong IncEdge list linking.");
     128    check(countIncEdges(G,n)==cnt,"Wrong IncEdge number.");
     129  }
     130
     131  template <class Graph>
     132  void checkGraphIncEdgeArcLists(const Graph &G, typename Graph::Node n,
     133                                 int cnt)
     134  {
     135    checkGraphIncEdgeList(G, n, cnt);
     136    checkGraphOutArcList(G, n, cnt);
     137    checkGraphInArcList(G, n, cnt);
     138  }
     139
     140  template<class Graph>
     141  void checkGraphBoundaryArcList(const Graph &G, typename Graph::Face n, int
     142    cnt)
     143  {
     144    typename Graph::CwBoundaryArcIt e(G,n);
     145    for(int i=0;i<cnt;i++) {
     146      check(e!=INVALID,"Wrong CwBoundaryArc list linking.");
     147      check(n==G.w1(e) || n==G.w2(e),"Wrong CwBoundaryArc list linking.");
     148      ++e;
     149    }
     150    check(e==INVALID,"Wrong BoundaryArc list linking.");
     151    check(countBoundaryArcs(G,n)==cnt,"Wrong IncEdge number.");
     152  }
     153
     154  template <class Graph>
     155  void checkGraphConArcList(const Graph &G, int cnt) {
     156    int i = 0;
     157    for (typename Graph::NodeIt u(G); u != INVALID; ++u) {
     158      for (typename Graph::NodeIt v(G); v != INVALID; ++v) {
     159        for (ConArcIt<Graph> a(G, u, v); a != INVALID; ++a) {
     160          check(G.source(a) == u, "Wrong iterator.");
     161          check(G.target(a) == v, "Wrong iterator.");
     162          ++i;
     163        }
     164      }
     165    }
     166    check(cnt == i, "Wrong iterator.");
     167  }
     168
     169  template <class Graph>
     170  void checkGraphConEdgeList(const Graph &G, int cnt) {
     171    int i = 0;
     172    for (typename Graph::NodeIt u(G); u != INVALID; ++u) {
     173      for (typename Graph::NodeIt v(G); v != INVALID; ++v) {
     174        for (ConEdgeIt<Graph> e(G, u, v); e != INVALID; ++e) {
     175          check((G.u(e) == u && G.v(e) == v) ||
     176                (G.u(e) == v && G.v(e) == u), "Wrong iterator.");
     177          i += u == v ? 2 : 1;
     178        }
     179      }
     180    }
     181    check(2 * cnt == i, "Wrong iterator.");
     182  }
     183
     184  template<class Graph>
     185  void checkGraphCommonNodeList(const Graph &G, typename Graph::Face f1,
     186      typename Graph::Face f2, int cnt)
     187  {
     188    typename Graph::CommonNodesIt e(G,f1,f2);
     189    for(int i=0;i<cnt;i++) {
     190      check(e!=INVALID,"Wrong CommonNodeIt.");
     191      ++e;
     192    }
     193    check(e==INVALID,"Wrong CommonNodeIt.");
     194  }
     195
     196  template<class Graph>
     197  void checkGraphCommonFaceList(const Graph &G, typename Graph::Node n1,
     198      typename Graph::Node n2, int cnt)
     199  {
     200    typename Graph::CommonFacesIt e(G,n1,n2);
     201    for(int i=0;i<cnt;i++) {
     202      check(e!=INVALID,"Wrong CommonNodeIt.");
     203      ++e;
     204    }
     205    check(e==INVALID,"Wrong CommonNodeIt.");
     206  }
     207
     208  template <typename Graph>
     209  void checkArcDirections(const Graph& G) {
     210    for (typename Graph::ArcIt a(G); a != INVALID; ++a) {
     211      check(G.source(a) == G.target(G.oppositeArc(a)), "Wrong direction");
     212      check(G.target(a) == G.source(G.oppositeArc(a)), "Wrong direction");
     213      check(G.direct(a, G.direction(a)) == a, "Wrong direction");
     214    }
     215  }
     216
     217  template <typename Graph>
     218  void checkNodeIds(const Graph& G) {
     219    std::set<int> values;
     220    for (typename Graph::NodeIt n(G); n != INVALID; ++n) {
     221      check(G.nodeFromId(G.id(n)) == n, "Wrong id");
     222      check(values.find(G.id(n)) == values.end(), "Wrong id");
     223      check(G.id(n) <= G.maxNodeId(), "Wrong maximum id");
     224      values.insert(G.id(n));
     225    }
     226  }
     227
     228  template <typename Graph>
     229  void checkArcIds(const Graph& G) {
     230    std::set<int> values;
     231    for (typename Graph::ArcIt a(G); a != INVALID; ++a) {
     232      check(G.arcFromId(G.id(a)) == a, "Wrong id");
     233      check(values.find(G.id(a)) == values.end(), "Wrong id");
     234      check(G.id(a) <= G.maxArcId(), "Wrong maximum id");
     235      values.insert(G.id(a));
     236    }
     237  }
     238
     239  template <typename Graph>
     240  void checkEdgeIds(const Graph& G) {
     241    std::set<int> values;
     242    for (typename Graph::EdgeIt e(G); e != INVALID; ++e) {
     243      check(G.edgeFromId(G.id(e)) == e, "Wrong id");
     244      check(values.find(G.id(e)) == values.end(), "Wrong id");
     245      check(G.id(e) <= G.maxEdgeId(), "Wrong maximum id");
     246      values.insert(G.id(e));
     247    }
     248  }
     249
     250  template <typename Graph>
     251  void checkFaceIds(const Graph& G) {
     252    std::set<int> values;
     253    for (typename Graph::FaceIt n(G); n != INVALID; ++n) {
     254      check(G.faceFromId(G.id(n)) == n, "Wrong id");
     255      check(values.find(G.id(n)) == values.end(), "Wrong id");
     256      check(G.id(n) <= G.maxFaceId(), "Wrong maximum id");
     257      values.insert(G.id(n));
     258    }
     259  }
     260
     261  template <typename Graph>
     262  void checkGraphNodeMap(const Graph& G) {
     263    typedef typename Graph::Node Node;
     264    typedef typename Graph::NodeIt NodeIt;
     265
     266    typedef typename Graph::template NodeMap<int> IntNodeMap;
     267    IntNodeMap map(G, 42);
     268    for (NodeIt it(G); it != INVALID; ++it) {
     269      check(map[it] == 42, "Wrong map constructor.");
     270    }
     271    int s = 0;
     272    for (NodeIt it(G); it != INVALID; ++it) {
     273      map[it] = 0;
     274      check(map[it] == 0, "Wrong operator[].");
     275      map.set(it, s);
     276      check(map[it] == s, "Wrong set.");
     277      ++s;
     278    }
     279    s = s * (s - 1) / 2;
     280    for (NodeIt it(G); it != INVALID; ++it) {
     281      s -= map[it];
     282    }
     283    check(s == 0, "Wrong sum.");
     284
     285    // map = constMap<Node>(12);
     286    // for (NodeIt it(G); it != INVALID; ++it) {
     287    //   check(map[it] == 12, "Wrong operator[].");
     288    // }
     289  }
     290
     291  template <typename Graph>
     292  void checkGraphArcMap(const Graph& G) {
     293    typedef typename Graph::Arc Arc;
     294    typedef typename Graph::ArcIt ArcIt;
     295
     296    typedef typename Graph::template ArcMap<int> IntArcMap;
     297    IntArcMap map(G, 42);
     298    for (ArcIt it(G); it != INVALID; ++it) {
     299      check(map[it] == 42, "Wrong map constructor.");
     300    }
     301    int s = 0;
     302    for (ArcIt it(G); it != INVALID; ++it) {
     303      map[it] = 0;
     304      check(map[it] == 0, "Wrong operator[].");
     305      map.set(it, s);
     306      check(map[it] == s, "Wrong set.");
     307      ++s;
     308    }
     309    s = s * (s - 1) / 2;
     310    for (ArcIt it(G); it != INVALID; ++it) {
     311      s -= map[it];
     312    }
     313    check(s == 0, "Wrong sum.");
     314
     315    // map = constMap<Arc>(12);
     316    // for (ArcIt it(G); it != INVALID; ++it) {
     317    //   check(map[it] == 12, "Wrong operator[].");
     318    // }
     319  }
     320
     321  template <typename Graph>
     322  void checkGraphEdgeMap(const Graph& G) {
     323    typedef typename Graph::Edge Edge;
     324    typedef typename Graph::EdgeIt EdgeIt;
     325
     326    typedef typename Graph::template EdgeMap<int> IntEdgeMap;
     327    IntEdgeMap map(G, 42);
     328    for (EdgeIt it(G); it != INVALID; ++it) {
     329      check(map[it] == 42, "Wrong map constructor.");
     330    }
     331    int s = 0;
     332    for (EdgeIt it(G); it != INVALID; ++it) {
     333      map[it] = 0;
     334      check(map[it] == 0, "Wrong operator[].");
     335      map.set(it, s);
     336      check(map[it] == s, "Wrong set.");
     337      ++s;
     338    }
     339    s = s * (s - 1) / 2;
     340    for (EdgeIt it(G); it != INVALID; ++it) {
     341      s -= map[it];
     342    }
     343    check(s == 0, "Wrong sum.");
     344
     345    // map = constMap<Edge>(12);
     346    // for (EdgeIt it(G); it != INVALID; ++it) {
     347    //   check(map[it] == 12, "Wrong operator[].");
     348    // }
     349  }
     350
     351
     352  template <typename Graph>
     353  void checkGraphFaceMap(const Graph& G) {
     354    typedef typename Graph::Face Face;
     355    typedef typename Graph::FaceIt FaceIt;
     356
     357    typedef typename Graph::template FaceMap<int> IntFaceMap;
     358    IntFaceMap map(G, 42);
     359    for (FaceIt it(G); it != INVALID; ++it) {
     360      check(map[it] == 42, "Wrong map constructor.");
     361    }
     362    int s = 0;
     363    for (FaceIt it(G); it != INVALID; ++it) {
     364      map[it] = 0;
     365      check(map[it] == 0, "Wrong operator[].");
     366      map.set(it, s);
     367      check(map[it] == s, "Wrong set.");
     368      ++s;
     369    }
     370    s = s * (s - 1) / 2;
     371    for (FaceIt it(G); it != INVALID; ++it) {
     372      s -= map[it];
     373    }
     374    check(s == 0, "Wrong sum.");
     375
     376    // map = constMap<Node>(12);
     377    // for (NodeIt it(G); it != INVALID; ++it) {
     378    //   check(map[it] == 12, "Wrong operator[].");
     379    // }
     380  }
     381
     382} //namespace lemon
     383
     384#endif