#define INFTY 1000000

/* A library of utility functions for representing graphs and templates
   drawn on surfaces and for 3-coloring them.  */

typedef vector<int> coloring;

/* Computes the winding number of a coloring.  */

static int
elwind (int u, int v)
{
  if (v == u + 1 || v == u - 2)
    return 1;
  else if (u == v + 1 || u == v - 2)
    return -1;
  else
    abort ();
}

static int
winding (coloring &c, int f, int t, int n)
{
  int wind = 0;

  for (int i = f; i != t; i = (i + 1) % n)
    wind += elwind (c[i], c[(i + 1) % n]);

  return wind;
}

/* Cyclic distance from F to T on an N-cycle. */

static int
cdist (int f, int t, int n)
{
  return (t - f + n) % n;
}

static int nedges;

/* Edge in a standard representation of an embedded graph.  Each edge
   is represented as a pair of opposite directed edge, for each edge
   we store a pointer to the previous and next edge in the cyclic order around its
   source vertex, the opposite edge, and the face to the left of the edge.  */

struct edge
{
  int to, left;
  int cid;
  edge *opp, *next, *prev;
  edge *copy;

  edge () : to(-1), left(-1), cid(-1), opp(NULL), next(NULL), prev(NULL), copy(NULL)
    {
      nedges++;
    }

  edge *face_next (void)
    {
      return opp->next;
    }
  edge *face_prev (void)
    {
      return prev->opp;
    }

  void link_after (edge *a)
    {
      edge *nx = a->next;

      a->next = this;
      nx->prev = this;
      prev = a;
      next = nx;
    }

  void unlink ()
    {
      edge *pv = prev;
      edge *nx = next;

      pv->next = nx;
      nx->prev = pv;

      prev = next = NULL;
    }

  ~edge ()
    {
      nedges--;
    }
};

/* Iterator for iterating over edges incident with a vertex.  */

struct vertex;
struct edge_iter_nbr
{
  edge *e, *end;
  bool lst;

  void init (edge *st)
    {
      if (!st)
	abort ();
      e = end = st;
      lst = false;
    }

  edge_iter_nbr (edge *st)
    {
      init (st);
    }
  edge_iter_nbr (vertex &v);

  edge *operator*()
    {
      return e;
    }

  edge_iter_nbr &operator++()
    {
      e = e->next;
      if (e == end)
	lst = true;
      return *this;
    }

  bool end_p ()
    {
      return lst;
    }
};

/* Iterator for iterating over edges incident with a face.  */

struct face;
struct edge_iter_face
{
  edge *e, *end;
  bool lst;

  void init (edge *st)
    {
      if (!st)
	abort ();
      e = end = st;
      lst = false;
    }
  edge_iter_face (edge *st)
    {
      init (st);
    }
  edge_iter_face (face &f);

  edge *operator*()
    {
      return e;
    }

  edge_iter_face &operator++()
    {
      e = e->face_next ();
      if (e == end)
	lst = true;
      return *this;
    }

  bool end_p ()
    {
      return lst;
    }
};

/* Create a new pair of opposite edges.  */

static edge *
edge_pair (int from, int to, int left, int right)
{
  edge *fw = new edge;
  edge *bw = new edge;

  fw->to = to;
  bw->to = from;
  fw->left = left;
  bw->left = right;
  fw->opp = bw;
  bw->opp = fw;

  return fw;
}

/* Is any of the edges in ES incident with V? */

static bool
touches (vector<edge *> &es, int v)
{
  if (es.empty ())
    return false;

  if (es.front ()->opp->to == v)
    return true;

  for (vector<edge *>::iterator e = es.begin (); e != es.end (); ++e)
    if ((*e)->to ==v)
      return true;

  return false;
}

/* A vertex in the representation of an embedded graph.  */

struct vertex
{
  edge *es;
  bool seen;
  int cid, deg;

  vertex ()
    {
      es = NULL;
      deg = 0;
    }

  vertex (edge *ea, int name, int fname)
    {
      es = edge_pair (name, ea->opp->to, fname, fname);
      es->prev = es->next = es;
      es->opp->link_after (ea);
      deg = 1;
    }

  edge *edge_to (int v)
    {
      for (edge_iter_nbr e(*this); !e.end_p (); ++e)
	if ((*e)->to == v)
      	  return *e;
        
      return NULL; 
    }

  int degree (void)
    {
      return deg;
    }

  void update_degree (void)
    {
      int dg = 0;
      for (edge_iter_nbr e(*this); !e.end_p (); ++e)
	dg++;
        
      deg = dg;
    }

  void print (FILE *f)
    {
      for (edge_iter_nbr e(es); !e.end_p (); ++e)
	fprintf (f, " %d", (*e)->to);
    }
};

edge_iter_nbr::edge_iter_nbr (vertex &v)
{
  init (v.es);
}

/* Face of an embedded graph or a template; FLOATING_FACES is
   the multiset assigned to the face in a template. */

typedef map<int,int> flfmap;
#define MAX_FLOATING 7

static void
add_flfmap (flfmap &ffs, flfmap &rg)
{
  for (flfmap::iterator i = rg.begin (); i != rg.end (); ++i)
    {
      if (ffs.count (i->first) > 0)
	ffs[i->first] += i->second;
      else
	ffs[i->first] = i->second;
    }
}

struct face
{
  flfmap floating_faces;
  edge *es;
  bool seen;
  int flength;

  face ()
    {
      es = NULL;
    }

  face (edge *_es)
    {
      es = _es;
    }

  int num_floating (int len) const
    {
      flfmap::const_iterator a = floating_faces.find (len);
      if (a == floating_faces.end ())
	return 0;

      return a->second;
    }

  void add_floating (int len, int cnt = 1)
    {
      if (len > MAX_FLOATING)
	abort ();
      if (cnt == 0)
	return;

      if (floating_faces.count (len) > 0)
	floating_faces[len] += cnt;
      else
	floating_faces[len] = cnt;
    }

  void remove_floating (int len, int cnt = 1)
    {
      if (cnt == 0)
	return;

      int a = floating_faces[len];
      if (a < cnt)
	abort ();
      else if (a == cnt)
	floating_faces.erase (len);
      else
	floating_faces[len] -= cnt;
    }

  void print (FILE *f)
    {
      for (edge_iter_face e(es); !e.end_p (); ++e)
	{
	  fprintf (f, " ");
	  fprintf (f, "%d", (*e)->to);
	}
      for (flfmap::iterator i = floating_faces.begin (); i != floating_faces.end (); ++i)
	{
	  fprintf (f, " fl ");
	  if (i->second != 1)
	    fprintf (f, "%d x ", i->second);
	  fprintf (f, "%d", i->first);
	}
    }

  edge *edge_to (int v)
    {
      for (edge_iter_face e(*this); !e.end_p (); ++e)
	if ((*e)->to == v)
      	  return *e;
      abort ();
    }

  void update_length ()
    {
      int len = 0;

      for (edge_iter_face e(*this); !e.end_p (); ++e)
	len++;

      flength = len;
    }

  int length () const
    {
      return flength;
    }
};

edge_iter_face::edge_iter_face (face &f)
{
  init (f.es);
}

/* List edges in the boundary of a face from F to T.  */

