/// \ingroup auxdat

   ///

   /// \brief A \e Union-Find data structure implementation which

   /// is able to store a priority for each item and retrieve the minimum of

   /// each class.

   ///

   /// A \e Union-Find data structure implementation which is able to

   /// store a priority for each item and retrieve the minimum of each

   /// class. In addition, it supports the joining and splitting the

   /// components. If you don't need this feature then you makes

   /// better to use the \ref UnionFind class which is more efficient.

   ///

   /// The union-find data strcuture based on a (2, 16)-tree with a

   /// tournament minimum selection on the internal nodes. The insert

   /// operation takes O(1), the find, set, decrease and increase takes

   /// O(log(n)), where n is the number of nodes in the current

   /// component.  The complexity of join and split is O(log(n)*k),

   /// where n is the sum of the number of the nodes and k is the

   /// number of joined components or the number of the components

   /// after the split.

   ///

   /// \pre You need to add all the elements by the \ref insert()

   /// method.

   ///

   template <typename _Value, typename _ItemIntMap, 

             typename _Comp = std::less<_Value> >

   class HeapUnionFind {

   public:

     typedef _Value Value;

     typedef typename _ItemIntMap::Key Item;

     typedef _ItemIntMap ItemIntMap;

     typedef _Comp Comp;

   private:

     static const int cmax = 3;

     ItemIntMap& index;

     struct ClassNode {

       int parent;

       int depth;

       int left, right;

       int next, prev;

     };

     int first_class;

     int first_free_class;

     std::vector<ClassNode> classes;

     int newClass() {

       if (first_free_class < 0) {

         int id = classes.size();

         classes.push_back(ClassNode());

         return id;

       } else {

         int id = first_free_class;

         first_free_class = classes[id].next;

         return id;

       }

     }

     void deleteClass(int id) {

       classes[id].next = first_free_class;

       first_free_class = id;

     }

     struct ItemNode {

       int parent;

       Item item;

       Value prio;

       int next, prev;

       int left, right;

       int size;

     };

     int first_free_node;

     std::vector<ItemNode> nodes;

     int newNode() {

       if (first_free_node < 0) {

         int id = nodes.size();

         nodes.push_back(ItemNode());

         return id;

       } else {

         int id = first_free_node;

         first_free_node = nodes[id].next;

         return id;

       }

     }

     void deleteNode(int id) {

       nodes[id].next = first_free_node;

       first_free_node = id;

     }

     Comp comp;

     int findClass(int id) const {

       int kd = id;

       while (kd >= 0) {

         kd = nodes[kd].parent;

       }

       return ~kd;

     }

     int leftNode(int id) const {

       int kd = ~(classes[id].parent);

       for (int i = 0; i < classes[id].depth; ++i) {

         kd = nodes[kd].left;

       }

       return kd;

     }

     int nextNode(int id) const {

       int depth = 0;

       while (id >= 0 && nodes[id].next == -1) {

         id = nodes[id].parent;

         ++depth;

       }

       if (id < 0) {

         return -1;

       }

       id = nodes[id].next;

       while (depth--) {

         id = nodes[id].left;      

       }

       return id;

     }

     void setPrio(int id) {

       int jd = nodes[id].left;

       nodes[id].prio = nodes[jd].prio;

       nodes[id].item = nodes[jd].item;

       jd = nodes[jd].next;

       while (jd != -1) {

         if (comp(nodes[jd].prio, nodes[id].prio)) {

           nodes[id].prio = nodes[jd].prio;

           nodes[id].item = nodes[jd].item;

         }

         jd = nodes[jd].next;

       }

     }

     void push(int id, int jd) {

       nodes[id].size = 1;

       nodes[id].left = nodes[id].right = jd;

       nodes[jd].next = nodes[jd].prev = -1;

       nodes[jd].parent = id;

     }

     void pushAfter(int id, int jd) {

       int kd = nodes[id].parent;

       if (nodes[id].next != -1) {

         nodes[nodes[id].next].prev = jd;

         if (kd >= 0) {

           nodes[kd].size += 1;

         }

       } else {

         if (kd >= 0) {

           nodes[kd].right = jd;

           nodes[kd].size += 1;

         }

       }

       nodes[jd].next = nodes[id].next;

       nodes[jd].prev = id;

       nodes[id].next = jd;

       nodes[jd].parent = kd;

     }

     void pushRight(int id, int jd) {

       nodes[id].size += 1;

       nodes[jd].prev = nodes[id].right;

       nodes[jd].next = -1;

       nodes[nodes[id].right].next = jd;

       nodes[id].right = jd;

       nodes[jd].parent = id;

     }

     void popRight(int id) {

       nodes[id].size -= 1;

       int jd = nodes[id].right;

       nodes[nodes[jd].prev].next = -1;

       nodes[id].right = nodes[jd].prev;

     }

     void splice(int id, int jd) {

       nodes[id].size += nodes[jd].size;

       nodes[nodes[id].right].next = nodes[jd].left;

       nodes[nodes[jd].left].prev = nodes[id].right;

       int kd = nodes[jd].left;

       while (kd != -1) {

         nodes[kd].parent = id;

         kd = nodes[kd].next;

       }

     }

     void split(int id, int jd) {

       int kd = nodes[id].parent;

       nodes[kd].right = nodes[id].prev;

       nodes[nodes[id].prev].next = -1;

       nodes[jd].left = id;

       nodes[id].prev = -1;

       int num = 0;

       while (id != -1) {

         nodes[id].parent = jd;

         nodes[jd].right = id;

         id = nodes[id].next;

         ++num;

       }      

       nodes[kd].size -= num;

       nodes[jd].size = num;

     }

     void pushLeft(int id, int jd) {

       nodes[id].size += 1;

       nodes[jd].next = nodes[id].left;

       nodes[jd].prev = -1;

       nodes[nodes[id].left].prev = jd;

       nodes[id].left = jd;

       nodes[jd].parent = id;

     }

     void popLeft(int id) {

       nodes[id].size -= 1;

       int jd = nodes[id].left;

       nodes[nodes[jd].next].prev = -1;

       nodes[id].left = nodes[jd].next;

     }

     void repairLeft(int id) {

       int jd = ~(classes[id].parent);

       while (nodes[jd].left != -1) {

 	int kd = nodes[jd].left;

 	if (nodes[jd].size == 1) {

 	  if (nodes[jd].parent < 0) {

 	    classes[id].parent = ~kd;

 	    classes[id].depth -= 1;

 	    nodes[kd].parent = ~id;

 	    deleteNode(jd);

 	    jd = kd;

 	  } else {

 	    int pd = nodes[jd].parent;

 	    if (nodes[nodes[jd].next].size < cmax) {

 	      pushLeft(nodes[jd].next, nodes[jd].left);

 	      if (less(nodes[jd].left, nodes[jd].next)) {

 		nodes[nodes[jd].next].prio = nodes[nodes[jd].left].prio;

 		nodes[nodes[jd].next].item = nodes[nodes[jd].left].item;

 	      } 

 	      popLeft(pd);

 	      deleteNode(jd);

 	      jd = pd;

 	    } else {

 	      int ld = nodes[nodes[jd].next].left;

 	      popLeft(nodes[jd].next);

 	      pushRight(jd, ld);

 	      if (less(ld, nodes[jd].left)) {

 		nodes[jd].item = nodes[ld].item;

 		nodes[jd].prio = nodes[jd].prio;

 	      }

 	      if (nodes[nodes[jd].next].item == nodes[ld].item) {

 		setPrio(nodes[jd].next);

 	      }

 	      jd = nodes[jd].left;

 	    }

 	  }

 	} else {

 	  jd = nodes[jd].left;

 	}

       }

     }    

     void repairRight(int id) {

       int jd = ~(classes[id].parent);

       while (nodes[jd].right != -1) {

 	int kd = nodes[jd].right;

 	if (nodes[jd].size == 1) {

 	  if (nodes[jd].parent < 0) {

 	    classes[id].parent = ~kd;

 	    classes[id].depth -= 1;

 	    nodes[kd].parent = ~id;

 	    deleteNode(jd);

 	    jd = kd;

 	  } else {

 	    int pd = nodes[jd].parent;

 	    if (nodes[nodes[jd].prev].size < cmax) {

 	      pushRight(nodes[jd].prev, nodes[jd].right);

 	      if (less(nodes[jd].right, nodes[jd].prev)) {

 		nodes[nodes[jd].prev].prio = nodes[nodes[jd].right].prio;

 		nodes[nodes[jd].prev].item = nodes[nodes[jd].right].item;

 	      } 

 	      popRight(pd);

 	      deleteNode(jd);

 	      jd = pd;

 	    } else {

 	      int ld = nodes[nodes[jd].prev].right;

 	      popRight(nodes[jd].prev);

 	      pushLeft(jd, ld);

 	      if (less(ld, nodes[jd].right)) {

 		nodes[jd].item = nodes[ld].item;

 		nodes[jd].prio = nodes[jd].prio;

 	      }

 	      if (nodes[nodes[jd].prev].item == nodes[ld].item) {

 		setPrio(nodes[jd].prev);

 	      }

 	      jd = nodes[jd].right;

 	    }

 	  }

 	} else {

 	  jd = nodes[jd].right;

 	}

       }

     }

     bool less(int id, int jd) const {

       return comp(nodes[id].prio, nodes[jd].prio);

     }

     bool equal(int id, int jd) const {

       return !less(id, jd) && !less(jd, id);

     }

   public:

     /// \brief Returns true when the given class is alive.

     ///

     /// Returns true when the given class is alive, ie. the class is

     /// not nested into other class.

     bool alive(int cls) const {

       return classes[cls].parent < 0;

     }

     /// \brief Returns true when the given class is trivial.

     ///

     /// Returns true when the given class is trivial, ie. the class

     /// contains just one item directly.

     bool trivial(int cls) const {

       return classes[cls].left == -1;

     }

     /// \brief Constructs the union-find.

     ///

     /// Constructs the union-find.  

     /// \brief _index The index map of the union-find. The data

     /// structure uses internally for store references.

     HeapUnionFind(ItemIntMap& _index) 

       : index(_index), first_class(-1), 

 	first_free_class(-1), first_free_node(-1) {}

     /// \brief Insert a new node into a new component.

     ///

     /// Insert a new node into a new component.

     /// \param item The item of the new node.

     /// \param prio The priority of the new node.

     /// \return The class id of the one-item-heap.

     int insert(const Item& item, const Value& prio) {

       int id = newNode();

       nodes[id].item = item;

       nodes[id].prio = prio;

       nodes[id].size = 0;

       nodes[id].prev = -1;

       nodes[id].next = -1;

       nodes[id].left = -1;

       nodes[id].right = -1;

       nodes[id].item = item;

       index[item] = id;

       int class_id = newClass();

       classes[class_id].parent = ~id;

       classes[class_id].depth = 0;

       classes[class_id].left = -1;

       classes[class_id].right = -1;

       if (first_class != -1) {

         classes[first_class].prev = class_id;

       }

       classes[class_id].next = first_class;

       classes[class_id].prev = -1;

       first_class = class_id;

       nodes[id].parent = ~class_id;

       return class_id;

     }

     /// \brief The class of the item.

     ///

     /// \return The alive class id of the item, which is not nested into

     /// other classes.

     ///

     /// The time complexity is O(log(n)).

     int find(const Item& item) const {

       return findClass(index[item]);

     }

     /// \brief Joins the classes.

     ///

     /// The current function joins the given classes. The parameter is

     /// an STL range which should be contains valid class ids. The

     /// time complexity is O(log(n)*k) where n is the overall number

     /// of the joined nodes and k is the number of classes.

     /// \return The class of the joined classes.

     /// \pre The range should contain at least two class ids.

     template <typename Iterator>

     int join(Iterator begin, Iterator end) {

       std::vector<int> cs;

       for (Iterator it = begin; it != end; ++it) {

         cs.push_back(*it);

       }

       int class_id = newClass();

       { // creation union-find

         if (first_class != -1) {

           classes[first_class].prev = class_id;

         }

         classes[class_id].next = first_class;

         classes[class_id].prev = -1;

         first_class = class_id;

         classes[class_id].depth = classes[cs[0]].depth;

         classes[class_id].parent = classes[cs[0]].parent;

         nodes[~(classes[class_id].parent)].parent = ~class_id;

         int l = cs[0];

         classes[class_id].left = l;

         classes[class_id].right = l;

         if (classes[l].next != -1) {

           classes[classes[l].next].prev = classes[l].prev;

         }

         classes[classes[l].prev].next = classes[l].next;

         classes[l].prev = -1;

         classes[l].next = -1;

         classes[l].depth = leftNode(l);

         classes[l].parent = class_id;

       }

       { // merging of heap

         int l = class_id;

         for (int ci = 1; ci < int(cs.size()); ++ci) {

           int r = cs[ci];

           int rln = leftNode(r);

           if (classes[l].depth > classes[r].depth) {

             int id = ~(classes[l].parent);

             for (int i = classes[r].depth + 1; i < classes[l].depth; ++i) {

               id = nodes[id].right;

             }

             while (id >= 0 && nodes[id].size == cmax) {

               int new_id = newNode();

               int right_id = nodes[id].right;

               popRight(id);

               if (nodes[id].item == nodes[right_id].item) {

                 setPrio(id);

               }

               push(new_id, right_id);

               pushRight(new_id, ~(classes[r].parent));

               setPrio(new_id);

               id = nodes[id].parent;

               classes[r].parent = ~new_id;

             }

             if (id < 0) {

               int new_parent = newNode();

               nodes[new_parent].next = -1;

               nodes[new_parent].prev = -1;

               nodes[new_parent].parent = ~l;

               push(new_parent, ~(classes[l].parent));

               pushRight(new_parent, ~(classes[r].parent));

               setPrio(new_parent);

               classes[l].parent = ~new_parent;

               classes[l].depth += 1;

             } else {

               pushRight(id, ~(classes[r].parent));

               while (id >= 0 && less(~(classes[r].parent), id)) {

                 nodes[id].prio = nodes[~(classes[r].parent)].prio;

                 nodes[id].item = nodes[~(classes[r].parent)].item;

                 id = nodes[id].parent;

               }

             }

           } else if (classes[r].depth > classes[l].depth) {

             int id = ~(classes[r].parent);

             for (int i = classes[l].depth + 1; i < classes[r].depth; ++i) {

               id = nodes[id].left;

             }

             while (id >= 0 && nodes[id].size == cmax) {

               int new_id = newNode();

               int left_id = nodes[id].left;

               popLeft(id);

               if (nodes[id].prio == nodes[left_id].prio) {

                 setPrio(id);

               }

               push(new_id, left_id);

               pushLeft(new_id, ~(classes[l].parent));

               setPrio(new_id);

               id = nodes[id].parent;

               classes[l].parent = ~new_id;

             }

             if (id < 0) {

               int new_parent = newNode();

               nodes[new_parent].next = -1;

               nodes[new_parent].prev = -1;

               nodes[new_parent].parent = ~l;

               push(new_parent, ~(classes[r].parent));

               pushLeft(new_parent, ~(classes[l].parent));

               setPrio(new_parent);

               classes[r].parent = ~new_parent;

               classes[r].depth += 1;

             } else {

               pushLeft(id, ~(classes[l].parent));

               while (id >= 0 && less(~(classes[l].parent), id)) {

                 nodes[id].prio = nodes[~(classes[l].parent)].prio;

                 nodes[id].item = nodes[~(classes[l].parent)].item;

                 id = nodes[id].parent;

               }

             }

             nodes[~(classes[r].parent)].parent = ~l;

             classes[l].parent = classes[r].parent;

             classes[l].depth = classes[r].depth;

           } else {

             if (classes[l].depth != 0 && 

                 nodes[~(classes[l].parent)].size + 

                 nodes[~(classes[r].parent)].size <= cmax) {

               splice(~(classes[l].parent), ~(classes[r].parent));

               deleteNode(~(classes[r].parent));

               if (less(~(classes[r].parent), ~(classes[l].parent))) {

                 nodes[~(classes[l].parent)].prio = 

                   nodes[~(classes[r].parent)].prio;

                 nodes[~(classes[l].parent)].item = 

                   nodes[~(classes[r].parent)].item;

               }

             } else {

               int new_parent = newNode();

               nodes[new_parent].next = nodes[new_parent].prev = -1;

               push(new_parent, ~(classes[l].parent));

               pushRight(new_parent, ~(classes[r].parent));

               setPrio(new_parent);

               classes[l].parent = ~new_parent;

               classes[l].depth += 1;

               nodes[new_parent].parent = ~l;

             }

           }

           if (classes[r].next != -1) {

             classes[classes[r].next].prev = classes[r].prev;

           }

           classes[classes[r].prev].next = classes[r].next;

           classes[r].prev = classes[l].right;

           classes[classes[l].right].next = r;

           classes[l].right = r;

           classes[r].parent = l;

           classes[r].next = -1;

           classes[r].depth = rln;

         }

       }

       return class_id;

     }

     /// \brief Split the class to subclasses.

     ///

     /// The current function splits the given class. The join, which

     /// made the current class, stored a reference to the

     /// subclasses. The \c splitClass() member restores the classes

     /// and creates the heaps. The parameter is an STL output iterator

     /// which will be filled with the subclass ids. The time

     /// complexity is O(log(n)*k) where n is the overall number of

     /// nodes in the splitted classes and k is the number of the

     /// classes.

     template <typename Iterator>

     void split(int cls, Iterator out) {

       std::vector<int> cs;

       { // splitting union-find

         int id = cls;

         int l = classes[id].left;

         classes[l].parent = classes[id].parent;

         classes[l].depth = classes[id].depth;

         nodes[~(classes[l].parent)].parent = ~l;

         *out++ = l;

         while (l != -1) {

           cs.push_back(l);

           l = classes[l].next;

         }

         classes[classes[id].right].next = first_class;

         classes[first_class].prev = classes[id].right;

         first_class = classes[id].left;

         if (classes[id].next != -1) {

           classes[classes[id].next].prev = classes[id].prev;

         }

         classes[classes[id].prev].next = classes[id].next;

         deleteClass(id);

       }

       {

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

           int l = classes[cs[i]].depth;

           while (nodes[nodes[l].parent].left == l) {

             l = nodes[l].parent;

           }

           int r = l;    

           while (nodes[l].parent >= 0) {

             l = nodes[l].parent;

             int new_node = newNode();

             nodes[new_node].prev = -1;

             nodes[new_node].next = -1;

             split(r, new_node);

             pushAfter(l, new_node);

             setPrio(l);

             setPrio(new_node);

             r = new_node;

           }

           classes[cs[i]].parent = ~r;

           classes[cs[i]].depth = classes[~(nodes[l].parent)].depth;

           nodes[r].parent = ~cs[i];

           nodes[l].next = -1;

           nodes[r].prev = -1;

           repairRight(~(nodes[l].parent));

           repairLeft(cs[i]);

           *out++ = cs[i];

         }

       }

     }

     /// \brief Gives back the priority of the current item.

     ///

     /// \return Gives back the priority of the current item.

     const Value& operator[](const Item& item) const {

       return nodes[index[item]].prio;

     }

     /// \brief Sets the priority of the current item.

     ///

     /// Sets the priority of the current item.

     void set(const Item& item, const Value& prio) {

       if (comp(prio, nodes[index[item]].prio)) {

         decrease(item, prio);

       } else if (!comp(prio, nodes[index[item]].prio)) {

         increase(item, prio);

       }

     }

     /// \brief Increase the priority of the current item.

     ///

     /// Increase the priority of the current item.

     void increase(const Item& item, const Value& prio) {

       int id = index[item];

       int kd = nodes[id].parent;

       nodes[id].prio = prio;

       while (kd >= 0 && nodes[kd].item == item) {

         setPrio(kd);

         kd = nodes[kd].parent;

       }

     }

     /// \brief Increase the priority of the current item.

     ///

     /// Increase the priority of the current item.

     void decrease(const Item& item, const Value& prio) {

       int id = index[item];

       int kd = nodes[id].parent;

       nodes[id].prio = prio;

       while (kd >= 0 && less(id, kd)) {

         nodes[kd].prio = prio;

         nodes[kd].item = item;

         kd = nodes[kd].parent;

       }

     }

     /// \brief Gives back the minimum priority of the class.

     ///

     /// \return Gives back the minimum priority of the class.

     const Value& classPrio(int cls) const {

       return nodes[~(classes[cls].parent)].prio;

     }

     /// \brief Gives back the minimum priority item of the class.

     ///

     /// \return Gives back the minimum priority item of the class.

     const Item& classTop(int cls) const {

       return nodes[~(classes[cls].parent)].item;

     }

     /// \brief Gives back a representant item of the class.

     /// 

     /// The representant is indpendent from the priorities of the

     /// items. 

     /// \return Gives back a representant item of the class.

     const Item& classRep(int id) const {

       int parent = classes[id].parent;

       return nodes[parent >= 0 ? classes[id].depth : leftNode(id)].item;

     }

     /// \brief Lemon style iterator for the items of a class.

     ///

     /// ClassIt is a lemon style iterator for the components. It iterates

     /// on the items of a class. By example if you want to iterate on

     /// each items of each classes then you may write the next code.

     ///\code

     /// for (ClassIt cit(huf); cit != INVALID; ++cit) {

     ///   std::cout << "Class: ";

     ///   for (ItemIt iit(huf, cit); iit != INVALID; ++iit) {

     ///     std::cout << toString(iit) << ' ' << std::endl;

     ///   }

     ///   std::cout << std::endl;

     /// }

     ///\endcode

     class ItemIt {

     private:

       const HeapUnionFind* _huf;

       int _id, _lid;

     public:

       /// \brief Default constructor 

       ///

       /// Default constructor 

       ItemIt() {}

       ItemIt(const HeapUnionFind& huf, int cls) : _huf(&huf) {

         int id = cls;

         int parent = _huf->classes[id].parent;

         if (parent >= 0) {

           _id = _huf->classes[id].depth;

           if (_huf->classes[id].next != -1) {

             _lid = _huf->classes[_huf->classes[id].next].depth;

           } else {

             _lid = -1;

           }

         } else {

           _id = _huf->leftNode(id);

           _lid = -1;

         } 

       }

       /// \brief Increment operator

       ///

       /// It steps to the next item in the class.

       ItemIt& operator++() {

         _id = _huf->nextNode(_id);

         return *this;

       }

       /// \brief Conversion operator

       ///

       /// It converts the iterator to the current item.

       operator const Item&() const {

         return _huf->nodes[_id].item;

       }

       /// \brief Equality operator

       ///

       /// Equality operator

       bool operator==(const ItemIt& i) { 

         return i._id == _id;

       }

       /// \brief Inequality operator

       ///

       /// Inequality operator

       bool operator!=(const ItemIt& i) { 

         return i._id != _id;

       }

       /// \brief Equality operator

       ///

       /// Equality operator

       bool operator==(Invalid) { 

         return _id == _lid;

       }

       /// \brief Inequality operator

       ///

       /// Inequality operator

       bool operator!=(Invalid) { 

         return _id != _lid;

       }

     };

     /// \brief Class iterator

     ///

     /// The iterator stores 

     class ClassIt {

     private:

       const HeapUnionFind* _huf;

       int _id;

     public:

       ClassIt(const HeapUnionFind& huf) 

         : _huf(&huf), _id(huf.first_class) {}

       ClassIt(const HeapUnionFind& huf, int cls) 

         : _huf(&huf), _id(huf.classes[cls].left) {}

       ClassIt(Invalid) : _huf(0), _id(-1) {}

       const ClassIt& operator++() {

         _id = _huf->classes[_id].next;

 	return *this;

       }

       /// \brief Equality operator

       ///

       /// Equality operator

       bool operator==(const ClassIt& i) { 

         return i._id == _id;

       }

       /// \brief Inequality operator

       ///

       /// Inequality operator

       bool operator!=(const ClassIt& i) { 

         return i._id != _id;

       }      

       operator int() const {

 	return _id;

       }

     };

   };