static void
face_arc (edge *f, edge *t, vector<edge *> &arc)
{
  arc.clear ();
  while (f != t)
    {
      arc.push_back (f);
      f = f->face_next ();
    }
  arc.push_back (f);
}

struct path
{
  int f, t, len;
  bool divides;

  path (int _f, int _t, int _len, bool _div)
    {
      f = _f;
      t = _t;
      len = _len;
      divides = _div;
    }
};

struct disk_graph;
struct elist_test
{
  virtual bool operator() (vector<edge *> &apath)
    {
      return true;
    };
};

struct vlist_test
{
  virtual bool operator() (vector<int> vlist) = 0;
};

static bool
smaller_code (vector<int> &c1, vector<int> &c2)
{
  if (c1.size () != c2.size ())
    abort ();

  return c1 < c2;
}

/* When copying a graph, we many need to also replace some pointers
   to the original graph by pointers to the new copy.  This is done
   using the following updaters (for single pointers, vectors of
   pointers, ...) */

struct eref_updater
{
  virtual void operator() (void) const
    {
    }
};

struct join_updater : public eref_updater
{
  const eref_updater *u1, *u2;

  join_updater (const eref_updater &_u1, const eref_updater &_u2) : u1(&_u1), u2(&_u2)
    {
    }

  virtual void operator() (void) const
    {
      (*u1)();
      (*u2)();
    }
};

struct evec_updater : public eref_updater
{
  vector<edge *> *transl;

  evec_updater (vector<edge *> &t)
    {
      transl = &t;
    }

  virtual void operator() (void) const
    {
      for (vector<edge *>::iterator t = transl->begin (); t != transl->end (); ++t)
	(*t) = (*t)->copy;
    }
};

struct single_updater : public eref_updater
{
  edge **transl;

  single_updater (edge *&t)
    {
      transl = &t;
    }

  virtual void operator() (void) const
    {
      (*transl) = (*transl)->copy;
    }
};

struct single_updater_assign : public eref_updater
{
  edge *transl;
  edge **tgt;

  single_updater_assign (edge *what, edge *&to)
    {
      transl = what;
      tgt = &to;
    }

  virtual void operator() (void) const
    {
      *tgt = transl->copy;
    }
};

/* Utilities to maintain an undo list for backtracking in a 3-coloring
   algorithm.  */

struct recorder
{
  virtual unsigned &operator() (unsigned &) = 0;
  virtual void undo () = 0;
};

struct dummy_recorder : public recorder
{
  dummy_recorder ()
    {
    }

  virtual unsigned &operator() (unsigned &val)
    {
      return val;
    }

  virtual void undo ()
    {
    }
};

struct vec_recorder : public recorder
{
  typedef pair<unsigned *,unsigned> rcd;
  vector<rcd> to_undo;

  vec_recorder ()
    {
    }

  virtual unsigned &operator() (unsigned &val)
    {
      to_undo.push_back (rcd (&val, val));
      return val;
    }

  virtual void undo ()
    {
      while (!to_undo.empty ())
	{
	  unsigned *wh = to_undo.back ().first;
	  unsigned val = to_undo.back ().second;

	  *wh = val;
	  to_undo.pop_back ();
	}
    }
};

struct coloring_callback
{
  virtual bool operator() (unsigned col[])
    {
      return true;
    }
};

struct avoid_test
{
  virtual bool avoid_edge (edge *e) const
    {
      return false;
    }
  virtual bool avoid_vertex (int v) const
    {
      return false;
    }
};

typedef pair<edge *,int> flfrec;

/* Embedded graph/template and various utilities to operate on it.
   More involved operations are finding a canonical code (for isomorphism
   testing and storage) and 3-coloring with template constraints.  */

struct graph
{
  vector<vertex> vs;
  vector<face> fs;
#ifdef _OPENMP
  omp_lock_t copy_lock;
#endif

  graph (void)
    {
#ifdef _OPENMP
      omp_init_lock (&copy_lock);
#endif
      vs.resize (2);
      fs.resize (1);
      edge *e1 = edge_pair (0, 1, 0, 0);
      edge *e2 = e1->opp;
      vs[0].es = e1->next = e1->prev = e1;
      vs[1].es = e2->next = e2->prev = e2;
      vs[0].deg = vs[1].deg = 1;
      fs[0].es = e1;
      fs[0].flength = 2;
    }

  graph (int n)
    {
#ifdef _OPENMP
      omp_init_lock (&copy_lock);
#endif
      vs.resize (n);
      fs.resize (2);

      edge *ae[n];
      for (int i = 0; i < n; i++)
	{
	  vs[i].es = ae[i] = edge_pair (i, (i + 1) % n, 0, 1);
	  vs[i].deg = 2;
	}
      for (int i = 0; i < n; i++)
	{
	  edge *oe = ae[(i + n - 1) % n]->opp;
	  ae[i]->next = ae[i]->prev = oe;
	  oe->prev = oe->next = ae[i];
	}

      fs[0].es = ae[0];
      fs[1].es = ae[0]->opp;
      fs[0].flength = fs[1].flength = n;
    }

  graph (graph &g, const eref_updater &update = eref_updater (), bool unlock = true)
    {
#ifdef _OPENMP
      omp_init_lock (&copy_lock);
      omp_set_lock (&g.copy_lock);
#endif
      vs = g.vs;
      fs = g.fs;

      edge *e;
      int v;

      for (v = 0; v < nvs(); v++)
	if (vs[v].es)
	  for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	    (*ae)->copy = new edge;
      for (v = 0; v < nvs(); v++)
	{
	  if (!vs[v].es)
	    continue;

	  for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	    {
	      e = *ae;
	      edge *cp = e->copy;

	      cp->to = e->to;
	      cp->left = e->left;
	      cp->opp = e->opp->copy;
	      cp->next = e->next->copy;
	      cp->prev = e->prev->copy;
	    }
	  vs[v].es = vs[v].es->copy;
	}
      for (int f = 0; f < nfs(); f++)
	if (fs[f].es)
	  fs[f].es = fs[f].es->copy;
      update ();

#ifdef _OPENMP
      if (unlock)
	omp_unset_lock (&g.copy_lock);
#endif
    }

  graph (vector<int> &code)
    {
#ifdef _OPENMP
      omp_init_lock (&copy_lock);
#endif
      int maxcid = code.size ();

      int cid_to_vno[maxcid];
      edge *cid_to_edge[maxcid];
      int cid_to_edge_src[maxcid];
      vector<flfrec> fles;
      vector<int> avertex;
      vector<edge *> apathe;

      for (int i = 0; i < maxcid; i++)
	{
	  cid_to_vno[i] = -1;
	  cid_to_edge[i] = NULL;
	}

      vs.push_back (vertex ());
      cid_to_vno[0] = 0;
      avertex.push_back (0);

      vector<int>::iterator p = code.begin ();

      while (p != code.end ())
	{
	  edge *e = new edge;
	  while (*p < -1)
	    {
	      int flen = 3 - *p;
	      fles.push_back (flfrec (e, flen));
	      ++p;
	    }
	  int v = avertex.back ();
	  int eno = *p;
	  ++p;

	  if (!vs[v].es)
    	    e->prev = e->next = e;
	  else
    	    e->link_after (vs[v].es);
	  vs[v].es = e;

	  edge *oe;
	  int tgt = -1;
	  if (eno == -1)
	    {
	      oe = apathe.back ();
	      apathe.pop_back ();
	      avertex.pop_back ();
	      tgt = avertex.back ();
	    }
	  else if (cid_to_edge[eno])
	    {
	      oe = cid_to_edge[eno];
	      tgt = cid_to_edge_src[eno];
	    }
	  else
	    oe = NULL;

	  if (oe)
	    {
	      e->to = tgt;
	      e->opp = oe;
	      oe->opp = e;
	      continue;
	    }

	  cid_to_edge[eno] = e;
	  cid_to_edge_src[eno] = v;

	  int tgt_cid = *p;
	  p++;

	  if (cid_to_vno[tgt_cid] != -1)
	    {
    	      tgt = cid_to_vno[tgt_cid];
	      e->to = tgt;
	      continue;
	    }

	  tgt = nvs ();
	  vs.push_back (vertex ());
	  cid_to_vno[tgt_cid] = tgt;
	  e->to = tgt;

	  avertex.push_back (tgt);
	  apathe.push_back (e);
	}

      for (int v = 0; v < nvs (); v++)
	{
	  vs[v].update_degree ();
	  for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	    {
	      if ((*ae)->left != -1)
		continue;

	      int f = nfs ();
	      fs.push_back (face (*ae));
	      for (edge_iter_face fe(*ae); !fe.end_p (); ++fe)
		(*fe)->left = f;
	      fs.back ().update_length ();
	    }
	}

      for (vector<flfrec>::iterator flf = fles.begin (); flf != fles.end (); ++flf)
	fs[flf->first->left].add_floating (flf->second);
    }

  ~graph()
    {
      for (int i = 0; i < nvs(); i++)
	{
	  vector<edge *> to_del;
	  if (!vs[i].es)
	    continue;

	  for (edge_iter_nbr e(vs[i]); !e.end_p (); ++e)
	    to_del.push_back (*e);

	  for (vector<edge *>::iterator r = to_del.begin (); r != to_del.end (); ++r)
	    delete (*r);
	}

#ifdef _OPENMP
      omp_destroy_lock (&copy_lock);
#endif
    }

  int nvs (void)
    {
      return vs.size ();
    }

  int nfs (void)
    {
      return fs.size ();
    }

  edge *add_pendant_vertex (edge *at)
    {
      int fno = at->left;
      int u = at->to;
      int v = nvs ();
      vs.push_back (vertex (at->opp, v, fno));
      vs[u].deg++;
      fs[fno].flength += 2;

      return vs[v].es->opp;
    }

  void split_face (edge *etou, edge *etov, int len, const flfmap &move_right, vector<edge *> &aes)
    {
      int fno = etou->left;
      int u = etou->to;
      int v = etov->to;
      edge *ue = etou->opp;
      edge *ve = etov->opp;

      for (int i = 0; i < len - 1; i++)
	{
	  int x = nvs();
	  vs.push_back (vertex (ue, x, fno));
	  vs[u].deg++;
	  ue = vs[x].es;
	  u = x;
	  aes.push_back (ue->opp);
	}

      int nwf = nfs();
      edge *nwe = edge_pair (u, v, fno, nwf);
      vs[u].deg++;
      vs[v].deg++;
      aes.push_back (nwe);

      nwe->link_after (ue);
      nwe->opp->link_after (ve);

      fs[fno].es = nwe;
      fs.push_back (face (nwe->opp));

      for (flfmap::const_iterator i = move_right.cbegin (); i != move_right.cend (); ++i)
	{
      	  fs[fno].remove_floating (i->first, i->second);
	  fs[nwf].add_floating (i->first, i->second);
	}

      for (edge_iter_face tr(fs[nwf]); !tr.end_p (); ++tr)
	(*tr)->left = nwf;

      fs[fno].update_length ();
      fs[nwf].update_length ();
    }

  void subdivide (edge *e)
    {
      int lf = e->left;
      int rf = e->opp->left;
      
      int x = nvs();
      vs.push_back (vertex (e->opp, x, lf));
      edge *xe = vs[x].es;
      xe->opp->left = rf;

      vs[e->to].es = xe->opp;
      e->opp->unlink ();
      e->opp->link_after (xe);
      e->to = x;
      vs[x].deg++;

      fs[lf].flength++;
      fs[rf].flength++;
    }

  void decontract_4face (edge *ef, edge *et)
    {
      int v = ef->opp->to;

      edge *bf = et->next;
      edge *bt = ef->prev;

      if (bf == ef)
	{
	  flfmap dummy;
	  vector<edge *> es;
	  split_face (ef, et->opp->face_prev (), 2, dummy, es);
	  return;
	}

      int w = nvs ();

      vs[v].es = ef;
      vs.push_back (vertex ());
      vs[w].es = bf;

      et->next = ef;
      ef->prev = et;
      bt->next = bf;
      bf->prev = bt;
      for (edge_iter_nbr e(vs[w]); !e.end_p (); ++e)
	(*e)->opp->to = w;

      int ff = nfs ();
      fs.push_back (face (ef));
      int fcl = ef->left;
      int fcr = bf->left;
      fs[fcl].es = bt->opp;
      fs[fcr].es = bf;

      edge *nf = edge_pair (w, ef->to, fcl, ff);
      edge *nt = edge_pair (w, et->to, ff, fcr);
      nf->link_after (bt);
      nt->link_after (nf);
      nf->opp->link_after (ef->opp);
      nt->opp->link_after (et->opp->prev);
      for (edge_iter_face e(fs[ff]); !e.end_p (); ++e)
	(*e)->left = ff;
      vs[v].update_degree ();
      vs[w].update_degree ();
      vs[ef->to].update_degree ();
      vs[et->to].update_degree ();
      fs[ff].update_length ();
      fs[fcl].update_length ();
      fs[fcr].update_length ();
    }

  void shear (edge *e, edge *into)
    {
      int v = e->opp->to;
      int f = into->left;

      edge *ep = edge_pair (v, e->to, e->left, f);
      fs[e->left].es = ep;
      e->left = f;
      ep->link_after (e->prev);
      ep->opp->link_after (e->opp);

      edge *l = into->prev;
      int w = nvs ();

      vs[v].es = ep;
      vs.push_back (vertex ());
      vs[w].es = l;

      ep->next = into;
      into->prev = ep;
      l->next = e;
      e->prev = l;
      for (edge_iter_nbr e(vs[w]); !e.end_p (); ++e)
	(*e)->opp->to = w;

      vs[v].update_degree ();
      vs[w].update_degree ();
      fs[f].update_length ();
    }

  void print (FILE *f)
    {
      for (int v = 0; v < nvs (); v++)
	if (vs[v].es)
	  {
	    fprintf (f, "  vertex %d:", v);
	    vs[v].print (f);
	    fprintf (f, "\n");
	  }
      for (int a = 0; a < nfs (); a++)
	if (fs[a].es)
	  {
	    fprintf (f, "  face %d:", a);
	    fs[a].print (f);
	    fprintf (f, "\n");
	  }
    }

  bool cycles_rec (int minv, int rem, vector<edge *> &cyc, elist_test &eval)
    {
      int v = cyc.back ()->to;

      for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	{
	  edge *e = *ae;

	  if (e->to == minv && e->opp != cyc.front ())
	    {
	      cyc.push_back (e);
	      if (eval (cyc))
		return true;
	      cyc.pop_back ();
	      continue;
	    }

	  if (rem == 1 || e->to < minv || touches (cyc, e->to))
	    continue;

	  cyc.push_back (e);
	  if (cycles_rec (minv, rem - 1, cyc, eval))
	    return true;
	  cyc.pop_back ();
	}

      return false;
    }

  bool cycles_withmin (int minv, int maxlen, elist_test &eval)
    {
      for (edge_iter_nbr ae(vs[minv]); !ae.end_p (); ++ae)
	{
	  edge *e = *ae;
	  vector<edge *> cyc;

	  if (e->to < minv)
	    continue;

	  cyc.push_back (e);
	  if (cycles_rec (minv, maxlen - 1, cyc, eval))
	    return true;
	}

      return false;
    }

  bool cycles (int maxlen, elist_test &eval)
    {
      int v;

      for (v = 0; v < nvs (); v++)
	if (vs[v].es && cycles_withmin (v, maxlen, eval))
	  return true;

      return false;
    }

  bool any_cycle (int maxlen)
    {
      elist_test any;
      return cycles (maxlen, any);
    }

  void get_code_dfs (edge *e, bool rev, vector<int> &code, int &acid)
    {
      int fno = rev ? e->opp->left : e->left;
      if (!fs[fno].seen)
	{
	  fs[fno].seen = true;
	  flfmap &ff = fs[fno].floating_faces;
	  for (flfmap::iterator j = ff.begin (); j != ff.end (); ++j)
	    for (int k = 0; k < j->second; k++)
	      code.push_back (3 - j->first);
	}

      if (e->cid != -1)
	{
	  /* Forward edge.  */
	  code.push_back (e->cid);
	  return;
	}
      code.push_back (acid);
      e->cid = e->opp->cid = acid++;

      if (vs[e->to].cid != -1)
	{
	  /* Back edge.  */
	  code.push_back (vs[e->to].cid);
	  return;
	}
      /* Tree edge.  */
      code.push_back (acid);
      vs[e->to].cid = acid++;

      edge *end = e->opp;
      for (edge *ae = rev ? end->prev : end->next; ae != end; ae = rev ? ae->prev : ae->next)
	get_code_dfs (ae, rev, code, acid);
      code.push_back (-1);
    }

  void get_code (edge *from, bool rev, vector<int> &code)
    {
      code.clear ();

      for (int v = 0; v < nvs(); v++)
	{
	  vs[v].cid = -1;
	  if (vs[v].es)
	    {
	      for (edge_iter_nbr e(vs[v]); !e.end_p (); ++e)
		(*e)->cid = -1;
	    }
	}
      for (int f = 0; f < nfs (); f++)
	fs[f].seen = false;

      int acid = 0;
      vs[from->opp->to].cid = acid++;

      edge *e = from;
      do
	{
	  get_code_dfs (e, rev, code, acid);
	  e = rev ? e->prev : e->next;
	} while (e != from);
    }

  void get_min_code (vector<int> &code)
    {
      bool fst = true;
      vector<int> acode;

      for (int v = 0; v < nvs(); v++)
	{
	  if (vs[v].es == NULL)
	    continue;

	  for (edge_iter_nbr e(vs[v]); !e.end_p (); ++e)
	    {
	      get_code (*e, false, acode);
	      if (fst || smaller_code (acode, code))
		code = acode;
	      fst = false;
	      get_code (*e, true, acode);
	      if (smaller_code (acode, code))
		code = acode;
	    }
	}
    }

  void remove_edge (edge *e)
    {
      int f1 = e->left, f2 = e->opp->left;

      if (f1 != f2)
	{
	  for (edge_iter_face ae(fs[f2]); !ae.end_p (); ++ae)
	    (*ae)->left = f1;
	  fs[f2].es = NULL;

	  flfmap &rem = fs[f2].floating_faces;
	  for (flfmap::iterator i = rem.begin (); i != rem.end (); ++i)
	    fs[f1].add_floating (i->first, i->second);
	}

      edge *fe = fs[f1].es;
      while (fe == e || fe == e->opp)
	fe = fe->face_next ();
      fs[f1].es = fe;

      int u = e->opp->to, v = e->to;

      if (vs[u].es == e)
	vs[u].es = e->next;
      if (vs[u].es == e)
	vs[u].es = NULL;

      if (vs[v].es == e->opp)
	vs[v].es = e->opp->next;
      if (vs[v].es == e->opp)
	vs[v].es = NULL;

      vs[u].deg--;
      vs[v].deg--;

      e->unlink ();
      e->opp->unlink ();
      delete e->opp;
      delete e;

      fs[f1].update_length ();
    }

  bool has_le3cycle (void)
    {
      int v;

      for (v = 0; v < nvs (); v++)
	{
	  if (!vs[v].es)
	    continue;

	  for (edge_iter_nbr e1(vs[v]); !e1.end_p (); ++e1)
	    {
	      if ((*e1)->to == v)
		return true;

	      for (edge_iter_nbr e2(vs[(*e1)->to]); !e2.end_p (); ++e2)
		{
		  if ((*e2)->to == v)
		    {
		      if (*e2 != (*e1)->opp)
			return true;

		      continue;
		    }

		  if (vs[v].edge_to ((*e2)->to))
		    return true;
		}
	    }
	}

      return false;
    }

  void set_color (unsigned color[], unsigned list[], int v, int col, recorder &rec, edge *without)
    {
      if (color[v])
	abort ();

      rec (color[v]) = col;

      for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	{
	  if (*ae == without || (*ae)->opp == without)
	    continue;

	  int u = (*ae)->to;
	  unsigned l = list[u];

	  if (l & (1 << col))
	    {
	      l &= ~(1 << col);
	      rec (list[u]) = l;
	    }
	}
    }

  bool is_colorable_rec (unsigned color[], unsigned list[], edge *without, coloring_callback &cb)
    {
      int best_v = -1, bests = 4;
      int nv = nvs ();

      for (int v = 0; v < nv; v++)
	{
	  if (vs[v].es == NULL || color[v] != 0)
	    continue;

	  unsigned l = list[v];
	  int sz = ((l & 2) != 0) + ((l & 4) != 0) + ((l & 8) != 0);

	  if (sz < bests)
	    {
	      bests = sz;
	      best_v = v;
	    }
	}

      if (best_v == -1)
	return cb (color);

      if (bests == 0)
	return false;

      unsigned l = list[best_v];

      for (int c = 1; c <= 3; c++)
	if ((l >> c) & 1)
	  {
	    vec_recorder rec;
	    set_color (color, list, best_v, c, rec, without);
	    if (is_colorable_rec (color, list, without, cb))
	      return true;

	    rec.undo ();
	  }

      return false;
    }

  bool find_coloring (unsigned color[], coloring_callback &cb, const coloring &precoloring, edge *without = NULL)
    {
      int nv = nvs ();
      unsigned list[nv];
      edge *e = NULL;

      for (int v = 0; v < nv; v++)
	{
	  color[v] = 0;
	  list[v] = 0xe;

	  if (vs[v].es)
	    e = vs[v].es;
	}
      if (!e)
	return cb (color);
      if (e == without || e->opp == without)
	e = e->next;
      if (e == without || e->opp == without)
	abort ();

      dummy_recorder dr;
      if (precoloring.empty ())
	{
	  set_color (color, list, e->to, 1, dr, without);
	  set_color (color, list, e->opp->to, 2, dr, without);
	}
      else
	{
	  int n = precoloring.size ();
	  for (int v = 0; v < n; v++)
	    if (vs[v].es)
	      {
		int c = precoloring[v];
		if (((list[v] >> c) & 1) == 0)
		  return false;
		set_color (color, list, v, c, dr, without);
	      }
	}

      return is_colorable_rec (color, list, without, cb);
    }

  bool find_coloring (unsigned color[], coloring_callback &cb, edge *without = NULL)
    {
      return find_coloring (color, cb, coloring{}, without);
    }

  int distance (int u, int v, const avoid_test &avoid = avoid_test ())
    {
      for (int i = 0; i < nvs(); i++)
	vs[i].seen = false;

      vector<int> que, nxt;
      int dist = 0;
      que.push_back (u);
      vs[u].seen = true;

      while (!que.empty ())
	{
	  while (!que.empty ())
	    {
	      int a = que.back ();
	      que.pop_back ();

	      for (edge_iter_nbr ae(vs[a]); !ae.end_p (); ++ae)
		{
		  edge *e = *ae;

		  if (e->to == v && !avoid.avoid_edge (e))
		    return dist + 1;

		  if (avoid.avoid_vertex (e->to))
		    continue;

		  if (!vs[e->to].seen)
		    {
		      vs[e->to].seen = true;
		      nxt.push_back (e->to);
		    }
		}
	    }

       	  que = nxt;
	  dist++;
	  nxt.clear ();
	}

      return INFTY;
    }

  void mark_reachable_faces (int fno)
    {
      if (fs[fno].seen)
	return;
      fs[fno].seen = true;

      for (edge_iter_face e(fs[fno]); !e.end_p (); ++e)
	if ((*e)->cid)
	  mark_reachable_faces ((*e)->opp->left);
    }

  int winding_number (int fno, unsigned color[], edge *without = NULL)
    {
      int w = 0;

      for (edge_iter_face e(fs[fno]); !e.end_p (); ++e)
	{
	  if (*e == without || (*e)->opp == without)
	    continue;

	  int u = color[(*e)->opp->to];
	  int v = color[(*e)->to];

	  w += elwind (u, v);
	}

      return w;
    }

  void verify (void) __attribute__ ((used))
    {
      for (int v = 0; v < nvs (); v++)
	if (vs[v].es)
	  {
	    int deg = vs[v].degree ();
	    vs[v].update_degree ();
	    if (deg != vs[v].degree ())
	      abort ();
	  }
      for (int f = 0; f < nfs (); f++)
	if (fs[f].es)
	  {
	    for (edge_iter_face e(fs[f]); !e.end_p (); ++e)
	      if ((*e)->left != f)
		abort ();

	    int len = fs[f].length ();
	    fs[f].update_length ();
	    if (len != fs[f].length ())
	      abort ();
	  }
    }

  void clear_marks (void)
    {
      for (int f = 0; f < nfs (); f++)
	fs[f].seen = false;
      for (int v = 0; v < nvs (); v++)
	{
	  vs[v].seen = false;
	  if (vs[v].es)
	    for (edge_iter_nbr e(vs[v]); !e.end_p (); ++e)
	      (*e)->cid = 1;
	}
    }
  void mark_interior (vector<edge *> &cycle)
    {
      clear_marks ();
      for (vector<edge *>::iterator e = cycle.begin (); e != cycle.end (); ++e)
	{
	  (*e)->cid = (*e)->opp->cid = 0;
	  vs[(*e)->to].seen = true;
	}

      for (vector<edge *>::iterator e = cycle.begin (); e != cycle.end (); ++e)
	mark_reachable_faces ((*e)->left);
    }

  bool paths_to_rec (vector<edge *> &apath, int v, int tov, int rem, elist_test &eval)
    {
      if (v == tov)
	return eval (apath);
      if (rem <= 0)
	return false;

      for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	{
	  edge *e = *ae;

	  if (!vs[e->to].seen)
	    {
	      vs[e->to].seen = true;
	      apath.push_back (e);
	      if (paths_to_rec (apath, e->to, tov, rem - 1, eval))
		return true;
	      apath.pop_back ();
	      vs[e->to].seen = false;
	    }
	}

      return false;
    }

  bool paths_to (int v, int vto, int max_len, elist_test &eval)
    {
      vector<edge *> apath;

      for (int a = 0; a < nvs(); a++)
	vs[a].seen = false;

      vs[v].seen = true;
      return paths_to_rec (apath, v, vto, max_len, eval);
    }

  bool twoec_components_rec (edge *from, int v, int lev, int &reach, vector<int> &open, vector<int> &levels, vlist_test &eval)
    {
      levels[v] = lev;
      open.push_back (v);

      for (edge_iter_nbr e(vs[v]); !e.end_p (); ++e)
	{
	  int u = (*e)->to;

	  if ((*e)->opp == from || vs[u].seen)
	    continue;

	  if (levels[u] == nvs ())
	    {
	      int ureach = nvs ();
	      if (twoec_components_rec (*e, u, lev + 1, ureach, open, levels, eval))
		return true;
	      reach = min (reach, ureach);
	    }
	  else
	    reach = min (reach, levels[u]);
	}

      if (reach >= lev)
	{
	  vector<int> comp;
	  int x;

	  do
	    {
	      x = open.back ();
	      comp.push_back (x);
	      open.pop_back ();
	      vs[x].seen = true;
	    }
	  while (x != v);

	  if (eval (comp))
	    return true;
	}
      return false;
    }

  bool twoec_components (int v, vlist_test &eval)
    {
      vector<int> levels (nvs (), nvs ());
      vector<int> open;
      int minlev = nvs ();

      return twoec_components_rec (NULL, v, 0, minlev, open, levels, eval);
    }

  edge *corresponding_edge (edge *oge)
    {
      return vs[oge->opp->to].edge_to (oge->to);
    }
};

struct avoid_external_test : public avoid_test
{
  int nout;
  avoid_external_test (int _nout)
    {
      nout = _nout;
    }

  virtual bool avoid_edge (edge *e) const
    {
      return e->left == 0 || e->opp->left == 0;
    }
  virtual bool avoid_vertex (int v) const
    {
      return v < nout;
    }
};

/* Specialization of the GRAPH data structure to the case of a graph embedded in
   a disk with the boundary traced by a NOUT-cycle, and some utility functions
   only relevant for this case.  */

struct disk_graph : public graph
{
  int nout;

  disk_graph (int n) : graph (n)
    {
      nout = n;
    }

  disk_graph (disk_graph &g, const eref_updater &update = eref_updater ()) : graph (g, update, false)
    {
      nout = g.nout;
#ifdef _OPENMP
      omp_unset_lock (&g.copy_lock);
#endif
    }

  void print (FILE *f)
    {
      fprintf (f, "Outer face 0 .. %d\n", nout - 1);
     
      graph::print (f);
    }

  bool ip_rec (vector<edge *> &apath, int v, edge *without, elist_test &eval)
    {
      for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	{
	  edge *e = *ae;

	  if (e == without || e->opp == without)
	    continue;

	  if (!vs[e->to].seen)
	    {
	      vs[e->to].seen = true;
	      apath.push_back (e);
	      if (ip_rec (apath, e->to, without, eval))
		return true;
	      apath.pop_back ();
	      vs[e->to].seen = false;
	    }

	  if (e->to >= nout)
	    continue;

	  if (apath.empty ())
	    {
	      if (e->left == 0
		  || e->opp->left == 0)
		continue;
	    }
	  else if (e->to == apath.front ()->opp->to)
	    continue;

	  apath.push_back (e);
	  if (eval (apath))
	    return true;
	  apath.pop_back ();
	}

      return false;
    }

  bool interior_paths_from (int v, edge *without, elist_test &eval)
    {
      int a;

      vector<edge *> apath;

      for (a = 0; a < nout; a++)
	vs[a].seen = true;
      for (; a < nvs(); a++)
	vs[a].seen = false;

      return ip_rec (apath, v, without, eval);
    }

  bool interior_paths (edge *without, elist_test &eval)
    {
      int v;

      for (v = 0; v < nout; v++)
	if (interior_paths_from (v, without, eval))
	  return true;

      return false;
    }

  int internal_distance (int u, int v)
    {
      return distance (u, v, avoid_external_test (nout));
    }

  bool forced_deg2 (int v)
    {
      edge *e1 = vs[v].es;
      edge *e2 = e1->next;

      if (e2->next != e1)
	return false;

      int fno = e1->left == 0 ? e2->left : e1->left;
      int len = fs[fno].length ();
      return len <= 5;
    }

  bool precoloring_extends (unsigned boundcol[], coloring_callback &cb, edge *without = NULL)
    {
      int v, nv = nvs ();
      unsigned color[nv];
      unsigned list[nv];

      for (v = 0; v < nv; v++)
	{
	  color[v] = 0;
	  list[v] = 0xe;
	}

      dummy_recorder dr;
      for (v = 0; v < nout; v++)
	{
	  int c = boundcol[v];
	  if (((list[v] >> c) & 1) == 0)
	    return false;
	  set_color (color, list, v, c, dr, without);
	}

      return is_colorable_rec (color, list, without, cb);
    }

  void get_min_code (vector<int> &code)
    {
      bool fst = true;
      vector<int> acode;

      for (edge_iter_face e(fs[0]); !e.end_p (); ++e)
	{
	  get_code (*e, false, acode);
	  if (fst || smaller_code (acode, code))
	    code = acode;
	  fst = false;
	  get_code ((*e)->opp, true, acode);
	  if (smaller_code (acode, code))
    	    code = acode;
	}
    }

  void floating_inside (vector<edge *> &cycle, flfmap &flin)
    {
      flin.clear ();

      mark_interior (cycle);

      if (fs[0].seen)
	abort ();
      for (int f = 1; f < nfs (); f++)
	if (fs[f].seen)
	  {
	    flfmap &ff = fs[f].floating_faces;

	    for (flfmap::iterator i = ff.begin (); i != ff.end (); ++i)
	      {
		if (flin.count (i->first) > 0 || i->second != 1)
		  abort ();
		flin[i->first] = 1;
	      }
	  }
    }

  bool inward_facing (vector<edge *> &cycle)
    {
      clear_marks ();

      for (vector<edge *>::iterator e = cycle.begin (); e != cycle.end (); ++e)
	(*e)->cid = 0;

      for (vector<edge *>::iterator e = cycle.begin (); e != cycle.end (); ++e)
	mark_reachable_faces ((*e)->left);

      return !fs[0].seen;
    }
};

typedef pair<edge *,edge *> mergedsc;

struct mvec_updater : public eref_updater
{
  vector<mergedsc> *transl;

  mvec_updater (vector<mergedsc> &t)
    {
      transl = &t;
    }

  virtual void operator() (void) const
    {
      for (vector<mergedsc>::iterator t = transl->begin (); t != transl->end (); ++t)
	{
	  t->first = t->first->copy;
	  t->second = t->second->copy;
	}
    }
};

static void
unite (int merv[], int n, int a, int b)
{
  if (a == b)
    return;
  if (a > b)
    swap (a, b);

  for (int i = 0; i < n; i++)
    if (merv[i] == b)
      merv[i] = a;
}

static bool
check_winding_number (flfmap &ffs, int wind, const vector<int> *flfaces = NULL)
{
  int n57 = 0, n6 = 0;
  for (flfmap::iterator i = ffs.begin (); i != ffs.end (); ++i)
    {
      int len = flfaces ? (*flfaces)[i->first] : i->first;

      if (len == 5 || len == 7)
	n57 += i->second;
      else if (len == 6)
	n6 += i->second;
      else
	abort ();
    }

  wind = abs (wind);
  assert (wind % 3 == 0);
  wind /= 3;

  if (wind % 2 != n57 % 2)
    return false;

  return wind <= 2 * n6 + n57;
}

struct surf_graph;
struct exists_coloring_check_wno : public coloring_callback
{
  surf_graph *g;
  edge *without;

  exists_coloring_check_wno (surf_graph &_g, edge *wo)
    {
      g = &_g;
      without = wo;
    }

  virtual bool operator() (unsigned col[]);
};

struct count_colorings_check_wno : public coloring_callback
{
  int ncol;
  surf_graph *g;
  edge *without;

  count_colorings_check_wno (surf_graph &_g, edge *wo)
    {
      ncol = 0;
      g = &_g;
      without = wo;
    }

  virtual bool operator() (unsigned col[]);
};

static int
mincrit_oface (int n5, int n6, int n7)
{
  if (n6 == 0 && n7 == 0)
    {
      if (n5 == 0)
	return 4;
      if (n5 == 1)
	return 5;
      if (n5 == 2)
	return 8;
      if (n5 == 3)
	return 9;
      if (n5 == 4)
	return 10;
      abort ();
    }

  if (n6 == 1 && n7 == 0)
    {
      if (n5 == 0)
	return 6;
      else if (n5 == 1)
	return 9;
      else if (n5 == 2)
	return 10;

      abort ();
    }

  if (n7 == 1)
    {
      if (n5 == 0)
	return 7;
      else if (n5 == 1)
	return 10;

      abort ();
    }

  abort ();
}

/* Further utilities to operate on graphs on surfaces, including in particular
   pasting of a disk_graph into a 2-cell face.  */

struct surf_graph : public graph
{
  surf_graph (vector<int> &code) : graph (code)
    {
    }

  surf_graph () : graph ()
    {
    }

  void merge_vs (edge *a, edge *b)
    {
      a = a->opp;
      edge *x = a->next;
      edge *y = b->prev;
      a->next = b;
      b->prev = a;
      y->next = x;
      x->prev = y;
    }

  void suppress_bigon (edge *a, int merv[])
    {
      edge *b = a->face_next ();
      edge *ao = a->opp;
      edge *bo = b->opp;
      int u = merv[b->to];
      int v = merv[a->to];

      if (vs[u].es == a)
	vs[u].es = a->next;
      if (vs[v].es == b)
	vs[v].es = b->next;

      a->unlink ();
      b->unlink ();
      ao->opp = bo;
      bo->opp = ao;

      vs[u].deg--;
      vs[v].deg--;

      delete a;
      delete b;
    }

  surf_graph (int n) : graph (n)
    {
    }

  surf_graph (surf_graph &g, const eref_updater &update = eref_updater (), bool unlock = true) : graph (g, update, unlock)
    {
    }

  surf_graph (disk_graph &g, vector<mergedsc> &to_merge) : graph (g, mvec_updater (to_merge))
    {
      for (edge_iter_face ae(fs[0]); !ae.end_p (); ++ae)
	(*ae)->cid = 0;

      int nout = g.nout;
      int merv[nout];
      for (int v = 0; v < nout; v++)
	merv[v] = v;

      for (vector<mergedsc>::iterator m = to_merge.begin (); m != to_merge.end (); ++m)
	{
	  edge *a = m->first;
	  edge *b = m->second;
	  a->cid = b->cid = 1;

	  unite (merv, nout, merv[a->to], merv[b->opp->to]);
	  unite (merv, nout, merv[b->to], merv[a->opp->to]);
	}
      for (int v = 0; v < nout; v++)
       	if (merv[v] != v)
	  vs[v].es = NULL;

      vector<edge *> rem_edges;
      for (edge_iter_face ae(fs[0]); !ae.end_p (); ++ae)
	if (!(*ae)->cid)
	  rem_edges.push_back (*ae);
      fs[0].es = NULL;

      for (vector<mergedsc>::iterator m = to_merge.begin (); m != to_merge.end (); ++m)
	{
	  edge *a = m->first;
	  edge *b = m->second;

	  merge_vs (a, b);
	  merge_vs (b, a);
	  suppress_bigon (a, merv);
	}
      for (int v = 0; v < nout; v++)
	{
	  if (merv[v] != v)
	    continue;

	  for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	    (*ae)->opp->to = v;

	  vs[v].update_degree ();
	}

      for (vector<edge *>::iterator ri = rem_edges.begin (); ri != rem_edges.end (); ++ri)
	{
	  edge *e = *ri;
	  if (e->left != 0)
	    continue;

	  int fno = nfs ();
	  fs.push_back (face (e));

	  for (edge_iter_face ae(fs[fno]); !ae.end_p (); ++ae)
	    (*ae)->left = fno;
	  fs[fno].update_length ();
	}
    }

  bool check_winding_numbers (unsigned color[], edge *without, const vector<int> *flfaces = NULL)
    {
      for (int f = 0; f < nfs (); f++)
	{
	  if (!fs[f].es)
	    continue;

	  if (without
	      && (without->left == f || without->opp->left == f))
	    continue;

	  int w = winding_number (f, color, NULL);

	  if (!check_winding_number (fs[f].floating_faces, w, flfaces))
	    return false;
	}

      if (without == NULL)
	return true;

      int l = without->left;
      int r = without->opp->left;

      if (l == r)
	abort ();

      int w = winding_number (l, color, without) + winding_number (r, color, without);
      flfmap ffs = fs[l].floating_faces;
      flfmap &rg = fs[r].floating_faces;
      add_flfmap (ffs, rg);

      return check_winding_number (ffs, w, flfaces);
    }

  bool is_colorable_wno (edge *without = NULL)
    {
      int nv = nvs ();
      unsigned color[nv];
      exists_coloring_check_wno cb(*this, without);

      return find_coloring (color, cb, without);
    }

  int number_of_colorings (void)
    {
      int nv = nvs ();
      unsigned color[nv];
      count_colorings_check_wno nc(*this, NULL);
      find_coloring (color, nc);
      return nc.ncol;
    }

  bool find_leaf (int &v)
    {
      for (v = 0; v < nvs (); v++)
	{
	  edge *e = vs[v].es;
	  if (e && e->next == e)
	    return true;
	}

      return false;
    }

  void eliminate_leaves (void)
    {
      int v;

      while (find_leaf (v))
	remove_edge (vs[v].es);
    }

  edge *redundant_edge_at (int v)
    {
      if (!vs[v].es)
	return NULL;

      for (edge_iter_nbr ae(vs[v]); !ae.end_p (); ++ae)
	{
	  edge *e = *ae;

	  if (e->left == e->opp->left)
	    continue;

	  if (!is_colorable_wno (*ae))
	    return *ae;
	}

      return NULL;
    }

  bool is_minimal (void)
    {
      for (int v = 0; v < nvs (); v++)
	if (redundant_edge_at (v))
	  return false;

      return true;
    }

  bool minimize ()
    {
      bool changed = false;
      for (int v = 0; v < nvs (); v++)
	{
	  edge *e;
	  while ((e = redundant_edge_at (v)) != NULL)
	    {
	      if (e->left == e->opp->left)
		abort ();
	      remove_edge (e);
	      eliminate_leaves ();
	      changed = true;
	    }
	}
      for (int f = 0; f < nfs (); f++)
	{
	  if (fs[f].es == NULL)
	    continue;
	  int len = fs[f].length ();
	  if (len == 4)
	    {
	      fs[f].floating_faces.clear ();
	      continue;
	    }
	  if (len == 5)
	    {
	      fs[f].floating_faces.clear ();
	      fs[f].add_floating (5);
	      continue;
	    }
  
	  int n57 = 0, n6 = 0;
	  flfmap &ffs = fs[f].floating_faces;
	  for (flfmap::iterator i = ffs.begin (); i != ffs.end (); ++i)
	    if (i->first == 5 || i->first == 7)
	      n57 += i->second;
	    else if (i->first == 6)
	      n6 += i->second;
	    else
	      abort ();

	  int mwind = 2 * n6 + n57;
	  int fwind = len / 3;
	  if (fwind % 2 != len % 2)
	    fwind--;

	  if (mwind >= fwind)
	    {
	      fs[f].floating_faces.clear ();
	      fs[f].add_floating (len);
	    }
	  }

      return changed;
    }

  void decontract_2path (edge *ef, edge *et, int &fcl, int &fcr)
    {
      int v = ef->opp->to;

      edge *bf = et->next;
      edge *bt = ef->prev;
      int w = nvs ();

      vs[v].es = ef;
      vs.push_back (vertex ());
      vs[w].es = bf;

      et->next = ef;
      ef->prev = et;
      bt->next = bf;
      bf->prev = bt;
      for (edge_iter_nbr e(vs[w]); !e.end_p (); ++e)
	(*e)->opp->to = w;

      fcl = ef->left;
      fcr = et->opp->left;

      int z = nvs ();
      vs.push_back (vertex ());

      edge *wz = edge_pair (w, z, fcl, fcr);
      edge *zv = edge_pair (z, v, fcl, fcr);
      wz->link_after (bt);
      vs[z].es = wz->opp->prev = wz->opp->next = wz->opp;
      zv->link_after (wz->opp);
      zv->opp->link_after (et);

      vs[v].update_degree ();
      vs[w].update_degree ();
      vs[z].deg = 2;
      fs[fcl].update_length ();
      fs[fcr].update_length ();
    }

  int vs_inside_cycle (vector<edge *> &cycle, flfmap *fins = NULL)
    {
      clear_marks ();

      for (vector<edge *>::iterator e = cycle.begin (); e != cycle.end (); ++e)
	(*e)->cid = 0;

      for (vector<edge *>::iterator e = cycle.begin (); e != cycle.end (); ++e)
	mark_reachable_faces ((*e)->left);

      for (vector<edge *>::iterator e = cycle.begin (); e != cycle.end (); ++e)
	if (fs[(*e)->opp->left].seen)
	  return -1;

      int faces_seen = 0;
      int edges_seen = cycle.size ();
      for (int f = 0; f < nfs (); f++)
	if (fs[f].seen)
	  {
	    faces_seen++;
	    for (edge_iter_face e(fs[f]); !e.end_p (); ++e)
	      {
		edges_seen++;
		vs[(*e)->to].seen = true;
	      }
	  }
      if (edges_seen % 2 != 0)
	abort ();
      edges_seen /= 2;

      int vertices_seen = 0, tot = 0;
      for (int v = 0; v < nvs (); v++)
	{
	  if (!vs[v].es)
	    continue;
	  tot++;
	  if (vs[v].seen)
	    vertices_seen++;
	}

      int genus_inside = edges_seen + 2 - vertices_seen - (faces_seen + 1);
      if (genus_inside == 0)
	{
	  if (fins)
	    for (int f = 0; f < nfs (); f++)
      	      if (fs[f].seen)
		add_flfmap (*fins, fs[f].floating_faces);
	  return vertices_seen - cycle.size ();
	}
      else if (genus_inside == 2)
	{
	  if (fins)
	    for (int f = 0; f < nfs (); f++)
      	      if (fs[f].es != NULL && !fs[f].seen)
		add_flfmap (*fins, fs[f].floating_faces);
	  return tot - vertices_seen;
	}
      else
	abort ();
    }

  void paste_disk (int fno, disk_graph *dg)
    {
#ifdef _OPENMP
      omp_set_lock (&dg->copy_lock);
#endif
      int vmap[dg->nvs ()];
      int fdelta = nfs () - 1;
      edge *matche[dg->nout];
      edge *e;
      int v;

      v = dg->nout - 1;
      for (edge_iter_face e(fs[fno]); !e.end_p (); ++e)
	{
	  vmap[v] = (*e)->to;
	  matche[v] = *e;
	  v--;
	}

      for (v = 0; v < dg->nvs(); v++)
	{
	  for (edge_iter_nbr ae(dg->vs[v]); !ae.end_p (); ++ae)
	    if ((*ae)->left == 0)
	      (*ae)->copy = matche[(*ae)->opp->to]->opp;
	    else if ((*ae)->opp->left == 0)
	      (*ae)->copy = matche[(*ae)->to];
	    else
	      (*ae)->copy = new edge;
	}

      for (v = dg->nout; v < dg->nvs(); v++)
	{
      	  vmap[v] = nvs ();
	  vs.push_back (vertex ());
	  vs[vmap[v]].es = dg->vs[v].es->copy;
	}

      for (int f = 1; f < dg->nfs (); f++)
	{
	  face nwf = face (dg->fs[f].es->copy);
	  nwf.floating_faces = dg->fs[f].floating_faces;
	  fs.push_back (nwf);
	}

      for (v = 0; v < dg->nvs(); v++)
	{
	  for (edge_iter_nbr ae(dg->vs[v]); !ae.end_p (); ++ae)
	    {
	      e = *ae;
	      edge *cp = e->copy;

	      if (e->left == 0)
		cp->next = e->next->copy;
	      else if (e->opp->left == 0)
		{
		  cp->left = e->left + fdelta;
		  cp->prev = e->prev->copy;
		}
	      else
		{
		  cp->to = vmap[e->to];
		  cp->left = e->left + fdelta;
		  cp->opp = e->opp->copy;
		  cp->next = e->next->copy;
		  cp->prev = e->prev->copy;
		}
	    }
	}

      for (v = 0; v < dg->nvs (); v++)
	vs[vmap[v]].update_degree ();
      fs[fno].es = NULL;
      for (int f = 1; f < dg->nfs (); f++)
	fs[nfs () - f].update_length ();

#ifdef _OPENMP
      omp_unset_lock (&dg->copy_lock);
#endif
    }

  void paste_disks (disk_graph *ds[])
    {
      int orig_nf = nfs ();
      for (int f = 0; f < orig_nf; f++)
	if (ds[f])
	  paste_disk (f, ds[f]);
    }

  bool is_colorable (edge *without = NULL)
    {
      int nv = nvs ();
      unsigned color[nv];
      coloring_callback cb;

      return find_coloring (color, cb, without);
    }

  bool critical ()
    {
      if (is_colorable ())
	return false;

      for (int v = 0; v < nvs (); v++)
	for (edge_iter_nbr e(vs[v]); !e.end_p (); ++e)
	  {
	    if ((*e)->to < (*e)->opp->to)
	      continue;
	    if (!is_colorable (*e))
	      return false;
	  }

      return true;
    }
};

bool exists_coloring_check_wno::operator() (unsigned col[])
{
  return g->check_winding_numbers (col, without);
}

bool count_colorings_check_wno::operator() (unsigned col[])
{
  if (g->check_winding_numbers (col, without))
    ncol++;
  return false;
}

static bool facial (vector<edge *> &cycle)
{
  int i, len = cycle.size ();

  for (i = 0; i < len; i++)
    if (cycle[(i + 1) % len] != cycle[i]->face_next ())
      break;

  if (i == len)
    return true;

  for (i = 0; i < len; i++)
    if (cycle[(i + 1) % len] != cycle[i]->opp->prev)
      break;

  if (i == len)
    return true;

  return false;
}

struct code_type
{
  vector<unsigned> code;
  int last_bits_used;

  void add (unsigned x, unsigned width)
    {
      if (last_bits_used == 0)
	code.push_back (0);

      unsigned rem = 32 - last_bits_used;
      unsigned &lst = code.back ();

      if (width <= rem)
	{
	  last_bits_used = (last_bits_used + width) % 32;
	  lst = (lst << width) | x;
	  return;
	}

      unsigned xl = x & ((1 << rem) - 1);
      lst = (lst << rem) | xl;
      code.push_back (x >> rem);
      last_bits_used = width - rem;
    }

  code_type (vector<int> &c)
    {
      unsigned width = 3;
      last_bits_used = 0;

      for (vector<int>::iterator a = c.begin (); a != c.end (); ++a)
	{
	  unsigned x = *a + MAX_FLOATING - 3;
	  add (x, width);
	  if (x > (1u << width) - 1)
	    abort ();
	  if (x == (1u << width) - 1)
	    width++;
	}
    }

  bool operator== (const code_type& c) const
    {
      return last_bits_used == c.last_bits_used && code == c.code;
    }
};

struct code_type_hash
{
  size_t operator() (const code_type& c) const
    {
      unsigned ret = c.last_bits_used;

      for (vector<unsigned>::const_iterator i = c.code.begin (); i != c.code.end (); ++i)
	ret = 181 * ret + *i;

      return ret;
    }
};

typedef unordered_set<code_type,code_type_hash> rct;

static bool
register_cd (rct &where, vector<int> &code)
{
  code_type ct(code);

  if (where.count (ct))
    return true;

  where.insert (ct);
  return false;
}

static bool
read_code (FILE *f, vector<int> &code)
{
  code.clear ();
  char *buf = NULL, *p;
  size_t sze = 0;
  int n, a;

  if (getline (&buf, &sze, f) < 0)
    {
      free (buf);
      return false;
    }
  p = buf;
  while (sscanf (p, "%d%n", &a, &n) > 0)
    {
      code.push_back (a);
      p += n;
    }

  free (buf);

  return !code.empty ();
}
