lemon/unionfind.h
author Alpar Juttner <alpar@cs.elte.hu>
Thu, 24 Apr 2008 20:26:14 +0100
changeset 151 4d79daa40e9b
child 209 765619b7cbb2
permissions -rw-r--r--
Turn on built in Doxygen STL support
alpar@103
     1
/* -*- C++ -*-
alpar@103
     2
 *
alpar@103
     3
 * This file is a part of LEMON, a generic C++ optimization library
alpar@103
     4
 *
alpar@103
     5
 * Copyright (C) 2003-2008
alpar@103
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
alpar@103
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
alpar@103
     8
 *
alpar@103
     9
 * Permission to use, modify and distribute this software is granted
alpar@103
    10
 * provided that this copyright notice appears in all copies. For
alpar@103
    11
 * precise terms see the accompanying LICENSE file.
alpar@103
    12
 *
alpar@103
    13
 * This software is provided "AS IS" with no warranty of any kind,
alpar@103
    14
 * express or implied, and with no claim as to its suitability for any
alpar@103
    15
 * purpose.
alpar@103
    16
 *
alpar@103
    17
 */
alpar@103
    18
alpar@103
    19
#ifndef LEMON_UNION_FIND_H
alpar@103
    20
#define LEMON_UNION_FIND_H
alpar@103
    21
alpar@103
    22
//!\ingroup auxdat
alpar@103
    23
//!\file
alpar@103
    24
//!\brief Union-Find data structures.
alpar@103
    25
//!
alpar@103
    26
alpar@103
    27
#include <vector>
alpar@103
    28
#include <list>
alpar@103
    29
#include <utility>
alpar@103
    30
#include <algorithm>
alpar@103
    31
#include <functional>
alpar@103
    32
alpar@103
    33
#include <lemon/bits/invalid.h>
alpar@103
    34
alpar@103
    35
namespace lemon {
alpar@103
    36
alpar@103
    37
  /// \ingroup auxdat
alpar@103
    38
  ///
alpar@103
    39
  /// \brief A \e Union-Find data structure implementation
alpar@103
    40
  ///
alpar@103
    41
  /// The class implements the \e Union-Find data structure. 
alpar@103
    42
  /// The union operation uses rank heuristic, while
alpar@103
    43
  /// the find operation uses path compression.
alpar@103
    44
  /// This is a very simple but efficient implementation, providing 
alpar@103
    45
  /// only four methods: join (union), find, insert and size.
alpar@103
    46
  /// For more features see the \ref UnionFindEnum class.
alpar@103
    47
  ///
alpar@103
    48
  /// It is primarily used in Kruskal algorithm for finding minimal
alpar@103
    49
  /// cost spanning tree in a graph.
alpar@103
    50
  /// \sa kruskal()
alpar@103
    51
  ///
alpar@103
    52
  /// \pre You need to add all the elements by the \ref insert()
alpar@103
    53
  /// method.  
alpar@103
    54
  template <typename _ItemIntMap>
alpar@103
    55
  class UnionFind {
alpar@103
    56
  public:
alpar@103
    57
alpar@103
    58
    typedef _ItemIntMap ItemIntMap;
alpar@103
    59
    typedef typename ItemIntMap::Key Item;
alpar@103
    60
alpar@103
    61
  private:
alpar@103
    62
    // If the items vector stores negative value for an item then
alpar@103
    63
    // that item is root item and it has -items[it] component size.
alpar@103
    64
    // Else the items[it] contains the index of the parent.
alpar@103
    65
    std::vector<int> items;
alpar@103
    66
    ItemIntMap& index;
alpar@103
    67
alpar@103
    68
    bool rep(int idx) const {
alpar@103
    69
      return items[idx] < 0;
alpar@103
    70
    }
alpar@103
    71
alpar@103
    72
    int repIndex(int idx) const {
alpar@103
    73
      int k = idx;
alpar@103
    74
      while (!rep(k)) {
alpar@103
    75
        k = items[k] ;
alpar@103
    76
      }
alpar@103
    77
      while (idx != k) {
alpar@103
    78
        int next = items[idx];
alpar@103
    79
        const_cast<int&>(items[idx]) = k;
alpar@103
    80
        idx = next;
alpar@103
    81
      }
alpar@103
    82
      return k;
alpar@103
    83
    }
alpar@103
    84
alpar@103
    85
  public:
alpar@103
    86
alpar@103
    87
    /// \brief Constructor
alpar@103
    88
    ///
alpar@103
    89
    /// Constructor of the UnionFind class. You should give an item to
alpar@103
    90
    /// integer map which will be used from the data structure. If you
alpar@103
    91
    /// modify directly this map that may cause segmentation fault,
alpar@103
    92
    /// invalid data structure, or infinite loop when you use again
alpar@103
    93
    /// the union-find.
alpar@103
    94
    UnionFind(ItemIntMap& m) : index(m) {}
alpar@103
    95
alpar@103
    96
    /// \brief Returns the index of the element's component.
alpar@103
    97
    ///
alpar@103
    98
    /// The method returns the index of the element's component.
alpar@103
    99
    /// This is an integer between zero and the number of inserted elements.
alpar@103
   100
    ///
alpar@103
   101
    int find(const Item& a) {
alpar@103
   102
      return repIndex(index[a]);
alpar@103
   103
    }
alpar@103
   104
alpar@103
   105
    /// \brief Clears the union-find data structure
alpar@103
   106
    ///
alpar@103
   107
    /// Erase each item from the data structure.
alpar@103
   108
    void clear() {
alpar@103
   109
      items.clear();
alpar@103
   110
    }
alpar@103
   111
alpar@103
   112
    /// \brief Inserts a new element into the structure.
alpar@103
   113
    ///
alpar@103
   114
    /// This method inserts a new element into the data structure. 
alpar@103
   115
    ///
alpar@103
   116
    /// The method returns the index of the new component.
alpar@103
   117
    int insert(const Item& a) {
alpar@103
   118
      int n = items.size();
alpar@103
   119
      items.push_back(-1);
alpar@103
   120
      index.set(a,n);
alpar@103
   121
      return n;
alpar@103
   122
    }
alpar@103
   123
alpar@103
   124
    /// \brief Joining the components of element \e a and element \e b.
alpar@103
   125
    ///
alpar@103
   126
    /// This is the \e union operation of the Union-Find structure. 
alpar@103
   127
    /// Joins the component of element \e a and component of
alpar@103
   128
    /// element \e b. If \e a and \e b are in the same component then
alpar@103
   129
    /// it returns false otherwise it returns true.
alpar@103
   130
    bool join(const Item& a, const Item& b) {
alpar@103
   131
      int ka = repIndex(index[a]);
alpar@103
   132
      int kb = repIndex(index[b]);
alpar@103
   133
alpar@103
   134
      if ( ka == kb ) 
alpar@103
   135
	return false;
alpar@103
   136
alpar@103
   137
      if (items[ka] < items[kb]) {
alpar@103
   138
	items[ka] += items[kb];
alpar@103
   139
	items[kb] = ka;
alpar@103
   140
      } else {
alpar@103
   141
	items[kb] += items[ka];
alpar@103
   142
	items[ka] = kb;
alpar@103
   143
      }
alpar@103
   144
      return true;
alpar@103
   145
    }
alpar@103
   146
alpar@103
   147
    /// \brief Returns the size of the component of element \e a.
alpar@103
   148
    ///
alpar@103
   149
    /// Returns the size of the component of element \e a.
alpar@103
   150
    int size(const Item& a) {
alpar@103
   151
      int k = repIndex(index[a]);
alpar@103
   152
      return - items[k];
alpar@103
   153
    }
alpar@103
   154
alpar@103
   155
  };
alpar@103
   156
alpar@103
   157
  /// \ingroup auxdat
alpar@103
   158
  ///
alpar@103
   159
  /// \brief A \e Union-Find data structure implementation which
alpar@103
   160
  /// is able to enumerate the components.
alpar@103
   161
  ///
alpar@103
   162
  /// The class implements a \e Union-Find data structure
alpar@103
   163
  /// which is able to enumerate the components and the items in
alpar@103
   164
  /// a component. If you don't need this feature then perhaps it's
alpar@103
   165
  /// better to use the \ref UnionFind class which is more efficient.
alpar@103
   166
  ///
alpar@103
   167
  /// The union operation uses rank heuristic, while
alpar@103
   168
  /// the find operation uses path compression.
alpar@103
   169
  ///
alpar@103
   170
  /// \pre You need to add all the elements by the \ref insert()
alpar@103
   171
  /// method.
alpar@103
   172
  ///
alpar@103
   173
  template <typename _ItemIntMap>
alpar@103
   174
  class UnionFindEnum {
alpar@103
   175
  public:
alpar@103
   176
    
alpar@103
   177
    typedef _ItemIntMap ItemIntMap;
alpar@103
   178
    typedef typename ItemIntMap::Key Item;
alpar@103
   179
alpar@103
   180
  private:
alpar@103
   181
    
alpar@103
   182
    ItemIntMap& index;
alpar@103
   183
alpar@103
   184
    // If the parent stores negative value for an item then that item
alpar@103
   185
    // is root item and it has ~(items[it].parent) component id.  Else
alpar@103
   186
    // the items[it].parent contains the index of the parent.
alpar@103
   187
    //
alpar@103
   188
    // The \c next and \c prev provides the double-linked
alpar@103
   189
    // cyclic list of one component's items.
alpar@103
   190
    struct ItemT {
alpar@103
   191
      int parent;
alpar@103
   192
      Item item;
alpar@103
   193
alpar@103
   194
      int next, prev;
alpar@103
   195
    };
alpar@103
   196
alpar@103
   197
    std::vector<ItemT> items;
alpar@103
   198
    int firstFreeItem;
alpar@103
   199
alpar@103
   200
    struct ClassT {
alpar@103
   201
      int size;
alpar@103
   202
      int firstItem;
alpar@103
   203
      int next, prev;
alpar@103
   204
    };
alpar@103
   205
    
alpar@103
   206
    std::vector<ClassT> classes;
alpar@103
   207
    int firstClass, firstFreeClass;
alpar@103
   208
alpar@103
   209
    int newClass() {
alpar@103
   210
      if (firstFreeClass == -1) {
alpar@103
   211
	int cdx = classes.size();
alpar@103
   212
	classes.push_back(ClassT());
alpar@103
   213
	return cdx;
alpar@103
   214
      } else {
alpar@103
   215
	int cdx = firstFreeClass;
alpar@103
   216
	firstFreeClass = classes[firstFreeClass].next;
alpar@103
   217
	return cdx;
alpar@103
   218
      }
alpar@103
   219
    }
alpar@103
   220
alpar@103
   221
    int newItem() {
alpar@103
   222
      if (firstFreeItem == -1) {
alpar@103
   223
	int idx = items.size();
alpar@103
   224
	items.push_back(ItemT());
alpar@103
   225
	return idx;
alpar@103
   226
      } else {
alpar@103
   227
	int idx = firstFreeItem;
alpar@103
   228
	firstFreeItem = items[firstFreeItem].next;
alpar@103
   229
	return idx;
alpar@103
   230
      }
alpar@103
   231
    }
alpar@103
   232
alpar@103
   233
alpar@103
   234
    bool rep(int idx) const {
alpar@103
   235
      return items[idx].parent < 0;
alpar@103
   236
    }
alpar@103
   237
alpar@103
   238
    int repIndex(int idx) const {
alpar@103
   239
      int k = idx;
alpar@103
   240
      while (!rep(k)) {
alpar@103
   241
        k = items[k].parent;
alpar@103
   242
      }
alpar@103
   243
      while (idx != k) {
alpar@103
   244
        int next = items[idx].parent;
alpar@103
   245
        const_cast<int&>(items[idx].parent) = k;
alpar@103
   246
        idx = next;
alpar@103
   247
      }
alpar@103
   248
      return k;
alpar@103
   249
    }
alpar@103
   250
alpar@103
   251
    int classIndex(int idx) const {
alpar@103
   252
      return ~(items[repIndex(idx)].parent);
alpar@103
   253
    }
alpar@103
   254
alpar@103
   255
    void singletonItem(int idx) {
alpar@103
   256
      items[idx].next = idx;
alpar@103
   257
      items[idx].prev = idx;
alpar@103
   258
    }
alpar@103
   259
alpar@103
   260
    void laceItem(int idx, int rdx) {
alpar@103
   261
      items[idx].prev = rdx;
alpar@103
   262
      items[idx].next = items[rdx].next;
alpar@103
   263
      items[items[rdx].next].prev = idx;
alpar@103
   264
      items[rdx].next = idx;
alpar@103
   265
    }
alpar@103
   266
alpar@103
   267
    void unlaceItem(int idx) {
alpar@103
   268
      items[items[idx].prev].next = items[idx].next;
alpar@103
   269
      items[items[idx].next].prev = items[idx].prev;
alpar@103
   270
      
alpar@103
   271
      items[idx].next = firstFreeItem;
alpar@103
   272
      firstFreeItem = idx;
alpar@103
   273
    }
alpar@103
   274
alpar@103
   275
    void spliceItems(int ak, int bk) {
alpar@103
   276
      items[items[ak].prev].next = bk;
alpar@103
   277
      items[items[bk].prev].next = ak;
alpar@103
   278
      int tmp = items[ak].prev;
alpar@103
   279
      items[ak].prev = items[bk].prev;
alpar@103
   280
      items[bk].prev = tmp;
alpar@103
   281
        
alpar@103
   282
    }
alpar@103
   283
alpar@103
   284
    void laceClass(int cls) {
alpar@103
   285
      if (firstClass != -1) {
alpar@103
   286
        classes[firstClass].prev = cls;
alpar@103
   287
      }
alpar@103
   288
      classes[cls].next = firstClass;
alpar@103
   289
      classes[cls].prev = -1;
alpar@103
   290
      firstClass = cls;
alpar@103
   291
    } 
alpar@103
   292
alpar@103
   293
    void unlaceClass(int cls) {
alpar@103
   294
      if (classes[cls].prev != -1) {
alpar@103
   295
        classes[classes[cls].prev].next = classes[cls].next;
alpar@103
   296
      } else {
alpar@103
   297
        firstClass = classes[cls].next;
alpar@103
   298
      }
alpar@103
   299
      if (classes[cls].next != -1) {
alpar@103
   300
        classes[classes[cls].next].prev = classes[cls].prev;
alpar@103
   301
      }
alpar@103
   302
      
alpar@103
   303
      classes[cls].next = firstFreeClass;
alpar@103
   304
      firstFreeClass = cls;
alpar@103
   305
    } 
alpar@103
   306
alpar@103
   307
  public:
alpar@103
   308
alpar@103
   309
    UnionFindEnum(ItemIntMap& _index) 
alpar@103
   310
      : index(_index), items(), firstFreeItem(-1), 
alpar@103
   311
	firstClass(-1), firstFreeClass(-1) {}
alpar@103
   312
    
alpar@103
   313
    /// \brief Inserts the given element into a new component.
alpar@103
   314
    ///
alpar@103
   315
    /// This method creates a new component consisting only of the
alpar@103
   316
    /// given element.
alpar@103
   317
    ///
alpar@103
   318
    int insert(const Item& item) {
alpar@103
   319
      int idx = newItem();
alpar@103
   320
alpar@103
   321
      index.set(item, idx);
alpar@103
   322
alpar@103
   323
      singletonItem(idx);
alpar@103
   324
      items[idx].item = item;
alpar@103
   325
alpar@103
   326
      int cdx = newClass();
alpar@103
   327
alpar@103
   328
      items[idx].parent = ~cdx;
alpar@103
   329
alpar@103
   330
      laceClass(cdx);
alpar@103
   331
      classes[cdx].size = 1;
alpar@103
   332
      classes[cdx].firstItem = idx;
alpar@103
   333
alpar@103
   334
      firstClass = cdx;
alpar@103
   335
      
alpar@103
   336
      return cdx;
alpar@103
   337
    }
alpar@103
   338
alpar@103
   339
    /// \brief Inserts the given element into the component of the others.
alpar@103
   340
    ///
alpar@103
   341
    /// This methods inserts the element \e a into the component of the
alpar@103
   342
    /// element \e comp. 
alpar@103
   343
    void insert(const Item& item, int cls) {
alpar@103
   344
      int rdx = classes[cls].firstItem;
alpar@103
   345
      int idx = newItem();
alpar@103
   346
alpar@103
   347
      index.set(item, idx);
alpar@103
   348
alpar@103
   349
      laceItem(idx, rdx);
alpar@103
   350
alpar@103
   351
      items[idx].item = item;
alpar@103
   352
      items[idx].parent = rdx;
alpar@103
   353
alpar@103
   354
      ++classes[~(items[rdx].parent)].size;
alpar@103
   355
    }
alpar@103
   356
alpar@103
   357
    /// \brief Clears the union-find data structure
alpar@103
   358
    ///
alpar@103
   359
    /// Erase each item from the data structure.
alpar@103
   360
    void clear() {
alpar@103
   361
      items.clear();
alpar@103
   362
      firstClass = -1;
alpar@103
   363
      firstFreeItem = -1;
alpar@103
   364
    }
alpar@103
   365
alpar@103
   366
    /// \brief Finds the component of the given element.
alpar@103
   367
    ///
alpar@103
   368
    /// The method returns the component id of the given element.
alpar@103
   369
    int find(const Item &item) const {
alpar@103
   370
      return ~(items[repIndex(index[item])].parent);
alpar@103
   371
    }
alpar@103
   372
alpar@103
   373
    /// \brief Joining the component of element \e a and element \e b.
alpar@103
   374
    ///
alpar@103
   375
    /// This is the \e union operation of the Union-Find structure. 
alpar@103
   376
    /// Joins the component of element \e a and component of
alpar@103
   377
    /// element \e b. If \e a and \e b are in the same component then
alpar@103
   378
    /// returns -1 else returns the remaining class.
alpar@103
   379
    int join(const Item& a, const Item& b) {
alpar@103
   380
alpar@103
   381
      int ak = repIndex(index[a]);
alpar@103
   382
      int bk = repIndex(index[b]);
alpar@103
   383
alpar@103
   384
      if (ak == bk) {
alpar@103
   385
	return -1;
alpar@103
   386
      }
alpar@103
   387
alpar@103
   388
      int acx = ~(items[ak].parent);
alpar@103
   389
      int bcx = ~(items[bk].parent);
alpar@103
   390
alpar@103
   391
      int rcx;
alpar@103
   392
alpar@103
   393
      if (classes[acx].size > classes[bcx].size) {
alpar@103
   394
	classes[acx].size += classes[bcx].size;
alpar@103
   395
	items[bk].parent = ak;
alpar@103
   396
        unlaceClass(bcx);
alpar@103
   397
	rcx = acx;
alpar@103
   398
      } else {
alpar@103
   399
	classes[bcx].size += classes[acx].size;
alpar@103
   400
	items[ak].parent = bk;
alpar@103
   401
        unlaceClass(acx);
alpar@103
   402
	rcx = bcx;
alpar@103
   403
      }
alpar@103
   404
      spliceItems(ak, bk);
alpar@103
   405
alpar@103
   406
      return rcx;
alpar@103
   407
    }
alpar@103
   408
alpar@103
   409
    /// \brief Returns the size of the class.
alpar@103
   410
    ///
alpar@103
   411
    /// Returns the size of the class.
alpar@103
   412
    int size(int cls) const {
alpar@103
   413
      return classes[cls].size;
alpar@103
   414
    }
alpar@103
   415
alpar@103
   416
    /// \brief Splits up the component. 
alpar@103
   417
    ///
alpar@103
   418
    /// Splitting the component into singleton components (component
alpar@103
   419
    /// of size one).
alpar@103
   420
    void split(int cls) {
alpar@103
   421
      int fdx = classes[cls].firstItem;
alpar@103
   422
      int idx = items[fdx].next;
alpar@103
   423
      while (idx != fdx) {
alpar@103
   424
        int next = items[idx].next;
alpar@103
   425
alpar@103
   426
	singletonItem(idx);
alpar@103
   427
alpar@103
   428
	int cdx = newClass();        
alpar@103
   429
        items[idx].parent = ~cdx;
alpar@103
   430
alpar@103
   431
	laceClass(cdx);
alpar@103
   432
	classes[cdx].size = 1;
alpar@103
   433
	classes[cdx].firstItem = idx;
alpar@103
   434
        
alpar@103
   435
        idx = next;
alpar@103
   436
      }
alpar@103
   437
alpar@103
   438
      items[idx].prev = idx;
alpar@103
   439
      items[idx].next = idx;
alpar@103
   440
alpar@103
   441
      classes[~(items[idx].parent)].size = 1;
alpar@103
   442
      
alpar@103
   443
    }
alpar@103
   444
alpar@103
   445
    /// \brief Removes the given element from the structure.
alpar@103
   446
    ///
alpar@103
   447
    /// Removes the element from its component and if the component becomes
alpar@103
   448
    /// empty then removes that component from the component list.
alpar@103
   449
    ///
alpar@103
   450
    /// \warning It is an error to remove an element which is not in
alpar@103
   451
    /// the structure.
alpar@103
   452
    /// \warning This running time of this operation is proportional to the
alpar@103
   453
    /// number of the items in this class.
alpar@103
   454
    void erase(const Item& item) {
alpar@103
   455
      int idx = index[item];
alpar@103
   456
      int fdx = items[idx].next;
alpar@103
   457
alpar@103
   458
      int cdx = classIndex(idx);
alpar@103
   459
      if (idx == fdx) {
alpar@103
   460
	unlaceClass(cdx);
alpar@103
   461
	items[idx].next = firstFreeItem;
alpar@103
   462
	firstFreeItem = idx;
alpar@103
   463
	return;
alpar@103
   464
      } else {
alpar@103
   465
	classes[cdx].firstItem = fdx;
alpar@103
   466
	--classes[cdx].size;
alpar@103
   467
	items[fdx].parent = ~cdx;
alpar@103
   468
alpar@103
   469
	unlaceItem(idx);
alpar@103
   470
	idx = items[fdx].next;
alpar@103
   471
	while (idx != fdx) {
alpar@103
   472
	  items[idx].parent = fdx;
alpar@103
   473
	  idx = items[idx].next;
alpar@103
   474
	}
alpar@103
   475
          
alpar@103
   476
      }
alpar@103
   477
alpar@103
   478
    }
alpar@103
   479
alpar@103
   480
    /// \brief Gives back a representant item of the component.
alpar@103
   481
    ///
alpar@103
   482
    /// Gives back a representant item of the component.
alpar@103
   483
    Item item(int cls) const {
alpar@103
   484
      return items[classes[cls].firstItem].item;
alpar@103
   485
    }
alpar@103
   486
alpar@103
   487
    /// \brief Removes the component of the given element from the structure.
alpar@103
   488
    ///
alpar@103
   489
    /// Removes the component of the given element from the structure.
alpar@103
   490
    ///
alpar@103
   491
    /// \warning It is an error to give an element which is not in the
alpar@103
   492
    /// structure.
alpar@103
   493
    void eraseClass(int cls) {
alpar@103
   494
      int fdx = classes[cls].firstItem;
alpar@103
   495
      unlaceClass(cls);
alpar@103
   496
      items[items[fdx].prev].next = firstFreeItem;
alpar@103
   497
      firstFreeItem = fdx;
alpar@103
   498
    }
alpar@103
   499
alpar@103
   500
    /// \brief Lemon style iterator for the representant items.
alpar@103
   501
    ///
alpar@103
   502
    /// ClassIt is a lemon style iterator for the components. It iterates
alpar@103
   503
    /// on the ids of the classes.
alpar@103
   504
    class ClassIt {
alpar@103
   505
    public:
alpar@103
   506
      /// \brief Constructor of the iterator
alpar@103
   507
      ///
alpar@103
   508
      /// Constructor of the iterator
alpar@103
   509
      ClassIt(const UnionFindEnum& ufe) : unionFind(&ufe) {
alpar@103
   510
        cdx = unionFind->firstClass;
alpar@103
   511
      }
alpar@103
   512
alpar@103
   513
      /// \brief Constructor to get invalid iterator
alpar@103
   514
      ///
alpar@103
   515
      /// Constructor to get invalid iterator
alpar@103
   516
      ClassIt(Invalid) : unionFind(0), cdx(-1) {}
alpar@103
   517
      
alpar@103
   518
      /// \brief Increment operator
alpar@103
   519
      ///
alpar@103
   520
      /// It steps to the next representant item.
alpar@103
   521
      ClassIt& operator++() {
alpar@103
   522
        cdx = unionFind->classes[cdx].next;
alpar@103
   523
        return *this;
alpar@103
   524
      }
alpar@103
   525
      
alpar@103
   526
      /// \brief Conversion operator
alpar@103
   527
      ///
alpar@103
   528
      /// It converts the iterator to the current representant item.
alpar@103
   529
      operator int() const {
alpar@103
   530
        return cdx;
alpar@103
   531
      }
alpar@103
   532
alpar@103
   533
      /// \brief Equality operator
alpar@103
   534
      ///
alpar@103
   535
      /// Equality operator
alpar@103
   536
      bool operator==(const ClassIt& i) { 
alpar@103
   537
        return i.cdx == cdx;
alpar@103
   538
      }
alpar@103
   539
alpar@103
   540
      /// \brief Inequality operator
alpar@103
   541
      ///
alpar@103
   542
      /// Inequality operator
alpar@103
   543
      bool operator!=(const ClassIt& i) { 
alpar@103
   544
        return i.cdx != cdx;
alpar@103
   545
      }
alpar@103
   546
      
alpar@103
   547
    private:
alpar@103
   548
      const UnionFindEnum* unionFind;
alpar@103
   549
      int cdx;
alpar@103
   550
    };
alpar@103
   551
alpar@103
   552
    /// \brief Lemon style iterator for the items of a component.
alpar@103
   553
    ///
alpar@103
   554
    /// ClassIt is a lemon style iterator for the components. It iterates
alpar@103
   555
    /// on the items of a class. By example if you want to iterate on
alpar@103
   556
    /// each items of each classes then you may write the next code.
alpar@103
   557
    ///\code
alpar@103
   558
    /// for (ClassIt cit(ufe); cit != INVALID; ++cit) {
alpar@103
   559
    ///   std::cout << "Class: ";
alpar@103
   560
    ///   for (ItemIt iit(ufe, cit); iit != INVALID; ++iit) {
alpar@103
   561
    ///     std::cout << toString(iit) << ' ' << std::endl;
alpar@103
   562
    ///   }
alpar@103
   563
    ///   std::cout << std::endl;
alpar@103
   564
    /// }
alpar@103
   565
    ///\endcode
alpar@103
   566
    class ItemIt {
alpar@103
   567
    public:
alpar@103
   568
      /// \brief Constructor of the iterator
alpar@103
   569
      ///
alpar@103
   570
      /// Constructor of the iterator. The iterator iterates
alpar@103
   571
      /// on the class of the \c item.
alpar@103
   572
      ItemIt(const UnionFindEnum& ufe, int cls) : unionFind(&ufe) {
alpar@103
   573
        fdx = idx = unionFind->classes[cls].firstItem;
alpar@103
   574
      }
alpar@103
   575
alpar@103
   576
      /// \brief Constructor to get invalid iterator
alpar@103
   577
      ///
alpar@103
   578
      /// Constructor to get invalid iterator
alpar@103
   579
      ItemIt(Invalid) : unionFind(0), idx(-1) {}
alpar@103
   580
      
alpar@103
   581
      /// \brief Increment operator
alpar@103
   582
      ///
alpar@103
   583
      /// It steps to the next item in the class.
alpar@103
   584
      ItemIt& operator++() {
alpar@103
   585
        idx = unionFind->items[idx].next;
alpar@103
   586
        if (idx == fdx) idx = -1;
alpar@103
   587
        return *this;
alpar@103
   588
      }
alpar@103
   589
      
alpar@103
   590
      /// \brief Conversion operator
alpar@103
   591
      ///
alpar@103
   592
      /// It converts the iterator to the current item.
alpar@103
   593
      operator const Item&() const {
alpar@103
   594
        return unionFind->items[idx].item;
alpar@103
   595
      }
alpar@103
   596
alpar@103
   597
      /// \brief Equality operator
alpar@103
   598
      ///
alpar@103
   599
      /// Equality operator
alpar@103
   600
      bool operator==(const ItemIt& i) { 
alpar@103
   601
        return i.idx == idx;
alpar@103
   602
      }
alpar@103
   603
alpar@103
   604
      /// \brief Inequality operator
alpar@103
   605
      ///
alpar@103
   606
      /// Inequality operator
alpar@103
   607
      bool operator!=(const ItemIt& i) { 
alpar@103
   608
        return i.idx != idx;
alpar@103
   609
      }
alpar@103
   610
      
alpar@103
   611
    private:
alpar@103
   612
      const UnionFindEnum* unionFind;
alpar@103
   613
      int idx, fdx;
alpar@103
   614
    };
alpar@103
   615
alpar@103
   616
  };
alpar@103
   617
alpar@103
   618
  /// \ingroup auxdat
alpar@103
   619
  ///
alpar@103
   620
  /// \brief A \e Extend-Find data structure implementation which
alpar@103
   621
  /// is able to enumerate the components.
alpar@103
   622
  ///
alpar@103
   623
  /// The class implements an \e Extend-Find data structure which is
alpar@103
   624
  /// able to enumerate the components and the items in a
alpar@103
   625
  /// component. The data structure is a simplification of the
alpar@103
   626
  /// Union-Find structure, and it does not allow to merge two components.
alpar@103
   627
  ///
alpar@103
   628
  /// \pre You need to add all the elements by the \ref insert()
alpar@103
   629
  /// method.
alpar@103
   630
  template <typename _ItemIntMap>
alpar@103
   631
  class ExtendFindEnum {
alpar@103
   632
  public:
alpar@103
   633
    
alpar@103
   634
    typedef _ItemIntMap ItemIntMap;
alpar@103
   635
    typedef typename ItemIntMap::Key Item;
alpar@103
   636
alpar@103
   637
  private:
alpar@103
   638
    
alpar@103
   639
    ItemIntMap& index;
alpar@103
   640
alpar@103
   641
    struct ItemT {
alpar@103
   642
      int cls;
alpar@103
   643
      Item item;
alpar@103
   644
      int next, prev;
alpar@103
   645
    };
alpar@103
   646
alpar@103
   647
    std::vector<ItemT> items;
alpar@103
   648
    int firstFreeItem;
alpar@103
   649
alpar@103
   650
    struct ClassT {
alpar@103
   651
      int firstItem;
alpar@103
   652
      int next, prev;
alpar@103
   653
    };
alpar@103
   654
alpar@103
   655
    std::vector<ClassT> classes;
alpar@103
   656
alpar@103
   657
    int firstClass, firstFreeClass;
alpar@103
   658
alpar@103
   659
    int newClass() {
alpar@103
   660
      if (firstFreeClass != -1) {
alpar@103
   661
	int cdx = firstFreeClass;
alpar@103
   662
	firstFreeClass = classes[cdx].next;
alpar@103
   663
	return cdx;
alpar@103
   664
      } else {
alpar@103
   665
	classes.push_back(ClassT());
alpar@103
   666
	return classes.size() - 1;
alpar@103
   667
      }
alpar@103
   668
    }
alpar@103
   669
alpar@103
   670
    int newItem() {
alpar@103
   671
      if (firstFreeItem != -1) {
alpar@103
   672
	int idx = firstFreeItem;
alpar@103
   673
	firstFreeItem = items[idx].next;
alpar@103
   674
	return idx;
alpar@103
   675
      } else {
alpar@103
   676
	items.push_back(ItemT());
alpar@103
   677
	return items.size() - 1;
alpar@103
   678
      }
alpar@103
   679
    }
alpar@103
   680
alpar@103
   681
  public:
alpar@103
   682
alpar@103
   683
    /// \brief Constructor
alpar@103
   684
    ExtendFindEnum(ItemIntMap& _index) 
alpar@103
   685
      : index(_index), items(), firstFreeItem(-1), 
alpar@103
   686
	classes(), firstClass(-1), firstFreeClass(-1) {}
alpar@103
   687
    
alpar@103
   688
    /// \brief Inserts the given element into a new component.
alpar@103
   689
    ///
alpar@103
   690
    /// This method creates a new component consisting only of the
alpar@103
   691
    /// given element.
alpar@103
   692
    int insert(const Item& item) {
alpar@103
   693
      int cdx = newClass();
alpar@103
   694
      classes[cdx].prev = -1;
alpar@103
   695
      classes[cdx].next = firstClass;
alpar@103
   696
      if (firstClass != -1) {
alpar@103
   697
	classes[firstClass].prev = cdx;
alpar@103
   698
      }
alpar@103
   699
      firstClass = cdx;
alpar@103
   700
      
alpar@103
   701
      int idx = newItem();
alpar@103
   702
      items[idx].item = item;
alpar@103
   703
      items[idx].cls = cdx;
alpar@103
   704
      items[idx].prev = idx;
alpar@103
   705
      items[idx].next = idx;
alpar@103
   706
alpar@103
   707
      classes[cdx].firstItem = idx;
alpar@103
   708
alpar@103
   709
      index.set(item, idx);
alpar@103
   710
      
alpar@103
   711
      return cdx;
alpar@103
   712
    }
alpar@103
   713
alpar@103
   714
    /// \brief Inserts the given element into the given component.
alpar@103
   715
    ///
alpar@103
   716
    /// This methods inserts the element \e item a into the \e cls class.
alpar@103
   717
    void insert(const Item& item, int cls) {
alpar@103
   718
      int idx = newItem();
alpar@103
   719
      int rdx = classes[cls].firstItem;
alpar@103
   720
      items[idx].item = item;
alpar@103
   721
      items[idx].cls = cls;
alpar@103
   722
alpar@103
   723
      items[idx].prev = rdx;
alpar@103
   724
      items[idx].next = items[rdx].next;
alpar@103
   725
      items[items[rdx].next].prev = idx;
alpar@103
   726
      items[rdx].next = idx;
alpar@103
   727
alpar@103
   728
      index.set(item, idx);
alpar@103
   729
    }
alpar@103
   730
alpar@103
   731
    /// \brief Clears the union-find data structure
alpar@103
   732
    ///
alpar@103
   733
    /// Erase each item from the data structure.
alpar@103
   734
    void clear() {
alpar@103
   735
      items.clear();
alpar@103
   736
      classes.clear;
alpar@103
   737
      firstClass = firstFreeClass = firstFreeItem = -1;
alpar@103
   738
    }
alpar@103
   739
alpar@103
   740
    /// \brief Gives back the class of the \e item.
alpar@103
   741
    ///
alpar@103
   742
    /// Gives back the class of the \e item.
alpar@103
   743
    int find(const Item &item) const {
alpar@103
   744
      return items[index[item]].cls;
alpar@103
   745
    }
alpar@103
   746
alpar@103
   747
    /// \brief Gives back a representant item of the component.
alpar@103
   748
    ///
alpar@103
   749
    /// Gives back a representant item of the component.
alpar@103
   750
    Item item(int cls) const {
alpar@103
   751
      return items[classes[cls].firstItem].item;
alpar@103
   752
    }
alpar@103
   753
    
alpar@103
   754
    /// \brief Removes the given element from the structure.
alpar@103
   755
    ///
alpar@103
   756
    /// Removes the element from its component and if the component becomes
alpar@103
   757
    /// empty then removes that component from the component list.
alpar@103
   758
    ///
alpar@103
   759
    /// \warning It is an error to remove an element which is not in
alpar@103
   760
    /// the structure.
alpar@103
   761
    void erase(const Item &item) {
alpar@103
   762
      int idx = index[item];
alpar@103
   763
      int cdx = items[idx].cls;
alpar@103
   764
      
alpar@103
   765
      if (idx == items[idx].next) {
alpar@103
   766
	if (classes[cdx].prev != -1) {
alpar@103
   767
	  classes[classes[cdx].prev].next = classes[cdx].next; 
alpar@103
   768
	} else {
alpar@103
   769
	  firstClass = classes[cdx].next;
alpar@103
   770
	}
alpar@103
   771
	if (classes[cdx].next != -1) {
alpar@103
   772
	  classes[classes[cdx].next].prev = classes[cdx].prev; 
alpar@103
   773
	}
alpar@103
   774
	classes[cdx].next = firstFreeClass;
alpar@103
   775
	firstFreeClass = cdx;
alpar@103
   776
      } else {
alpar@103
   777
	classes[cdx].firstItem = items[idx].next;
alpar@103
   778
	items[items[idx].next].prev = items[idx].prev;
alpar@103
   779
	items[items[idx].prev].next = items[idx].next;
alpar@103
   780
      }
alpar@103
   781
      items[idx].next = firstFreeItem;
alpar@103
   782
      firstFreeItem = idx;
alpar@103
   783
	
alpar@103
   784
    }    
alpar@103
   785
alpar@103
   786
    
alpar@103
   787
    /// \brief Removes the component of the given element from the structure.
alpar@103
   788
    ///
alpar@103
   789
    /// Removes the component of the given element from the structure.
alpar@103
   790
    ///
alpar@103
   791
    /// \warning It is an error to give an element which is not in the
alpar@103
   792
    /// structure.
alpar@103
   793
    void eraseClass(int cdx) {
alpar@103
   794
      int idx = classes[cdx].firstItem;
alpar@103
   795
      items[items[idx].prev].next = firstFreeItem;
alpar@103
   796
      firstFreeItem = idx;
alpar@103
   797
alpar@103
   798
      if (classes[cdx].prev != -1) {
alpar@103
   799
	classes[classes[cdx].prev].next = classes[cdx].next; 
alpar@103
   800
      } else {
alpar@103
   801
	firstClass = classes[cdx].next;
alpar@103
   802
      }
alpar@103
   803
      if (classes[cdx].next != -1) {
alpar@103
   804
	classes[classes[cdx].next].prev = classes[cdx].prev; 
alpar@103
   805
      }
alpar@103
   806
      classes[cdx].next = firstFreeClass;
alpar@103
   807
      firstFreeClass = cdx;
alpar@103
   808
    }
alpar@103
   809
alpar@103
   810
    /// \brief Lemon style iterator for the classes.
alpar@103
   811
    ///
alpar@103
   812
    /// ClassIt is a lemon style iterator for the components. It iterates
alpar@103
   813
    /// on the ids of classes.
alpar@103
   814
    class ClassIt {
alpar@103
   815
    public:
alpar@103
   816
      /// \brief Constructor of the iterator
alpar@103
   817
      ///
alpar@103
   818
      /// Constructor of the iterator
alpar@103
   819
      ClassIt(const ExtendFindEnum& ufe) : extendFind(&ufe) {
alpar@103
   820
        cdx = extendFind->firstClass;
alpar@103
   821
      }
alpar@103
   822
alpar@103
   823
      /// \brief Constructor to get invalid iterator
alpar@103
   824
      ///
alpar@103
   825
      /// Constructor to get invalid iterator
alpar@103
   826
      ClassIt(Invalid) : extendFind(0), cdx(-1) {}
alpar@103
   827
      
alpar@103
   828
      /// \brief Increment operator
alpar@103
   829
      ///
alpar@103
   830
      /// It steps to the next representant item.
alpar@103
   831
      ClassIt& operator++() {
alpar@103
   832
        cdx = extendFind->classes[cdx].next;
alpar@103
   833
        return *this;
alpar@103
   834
      }
alpar@103
   835
      
alpar@103
   836
      /// \brief Conversion operator
alpar@103
   837
      ///
alpar@103
   838
      /// It converts the iterator to the current class id.
alpar@103
   839
      operator int() const {
alpar@103
   840
        return cdx;
alpar@103
   841
      }
alpar@103
   842
alpar@103
   843
      /// \brief Equality operator
alpar@103
   844
      ///
alpar@103
   845
      /// Equality operator
alpar@103
   846
      bool operator==(const ClassIt& i) { 
alpar@103
   847
        return i.cdx == cdx;
alpar@103
   848
      }
alpar@103
   849
alpar@103
   850
      /// \brief Inequality operator
alpar@103
   851
      ///
alpar@103
   852
      /// Inequality operator
alpar@103
   853
      bool operator!=(const ClassIt& i) { 
alpar@103
   854
        return i.cdx != cdx;
alpar@103
   855
      }
alpar@103
   856
      
alpar@103
   857
    private:
alpar@103
   858
      const ExtendFindEnum* extendFind;
alpar@103
   859
      int cdx;
alpar@103
   860
    };
alpar@103
   861
alpar@103
   862
    /// \brief Lemon style iterator for the items of a component.
alpar@103
   863
    ///
alpar@103
   864
    /// ClassIt is a lemon style iterator for the components. It iterates
alpar@103
   865
    /// on the items of a class. By example if you want to iterate on
alpar@103
   866
    /// each items of each classes then you may write the next code.
alpar@103
   867
    ///\code
alpar@103
   868
    /// for (ClassIt cit(ufe); cit != INVALID; ++cit) {
alpar@103
   869
    ///   std::cout << "Class: ";
alpar@103
   870
    ///   for (ItemIt iit(ufe, cit); iit != INVALID; ++iit) {
alpar@103
   871
    ///     std::cout << toString(iit) << ' ' << std::endl;
alpar@103
   872
    ///   }
alpar@103
   873
    ///   std::cout << std::endl;
alpar@103
   874
    /// }
alpar@103
   875
    ///\endcode
alpar@103
   876
    class ItemIt {
alpar@103
   877
    public:
alpar@103
   878
      /// \brief Constructor of the iterator
alpar@103
   879
      ///
alpar@103
   880
      /// Constructor of the iterator. The iterator iterates
alpar@103
   881
      /// on the class of the \c item.
alpar@103
   882
      ItemIt(const ExtendFindEnum& ufe, int cls) : extendFind(&ufe) {
alpar@103
   883
        fdx = idx = extendFind->classes[cls].firstItem;
alpar@103
   884
      }
alpar@103
   885
alpar@103
   886
      /// \brief Constructor to get invalid iterator
alpar@103
   887
      ///
alpar@103
   888
      /// Constructor to get invalid iterator
alpar@103
   889
      ItemIt(Invalid) : extendFind(0), idx(-1) {}
alpar@103
   890
      
alpar@103
   891
      /// \brief Increment operator
alpar@103
   892
      ///
alpar@103
   893
      /// It steps to the next item in the class.
alpar@103
   894
      ItemIt& operator++() {
alpar@103
   895
        idx = extendFind->items[idx].next;
alpar@103
   896
	if (fdx == idx) idx = -1;
alpar@103
   897
        return *this;
alpar@103
   898
      }
alpar@103
   899
      
alpar@103
   900
      /// \brief Conversion operator
alpar@103
   901
      ///
alpar@103
   902
      /// It converts the iterator to the current item.
alpar@103
   903
      operator const Item&() const {
alpar@103
   904
        return extendFind->items[idx].item;
alpar@103
   905
      }
alpar@103
   906
alpar@103
   907
      /// \brief Equality operator
alpar@103
   908
      ///
alpar@103
   909
      /// Equality operator
alpar@103
   910
      bool operator==(const ItemIt& i) { 
alpar@103
   911
        return i.idx == idx;
alpar@103
   912
      }
alpar@103
   913
alpar@103
   914
      /// \brief Inequality operator
alpar@103
   915
      ///
alpar@103
   916
      /// Inequality operator
alpar@103
   917
      bool operator!=(const ItemIt& i) { 
alpar@103
   918
        return i.idx != idx;
alpar@103
   919
      }
alpar@103
   920
      
alpar@103
   921
    private:
alpar@103
   922
      const ExtendFindEnum* extendFind;
alpar@103
   923
      int idx, fdx;
alpar@103
   924
    };
alpar@103
   925
alpar@103
   926
  };
alpar@103
   927
alpar@103
   928
  /// \ingroup auxdat
alpar@103
   929
  ///
alpar@103
   930
  /// \brief A \e Union-Find data structure implementation which
alpar@103
   931
  /// is able to store a priority for each item and retrieve the minimum of
alpar@103
   932
  /// each class.
alpar@103
   933
  ///
alpar@103
   934
  /// A \e Union-Find data structure implementation which is able to
alpar@103
   935
  /// store a priority for each item and retrieve the minimum of each
alpar@103
   936
  /// class. In addition, it supports the joining and splitting the
alpar@103
   937
  /// components. If you don't need this feature then you makes
alpar@103
   938
  /// better to use the \ref UnionFind class which is more efficient.
alpar@103
   939
  ///
alpar@103
   940
  /// The union-find data strcuture based on a (2, 16)-tree with a
alpar@103
   941
  /// tournament minimum selection on the internal nodes. The insert
alpar@103
   942
  /// operation takes O(1), the find, set, decrease and increase takes
alpar@103
   943
  /// O(log(n)), where n is the number of nodes in the current
alpar@103
   944
  /// component.  The complexity of join and split is O(log(n)*k),
alpar@103
   945
  /// where n is the sum of the number of the nodes and k is the
alpar@103
   946
  /// number of joined components or the number of the components
alpar@103
   947
  /// after the split.
alpar@103
   948
  ///
alpar@103
   949
  /// \pre You need to add all the elements by the \ref insert()
alpar@103
   950
  /// method.
alpar@103
   951
  ///
alpar@103
   952
  template <typename _Value, typename _ItemIntMap, 
alpar@103
   953
            typename _Comp = std::less<_Value> >
alpar@103
   954
  class HeapUnionFind {
alpar@103
   955
  public:
alpar@103
   956
    
alpar@103
   957
    typedef _Value Value;
alpar@103
   958
    typedef typename _ItemIntMap::Key Item;
alpar@103
   959
alpar@103
   960
    typedef _ItemIntMap ItemIntMap;
alpar@103
   961
alpar@103
   962
    typedef _Comp Comp;
alpar@103
   963
alpar@103
   964
  private:
alpar@103
   965
alpar@103
   966
    static const int cmax = 16;
alpar@103
   967
alpar@103
   968
    ItemIntMap& index;
alpar@103
   969
alpar@103
   970
    struct ClassNode {
alpar@103
   971
      int parent;
alpar@103
   972
      int depth;
alpar@103
   973
alpar@103
   974
      int left, right;
alpar@103
   975
      int next, prev;
alpar@103
   976
    };
alpar@103
   977
alpar@103
   978
    int first_class;
alpar@103
   979
    int first_free_class;
alpar@103
   980
    std::vector<ClassNode> classes;
alpar@103
   981
alpar@103
   982
    int newClass() {
alpar@103
   983
      if (first_free_class < 0) {
alpar@103
   984
        int id = classes.size();
alpar@103
   985
        classes.push_back(ClassNode());
alpar@103
   986
        return id;
alpar@103
   987
      } else {
alpar@103
   988
        int id = first_free_class;
alpar@103
   989
        first_free_class = classes[id].next;
alpar@103
   990
        return id;
alpar@103
   991
      }
alpar@103
   992
    }
alpar@103
   993
alpar@103
   994
    void deleteClass(int id) {
alpar@103
   995
      classes[id].next = first_free_class;
alpar@103
   996
      first_free_class = id;
alpar@103
   997
    }
alpar@103
   998
alpar@103
   999
    struct ItemNode {
alpar@103
  1000
      int parent;
alpar@103
  1001
      Item item;
alpar@103
  1002
      Value prio;
alpar@103
  1003
      int next, prev;
alpar@103
  1004
      int left, right;
alpar@103
  1005
      int size;
alpar@103
  1006
    };
alpar@103
  1007
alpar@103
  1008
    int first_free_node;
alpar@103
  1009
    std::vector<ItemNode> nodes;
alpar@103
  1010
alpar@103
  1011
    int newNode() {
alpar@103
  1012
      if (first_free_node < 0) {
alpar@103
  1013
        int id = nodes.size();
alpar@103
  1014
        nodes.push_back(ItemNode());
alpar@103
  1015
        return id;
alpar@103
  1016
      } else {
alpar@103
  1017
        int id = first_free_node;
alpar@103
  1018
        first_free_node = nodes[id].next;
alpar@103
  1019
        return id;
alpar@103
  1020
      }
alpar@103
  1021
    }
alpar@103
  1022
alpar@103
  1023
    void deleteNode(int id) {
alpar@103
  1024
      nodes[id].next = first_free_node;
alpar@103
  1025
      first_free_node = id;
alpar@103
  1026
    }
alpar@103
  1027
alpar@103
  1028
    Comp comp;
alpar@103
  1029
alpar@103
  1030
    int findClass(int id) const {
alpar@103
  1031
      int kd = id;
alpar@103
  1032
      while (kd >= 0) {
alpar@103
  1033
        kd = nodes[kd].parent;
alpar@103
  1034
      }
alpar@103
  1035
      return ~kd;
alpar@103
  1036
    }
alpar@103
  1037
alpar@103
  1038
    int leftNode(int id) const {
alpar@103
  1039
      int kd = ~(classes[id].parent);
alpar@103
  1040
      for (int i = 0; i < classes[id].depth; ++i) {
alpar@103
  1041
        kd = nodes[kd].left;
alpar@103
  1042
      }
alpar@103
  1043
      return kd;
alpar@103
  1044
    }
alpar@103
  1045
alpar@103
  1046
    int nextNode(int id) const {
alpar@103
  1047
      int depth = 0;
alpar@103
  1048
      while (id >= 0 && nodes[id].next == -1) {
alpar@103
  1049
        id = nodes[id].parent;
alpar@103
  1050
        ++depth;
alpar@103
  1051
      }
alpar@103
  1052
      if (id < 0) {
alpar@103
  1053
        return -1;
alpar@103
  1054
      }
alpar@103
  1055
      id = nodes[id].next;
alpar@103
  1056
      while (depth--) {
alpar@103
  1057
        id = nodes[id].left;      
alpar@103
  1058
      }
alpar@103
  1059
      return id;
alpar@103
  1060
    }
alpar@103
  1061
alpar@103
  1062
alpar@103
  1063
    void setPrio(int id) {
alpar@103
  1064
      int jd = nodes[id].left;
alpar@103
  1065
      nodes[id].prio = nodes[jd].prio;
alpar@103
  1066
      nodes[id].item = nodes[jd].item;
alpar@103
  1067
      jd = nodes[jd].next;
alpar@103
  1068
      while (jd != -1) {
alpar@103
  1069
        if (comp(nodes[jd].prio, nodes[id].prio)) {
alpar@103
  1070
          nodes[id].prio = nodes[jd].prio;
alpar@103
  1071
          nodes[id].item = nodes[jd].item;
alpar@103
  1072
        }
alpar@103
  1073
        jd = nodes[jd].next;
alpar@103
  1074
      }
alpar@103
  1075
    }
alpar@103
  1076
alpar@103
  1077
    void push(int id, int jd) {
alpar@103
  1078
      nodes[id].size = 1;
alpar@103
  1079
      nodes[id].left = nodes[id].right = jd;
alpar@103
  1080
      nodes[jd].next = nodes[jd].prev = -1;
alpar@103
  1081
      nodes[jd].parent = id;
alpar@103
  1082
    }
alpar@103
  1083
alpar@103
  1084
    void pushAfter(int id, int jd) {
alpar@103
  1085
      int kd = nodes[id].parent;
alpar@103
  1086
      if (nodes[id].next != -1) {
alpar@103
  1087
        nodes[nodes[id].next].prev = jd;
alpar@103
  1088
        if (kd >= 0) {
alpar@103
  1089
          nodes[kd].size += 1;
alpar@103
  1090
        }
alpar@103
  1091
      } else {
alpar@103
  1092
        if (kd >= 0) {
alpar@103
  1093
          nodes[kd].right = jd;
alpar@103
  1094
          nodes[kd].size += 1;
alpar@103
  1095
        }
alpar@103
  1096
      }
alpar@103
  1097
      nodes[jd].next = nodes[id].next;
alpar@103
  1098
      nodes[jd].prev = id;
alpar@103
  1099
      nodes[id].next = jd;
alpar@103
  1100
      nodes[jd].parent = kd;
alpar@103
  1101
    }
alpar@103
  1102
alpar@103
  1103
    void pushRight(int id, int jd) {
alpar@103
  1104
      nodes[id].size += 1;
alpar@103
  1105
      nodes[jd].prev = nodes[id].right;
alpar@103
  1106
      nodes[jd].next = -1;
alpar@103
  1107
      nodes[nodes[id].right].next = jd;
alpar@103
  1108
      nodes[id].right = jd;
alpar@103
  1109
      nodes[jd].parent = id;
alpar@103
  1110
    }
alpar@103
  1111
alpar@103
  1112
    void popRight(int id) {
alpar@103
  1113
      nodes[id].size -= 1;
alpar@103
  1114
      int jd = nodes[id].right;
alpar@103
  1115
      nodes[nodes[jd].prev].next = -1;
alpar@103
  1116
      nodes[id].right = nodes[jd].prev;
alpar@103
  1117
    }
alpar@103
  1118
alpar@103
  1119
    void splice(int id, int jd) {
alpar@103
  1120
      nodes[id].size += nodes[jd].size;
alpar@103
  1121
      nodes[nodes[id].right].next = nodes[jd].left;
alpar@103
  1122
      nodes[nodes[jd].left].prev = nodes[id].right;
alpar@103
  1123
      int kd = nodes[jd].left;
alpar@103
  1124
      while (kd != -1) {
alpar@103
  1125
        nodes[kd].parent = id;
alpar@103
  1126
        kd = nodes[kd].next;
alpar@103
  1127
      }
alpar@103
  1128
      nodes[id].right = nodes[jd].right;
alpar@103
  1129
    }
alpar@103
  1130
alpar@103
  1131
    void split(int id, int jd) {
alpar@103
  1132
      int kd = nodes[id].parent;
alpar@103
  1133
      nodes[kd].right = nodes[id].prev;
alpar@103
  1134
      nodes[nodes[id].prev].next = -1;
alpar@103
  1135
      
alpar@103
  1136
      nodes[jd].left = id;
alpar@103
  1137
      nodes[id].prev = -1;
alpar@103
  1138
      int num = 0;
alpar@103
  1139
      while (id != -1) {
alpar@103
  1140
        nodes[id].parent = jd;
alpar@103
  1141
        nodes[jd].right = id;
alpar@103
  1142
        id = nodes[id].next;
alpar@103
  1143
        ++num;
alpar@103
  1144
      }      
alpar@103
  1145
      nodes[kd].size -= num;
alpar@103
  1146
      nodes[jd].size = num;
alpar@103
  1147
    }
alpar@103
  1148
alpar@103
  1149
    void pushLeft(int id, int jd) {
alpar@103
  1150
      nodes[id].size += 1;
alpar@103
  1151
      nodes[jd].next = nodes[id].left;
alpar@103
  1152
      nodes[jd].prev = -1;
alpar@103
  1153
      nodes[nodes[id].left].prev = jd;
alpar@103
  1154
      nodes[id].left = jd;
alpar@103
  1155
      nodes[jd].parent = id;
alpar@103
  1156
    }
alpar@103
  1157
alpar@103
  1158
    void popLeft(int id) {
alpar@103
  1159
      nodes[id].size -= 1;
alpar@103
  1160
      int jd = nodes[id].left;
alpar@103
  1161
      nodes[nodes[jd].next].prev = -1;
alpar@103
  1162
      nodes[id].left = nodes[jd].next;
alpar@103
  1163
    }
alpar@103
  1164
alpar@103
  1165
    void repairLeft(int id) {
alpar@103
  1166
      int jd = ~(classes[id].parent);
alpar@103
  1167
      while (nodes[jd].left != -1) {
alpar@103
  1168
	int kd = nodes[jd].left;
alpar@103
  1169
	if (nodes[jd].size == 1) {
alpar@103
  1170
	  if (nodes[jd].parent < 0) {
alpar@103
  1171
	    classes[id].parent = ~kd;
alpar@103
  1172
	    classes[id].depth -= 1;
alpar@103
  1173
	    nodes[kd].parent = ~id;
alpar@103
  1174
	    deleteNode(jd);
alpar@103
  1175
	    jd = kd;
alpar@103
  1176
	  } else {
alpar@103
  1177
	    int pd = nodes[jd].parent;
alpar@103
  1178
	    if (nodes[nodes[jd].next].size < cmax) {
alpar@103
  1179
	      pushLeft(nodes[jd].next, nodes[jd].left);
alpar@103
  1180
	      if (less(nodes[jd].left, nodes[jd].next)) {
alpar@103
  1181
		nodes[nodes[jd].next].prio = nodes[nodes[jd].left].prio;
alpar@103
  1182
		nodes[nodes[jd].next].item = nodes[nodes[jd].left].item;
alpar@103
  1183
	      } 
alpar@103
  1184
	      popLeft(pd);
alpar@103
  1185
	      deleteNode(jd);
alpar@103
  1186
	      jd = pd;
alpar@103
  1187
	    } else {
alpar@103
  1188
	      int ld = nodes[nodes[jd].next].left;
alpar@103
  1189
	      popLeft(nodes[jd].next);
alpar@103
  1190
	      pushRight(jd, ld);
alpar@103
  1191
	      if (less(ld, nodes[jd].left)) {
alpar@103
  1192
		nodes[jd].item = nodes[ld].item;
alpar@103
  1193
		nodes[jd].prio = nodes[jd].prio;
alpar@103
  1194
	      }
alpar@103
  1195
	      if (nodes[nodes[jd].next].item == nodes[ld].item) {
alpar@103
  1196
		setPrio(nodes[jd].next);
alpar@103
  1197
	      }
alpar@103
  1198
	      jd = nodes[jd].left;
alpar@103
  1199
	    }
alpar@103
  1200
	  }
alpar@103
  1201
	} else {
alpar@103
  1202
	  jd = nodes[jd].left;
alpar@103
  1203
	}
alpar@103
  1204
      }
alpar@103
  1205
    }    
alpar@103
  1206
alpar@103
  1207
    void repairRight(int id) {
alpar@103
  1208
      int jd = ~(classes[id].parent);
alpar@103
  1209
      while (nodes[jd].right != -1) {
alpar@103
  1210
	int kd = nodes[jd].right;
alpar@103
  1211
	if (nodes[jd].size == 1) {
alpar@103
  1212
	  if (nodes[jd].parent < 0) {
alpar@103
  1213
	    classes[id].parent = ~kd;
alpar@103
  1214
	    classes[id].depth -= 1;
alpar@103
  1215
	    nodes[kd].parent = ~id;
alpar@103
  1216
	    deleteNode(jd);
alpar@103
  1217
	    jd = kd;
alpar@103
  1218
	  } else {
alpar@103
  1219
	    int pd = nodes[jd].parent;
alpar@103
  1220
	    if (nodes[nodes[jd].prev].size < cmax) {
alpar@103
  1221
	      pushRight(nodes[jd].prev, nodes[jd].right);
alpar@103
  1222
	      if (less(nodes[jd].right, nodes[jd].prev)) {
alpar@103
  1223
		nodes[nodes[jd].prev].prio = nodes[nodes[jd].right].prio;
alpar@103
  1224
		nodes[nodes[jd].prev].item = nodes[nodes[jd].right].item;
alpar@103
  1225
	      } 
alpar@103
  1226
	      popRight(pd);
alpar@103
  1227
	      deleteNode(jd);
alpar@103
  1228
	      jd = pd;
alpar@103
  1229
	    } else {
alpar@103
  1230
	      int ld = nodes[nodes[jd].prev].right;
alpar@103
  1231
	      popRight(nodes[jd].prev);
alpar@103
  1232
	      pushLeft(jd, ld);
alpar@103
  1233
	      if (less(ld, nodes[jd].right)) {
alpar@103
  1234
		nodes[jd].item = nodes[ld].item;
alpar@103
  1235
		nodes[jd].prio = nodes[jd].prio;
alpar@103
  1236
	      }
alpar@103
  1237
	      if (nodes[nodes[jd].prev].item == nodes[ld].item) {
alpar@103
  1238
		setPrio(nodes[jd].prev);
alpar@103
  1239
	      }
alpar@103
  1240
	      jd = nodes[jd].right;
alpar@103
  1241
	    }
alpar@103
  1242
	  }
alpar@103
  1243
	} else {
alpar@103
  1244
	  jd = nodes[jd].right;
alpar@103
  1245
	}
alpar@103
  1246
      }
alpar@103
  1247
    }
alpar@103
  1248
alpar@103
  1249
alpar@103
  1250
    bool less(int id, int jd) const {
alpar@103
  1251
      return comp(nodes[id].prio, nodes[jd].prio);
alpar@103
  1252
    }
alpar@103
  1253
alpar@103
  1254
    bool equal(int id, int jd) const {
alpar@103
  1255
      return !less(id, jd) && !less(jd, id);
alpar@103
  1256
    }
alpar@103
  1257
alpar@103
  1258
alpar@103
  1259
  public:
alpar@103
  1260
alpar@103
  1261
    /// \brief Returns true when the given class is alive.
alpar@103
  1262
    ///
alpar@103
  1263
    /// Returns true when the given class is alive, ie. the class is
alpar@103
  1264
    /// not nested into other class.
alpar@103
  1265
    bool alive(int cls) const {
alpar@103
  1266
      return classes[cls].parent < 0;
alpar@103
  1267
    }
alpar@103
  1268
alpar@103
  1269
    /// \brief Returns true when the given class is trivial.
alpar@103
  1270
    ///
alpar@103
  1271
    /// Returns true when the given class is trivial, ie. the class
alpar@103
  1272
    /// contains just one item directly.
alpar@103
  1273
    bool trivial(int cls) const {
alpar@103
  1274
      return classes[cls].left == -1;
alpar@103
  1275
    }
alpar@103
  1276
alpar@103
  1277
    /// \brief Constructs the union-find.
alpar@103
  1278
    ///
alpar@103
  1279
    /// Constructs the union-find.  
alpar@103
  1280
    /// \brief _index The index map of the union-find. The data
alpar@103
  1281
    /// structure uses internally for store references.
alpar@103
  1282
    HeapUnionFind(ItemIntMap& _index) 
alpar@103
  1283
      : index(_index), first_class(-1), 
alpar@103
  1284
	first_free_class(-1), first_free_node(-1) {}
alpar@103
  1285
alpar@103
  1286
    /// \brief Insert a new node into a new component.
alpar@103
  1287
    ///
alpar@103
  1288
    /// Insert a new node into a new component.
alpar@103
  1289
    /// \param item The item of the new node.
alpar@103
  1290
    /// \param prio The priority of the new node.
alpar@103
  1291
    /// \return The class id of the one-item-heap.
alpar@103
  1292
    int insert(const Item& item, const Value& prio) {
alpar@103
  1293
      int id = newNode();
alpar@103
  1294
      nodes[id].item = item;
alpar@103
  1295
      nodes[id].prio = prio;
alpar@103
  1296
      nodes[id].size = 0;
alpar@103
  1297
alpar@103
  1298
      nodes[id].prev = -1;
alpar@103
  1299
      nodes[id].next = -1;
alpar@103
  1300
alpar@103
  1301
      nodes[id].left = -1;
alpar@103
  1302
      nodes[id].right = -1;
alpar@103
  1303
alpar@103
  1304
      nodes[id].item = item;
alpar@103
  1305
      index[item] = id;
alpar@103
  1306
      
alpar@103
  1307
      int class_id = newClass();
alpar@103
  1308
      classes[class_id].parent = ~id;
alpar@103
  1309
      classes[class_id].depth = 0;
alpar@103
  1310
alpar@103
  1311
      classes[class_id].left = -1;
alpar@103
  1312
      classes[class_id].right = -1;
alpar@103
  1313
      
alpar@103
  1314
      if (first_class != -1) {
alpar@103
  1315
        classes[first_class].prev = class_id;
alpar@103
  1316
      }
alpar@103
  1317
      classes[class_id].next = first_class;
alpar@103
  1318
      classes[class_id].prev = -1;
alpar@103
  1319
      first_class = class_id;
alpar@103
  1320
alpar@103
  1321
      nodes[id].parent = ~class_id;
alpar@103
  1322
      
alpar@103
  1323
      return class_id;
alpar@103
  1324
    }
alpar@103
  1325
alpar@103
  1326
    /// \brief The class of the item.
alpar@103
  1327
    ///
alpar@103
  1328
    /// \return The alive class id of the item, which is not nested into
alpar@103
  1329
    /// other classes.
alpar@103
  1330
    ///
alpar@103
  1331
    /// The time complexity is O(log(n)).
alpar@103
  1332
    int find(const Item& item) const {
alpar@103
  1333
      return findClass(index[item]);
alpar@103
  1334
    }
alpar@103
  1335
    
alpar@103
  1336
    /// \brief Joins the classes.
alpar@103
  1337
    ///
alpar@103
  1338
    /// The current function joins the given classes. The parameter is
alpar@103
  1339
    /// an STL range which should be contains valid class ids. The
alpar@103
  1340
    /// time complexity is O(log(n)*k) where n is the overall number
alpar@103
  1341
    /// of the joined nodes and k is the number of classes.
alpar@103
  1342
    /// \return The class of the joined classes.
alpar@103
  1343
    /// \pre The range should contain at least two class ids.
alpar@103
  1344
    template <typename Iterator>
alpar@103
  1345
    int join(Iterator begin, Iterator end) {
alpar@103
  1346
      std::vector<int> cs;
alpar@103
  1347
      for (Iterator it = begin; it != end; ++it) {
alpar@103
  1348
        cs.push_back(*it);
alpar@103
  1349
      }
alpar@103
  1350
alpar@103
  1351
      int class_id = newClass();
alpar@103
  1352
      { // creation union-find
alpar@103
  1353
alpar@103
  1354
        if (first_class != -1) {
alpar@103
  1355
          classes[first_class].prev = class_id;
alpar@103
  1356
        }
alpar@103
  1357
        classes[class_id].next = first_class;
alpar@103
  1358
        classes[class_id].prev = -1;
alpar@103
  1359
        first_class = class_id;
alpar@103
  1360
alpar@103
  1361
        classes[class_id].depth = classes[cs[0]].depth;
alpar@103
  1362
        classes[class_id].parent = classes[cs[0]].parent;
alpar@103
  1363
        nodes[~(classes[class_id].parent)].parent = ~class_id;
alpar@103
  1364
alpar@103
  1365
        int l = cs[0];
alpar@103
  1366
alpar@103
  1367
        classes[class_id].left = l;
alpar@103
  1368
        classes[class_id].right = l;
alpar@103
  1369
alpar@103
  1370
        if (classes[l].next != -1) {
alpar@103
  1371
          classes[classes[l].next].prev = classes[l].prev;
alpar@103
  1372
        }
alpar@103
  1373
        classes[classes[l].prev].next = classes[l].next;
alpar@103
  1374
        
alpar@103
  1375
        classes[l].prev = -1;
alpar@103
  1376
        classes[l].next = -1;
alpar@103
  1377
alpar@103
  1378
        classes[l].depth = leftNode(l);
alpar@103
  1379
        classes[l].parent = class_id;
alpar@103
  1380
        
alpar@103
  1381
      }
alpar@103
  1382
alpar@103
  1383
      { // merging of heap
alpar@103
  1384
        int l = class_id;
alpar@103
  1385
        for (int ci = 1; ci < int(cs.size()); ++ci) {
alpar@103
  1386
          int r = cs[ci];
alpar@103
  1387
          int rln = leftNode(r);
alpar@103
  1388
          if (classes[l].depth > classes[r].depth) {
alpar@103
  1389
            int id = ~(classes[l].parent);
alpar@103
  1390
            for (int i = classes[r].depth + 1; i < classes[l].depth; ++i) {
alpar@103
  1391
              id = nodes[id].right;
alpar@103
  1392
            }
alpar@103
  1393
            while (id >= 0 && nodes[id].size == cmax) {
alpar@103
  1394
              int new_id = newNode();
alpar@103
  1395
              int right_id = nodes[id].right;
alpar@103
  1396
alpar@103
  1397
              popRight(id);
alpar@103
  1398
              if (nodes[id].item == nodes[right_id].item) {
alpar@103
  1399
                setPrio(id);
alpar@103
  1400
              }
alpar@103
  1401
              push(new_id, right_id);
alpar@103
  1402
              pushRight(new_id, ~(classes[r].parent));
alpar@103
  1403
              setPrio(new_id);
alpar@103
  1404
alpar@103
  1405
              id = nodes[id].parent;
alpar@103
  1406
              classes[r].parent = ~new_id;
alpar@103
  1407
            }
alpar@103
  1408
            if (id < 0) {
alpar@103
  1409
              int new_parent = newNode();
alpar@103
  1410
              nodes[new_parent].next = -1;
alpar@103
  1411
              nodes[new_parent].prev = -1;
alpar@103
  1412
              nodes[new_parent].parent = ~l;
alpar@103
  1413
alpar@103
  1414
              push(new_parent, ~(classes[l].parent));
alpar@103
  1415
              pushRight(new_parent, ~(classes[r].parent));
alpar@103
  1416
              setPrio(new_parent);
alpar@103
  1417
alpar@103
  1418
              classes[l].parent = ~new_parent;
alpar@103
  1419
              classes[l].depth += 1;
alpar@103
  1420
            } else {
alpar@103
  1421
              pushRight(id, ~(classes[r].parent));
alpar@103
  1422
              while (id >= 0 && less(~(classes[r].parent), id)) {
alpar@103
  1423
                nodes[id].prio = nodes[~(classes[r].parent)].prio;
alpar@103
  1424
                nodes[id].item = nodes[~(classes[r].parent)].item;
alpar@103
  1425
                id = nodes[id].parent;
alpar@103
  1426
              }
alpar@103
  1427
            }
alpar@103
  1428
          } else if (classes[r].depth > classes[l].depth) {
alpar@103
  1429
            int id = ~(classes[r].parent);
alpar@103
  1430
            for (int i = classes[l].depth + 1; i < classes[r].depth; ++i) {
alpar@103
  1431
              id = nodes[id].left;
alpar@103
  1432
            }
alpar@103
  1433
            while (id >= 0 && nodes[id].size == cmax) {
alpar@103
  1434
              int new_id = newNode();
alpar@103
  1435
              int left_id = nodes[id].left;
alpar@103
  1436
alpar@103
  1437
              popLeft(id);
alpar@103
  1438
              if (nodes[id].prio == nodes[left_id].prio) {
alpar@103
  1439
                setPrio(id);
alpar@103
  1440
              }
alpar@103
  1441
              push(new_id, left_id);
alpar@103
  1442
              pushLeft(new_id, ~(classes[l].parent));
alpar@103
  1443
              setPrio(new_id);
alpar@103
  1444
alpar@103
  1445
              id = nodes[id].parent;
alpar@103
  1446
              classes[l].parent = ~new_id;
alpar@103
  1447
alpar@103
  1448
            }
alpar@103
  1449
            if (id < 0) {
alpar@103
  1450
              int new_parent = newNode();
alpar@103
  1451
              nodes[new_parent].next = -1;
alpar@103
  1452
              nodes[new_parent].prev = -1;
alpar@103
  1453
              nodes[new_parent].parent = ~l;
alpar@103
  1454
alpar@103
  1455
              push(new_parent, ~(classes[r].parent));
alpar@103
  1456
              pushLeft(new_parent, ~(classes[l].parent));
alpar@103
  1457
              setPrio(new_parent);
alpar@103
  1458
              
alpar@103
  1459
              classes[r].parent = ~new_parent;
alpar@103
  1460
              classes[r].depth += 1;
alpar@103
  1461
            } else {
alpar@103
  1462
              pushLeft(id, ~(classes[l].parent));
alpar@103
  1463
              while (id >= 0 && less(~(classes[l].parent), id)) {
alpar@103
  1464
                nodes[id].prio = nodes[~(classes[l].parent)].prio;
alpar@103
  1465
                nodes[id].item = nodes[~(classes[l].parent)].item;
alpar@103
  1466
                id = nodes[id].parent;
alpar@103
  1467
              }
alpar@103
  1468
            }
alpar@103
  1469
            nodes[~(classes[r].parent)].parent = ~l;
alpar@103
  1470
            classes[l].parent = classes[r].parent;
alpar@103
  1471
            classes[l].depth = classes[r].depth;
alpar@103
  1472
          } else {
alpar@103
  1473
            if (classes[l].depth != 0 && 
alpar@103
  1474
                nodes[~(classes[l].parent)].size + 
alpar@103
  1475
                nodes[~(classes[r].parent)].size <= cmax) {
alpar@103
  1476
              splice(~(classes[l].parent), ~(classes[r].parent));
alpar@103
  1477
              deleteNode(~(classes[r].parent));
alpar@103
  1478
              if (less(~(classes[r].parent), ~(classes[l].parent))) {
alpar@103
  1479
                nodes[~(classes[l].parent)].prio = 
alpar@103
  1480
                  nodes[~(classes[r].parent)].prio;
alpar@103
  1481
                nodes[~(classes[l].parent)].item = 
alpar@103
  1482
                  nodes[~(classes[r].parent)].item;
alpar@103
  1483
              }
alpar@103
  1484
            } else {
alpar@103
  1485
              int new_parent = newNode();
alpar@103
  1486
              nodes[new_parent].next = nodes[new_parent].prev = -1;
alpar@103
  1487
              push(new_parent, ~(classes[l].parent));
alpar@103
  1488
              pushRight(new_parent, ~(classes[r].parent));
alpar@103
  1489
              setPrio(new_parent);
alpar@103
  1490
            
alpar@103
  1491
              classes[l].parent = ~new_parent;
alpar@103
  1492
              classes[l].depth += 1;
alpar@103
  1493
              nodes[new_parent].parent = ~l;
alpar@103
  1494
            }
alpar@103
  1495
          }
alpar@103
  1496
          if (classes[r].next != -1) {
alpar@103
  1497
            classes[classes[r].next].prev = classes[r].prev;
alpar@103
  1498
          }
alpar@103
  1499
          classes[classes[r].prev].next = classes[r].next;
alpar@103
  1500
alpar@103
  1501
          classes[r].prev = classes[l].right;
alpar@103
  1502
          classes[classes[l].right].next = r;
alpar@103
  1503
          classes[l].right = r;
alpar@103
  1504
          classes[r].parent = l;
alpar@103
  1505
alpar@103
  1506
          classes[r].next = -1;
alpar@103
  1507
          classes[r].depth = rln;
alpar@103
  1508
        }
alpar@103
  1509
      }
alpar@103
  1510
      return class_id;
alpar@103
  1511
    }
alpar@103
  1512
alpar@103
  1513
    /// \brief Split the class to subclasses.
alpar@103
  1514
    ///
alpar@103
  1515
    /// The current function splits the given class. The join, which
alpar@103
  1516
    /// made the current class, stored a reference to the
alpar@103
  1517
    /// subclasses. The \c splitClass() member restores the classes
alpar@103
  1518
    /// and creates the heaps. The parameter is an STL output iterator
alpar@103
  1519
    /// which will be filled with the subclass ids. The time
alpar@103
  1520
    /// complexity is O(log(n)*k) where n is the overall number of
alpar@103
  1521
    /// nodes in the splitted classes and k is the number of the
alpar@103
  1522
    /// classes.
alpar@103
  1523
    template <typename Iterator>
alpar@103
  1524
    void split(int cls, Iterator out) {
alpar@103
  1525
      std::vector<int> cs;
alpar@103
  1526
      { // splitting union-find
alpar@103
  1527
        int id = cls;
alpar@103
  1528
        int l = classes[id].left;
alpar@103
  1529
alpar@103
  1530
        classes[l].parent = classes[id].parent;
alpar@103
  1531
        classes[l].depth = classes[id].depth;
alpar@103
  1532
alpar@103
  1533
        nodes[~(classes[l].parent)].parent = ~l;
alpar@103
  1534
alpar@103
  1535
        *out++ = l;
alpar@103
  1536
alpar@103
  1537
        while (l != -1) {
alpar@103
  1538
          cs.push_back(l);
alpar@103
  1539
          l = classes[l].next;
alpar@103
  1540
        }
alpar@103
  1541
alpar@103
  1542
        classes[classes[id].right].next = first_class;
alpar@103
  1543
        classes[first_class].prev = classes[id].right;
alpar@103
  1544
        first_class = classes[id].left;
alpar@103
  1545
        
alpar@103
  1546
        if (classes[id].next != -1) {
alpar@103
  1547
          classes[classes[id].next].prev = classes[id].prev;
alpar@103
  1548
        }
alpar@103
  1549
        classes[classes[id].prev].next = classes[id].next;
alpar@103
  1550
        
alpar@103
  1551
        deleteClass(id);
alpar@103
  1552
      }
alpar@103
  1553
alpar@103
  1554
      {
alpar@103
  1555
        for (int i = 1; i < int(cs.size()); ++i) {
alpar@103
  1556
          int l = classes[cs[i]].depth;
alpar@103
  1557
          while (nodes[nodes[l].parent].left == l) {
alpar@103
  1558
            l = nodes[l].parent;
alpar@103
  1559
          }
alpar@103
  1560
          int r = l;    
alpar@103
  1561
          while (nodes[l].parent >= 0) {
alpar@103
  1562
            l = nodes[l].parent;
alpar@103
  1563
            int new_node = newNode();
alpar@103
  1564
alpar@103
  1565
            nodes[new_node].prev = -1;
alpar@103
  1566
            nodes[new_node].next = -1;
alpar@103
  1567
alpar@103
  1568
            split(r, new_node);
alpar@103
  1569
            pushAfter(l, new_node);
alpar@103
  1570
            setPrio(l);
alpar@103
  1571
            setPrio(new_node);
alpar@103
  1572
            r = new_node;
alpar@103
  1573
          }
alpar@103
  1574
          classes[cs[i]].parent = ~r;
alpar@103
  1575
          classes[cs[i]].depth = classes[~(nodes[l].parent)].depth;
alpar@103
  1576
          nodes[r].parent = ~cs[i];
alpar@103
  1577
alpar@103
  1578
          nodes[l].next = -1;
alpar@103
  1579
          nodes[r].prev = -1;
alpar@103
  1580
alpar@103
  1581
          repairRight(~(nodes[l].parent));
alpar@103
  1582
          repairLeft(cs[i]);
alpar@103
  1583
          
alpar@103
  1584
          *out++ = cs[i];
alpar@103
  1585
        }
alpar@103
  1586
      }
alpar@103
  1587
    }
alpar@103
  1588
alpar@103
  1589
    /// \brief Gives back the priority of the current item.
alpar@103
  1590
    ///
alpar@103
  1591
    /// \return Gives back the priority of the current item.
alpar@103
  1592
    const Value& operator[](const Item& item) const {
alpar@103
  1593
      return nodes[index[item]].prio;
alpar@103
  1594
    }
alpar@103
  1595
alpar@103
  1596
    /// \brief Sets the priority of the current item.
alpar@103
  1597
    ///
alpar@103
  1598
    /// Sets the priority of the current item.
alpar@103
  1599
    void set(const Item& item, const Value& prio) {
alpar@103
  1600
      if (comp(prio, nodes[index[item]].prio)) {
alpar@103
  1601
        decrease(item, prio);
alpar@103
  1602
      } else if (!comp(prio, nodes[index[item]].prio)) {
alpar@103
  1603
        increase(item, prio);
alpar@103
  1604
      }
alpar@103
  1605
    }
alpar@103
  1606
      
alpar@103
  1607
    /// \brief Increase the priority of the current item.
alpar@103
  1608
    ///
alpar@103
  1609
    /// Increase the priority of the current item.
alpar@103
  1610
    void increase(const Item& item, const Value& prio) {
alpar@103
  1611
      int id = index[item];
alpar@103
  1612
      int kd = nodes[id].parent;
alpar@103
  1613
      nodes[id].prio = prio;
alpar@103
  1614
      while (kd >= 0 && nodes[kd].item == item) {
alpar@103
  1615
        setPrio(kd);
alpar@103
  1616
        kd = nodes[kd].parent;
alpar@103
  1617
      }
alpar@103
  1618
    }
alpar@103
  1619
alpar@103
  1620
    /// \brief Increase the priority of the current item.
alpar@103
  1621
    ///
alpar@103
  1622
    /// Increase the priority of the current item.
alpar@103
  1623
    void decrease(const Item& item, const Value& prio) {
alpar@103
  1624
      int id = index[item];
alpar@103
  1625
      int kd = nodes[id].parent;
alpar@103
  1626
      nodes[id].prio = prio;
alpar@103
  1627
      while (kd >= 0 && less(id, kd)) {
alpar@103
  1628
        nodes[kd].prio = prio;
alpar@103
  1629
        nodes[kd].item = item;
alpar@103
  1630
        kd = nodes[kd].parent;
alpar@103
  1631
      }
alpar@103
  1632
    }
alpar@103
  1633
    
alpar@103
  1634
    /// \brief Gives back the minimum priority of the class.
alpar@103
  1635
    ///
alpar@103
  1636
    /// \return Gives back the minimum priority of the class.
alpar@103
  1637
    const Value& classPrio(int cls) const {
alpar@103
  1638
      return nodes[~(classes[cls].parent)].prio;
alpar@103
  1639
    }
alpar@103
  1640
alpar@103
  1641
    /// \brief Gives back the minimum priority item of the class.
alpar@103
  1642
    ///
alpar@103
  1643
    /// \return Gives back the minimum priority item of the class.
alpar@103
  1644
    const Item& classTop(int cls) const {
alpar@103
  1645
      return nodes[~(classes[cls].parent)].item;
alpar@103
  1646
    }
alpar@103
  1647
alpar@103
  1648
    /// \brief Gives back a representant item of the class.
alpar@103
  1649
    /// 
alpar@103
  1650
    /// The representant is indpendent from the priorities of the
alpar@103
  1651
    /// items. 
alpar@103
  1652
    /// \return Gives back a representant item of the class.
alpar@103
  1653
    const Item& classRep(int id) const {
alpar@103
  1654
      int parent = classes[id].parent;
alpar@103
  1655
      return nodes[parent >= 0 ? classes[id].depth : leftNode(id)].item;
alpar@103
  1656
    }
alpar@103
  1657
alpar@103
  1658
    /// \brief Lemon style iterator for the items of a class.
alpar@103
  1659
    ///
alpar@103
  1660
    /// ClassIt is a lemon style iterator for the components. It iterates
alpar@103
  1661
    /// on the items of a class. By example if you want to iterate on
alpar@103
  1662
    /// each items of each classes then you may write the next code.
alpar@103
  1663
    ///\code
alpar@103
  1664
    /// for (ClassIt cit(huf); cit != INVALID; ++cit) {
alpar@103
  1665
    ///   std::cout << "Class: ";
alpar@103
  1666
    ///   for (ItemIt iit(huf, cit); iit != INVALID; ++iit) {
alpar@103
  1667
    ///     std::cout << toString(iit) << ' ' << std::endl;
alpar@103
  1668
    ///   }
alpar@103
  1669
    ///   std::cout << std::endl;
alpar@103
  1670
    /// }
alpar@103
  1671
    ///\endcode
alpar@103
  1672
    class ItemIt {
alpar@103
  1673
    private:
alpar@103
  1674
alpar@103
  1675
      const HeapUnionFind* _huf;
alpar@103
  1676
      int _id, _lid;
alpar@103
  1677
      
alpar@103
  1678
    public:
alpar@103
  1679
alpar@103
  1680
      /// \brief Default constructor 
alpar@103
  1681
      ///
alpar@103
  1682
      /// Default constructor 
alpar@103
  1683
      ItemIt() {}
alpar@103
  1684
alpar@103
  1685
      ItemIt(const HeapUnionFind& huf, int cls) : _huf(&huf) {
alpar@103
  1686
        int id = cls;
alpar@103
  1687
        int parent = _huf->classes[id].parent;
alpar@103
  1688
        if (parent >= 0) {
alpar@103
  1689
          _id = _huf->classes[id].depth;
alpar@103
  1690
          if (_huf->classes[id].next != -1) {
alpar@103
  1691
            _lid = _huf->classes[_huf->classes[id].next].depth;
alpar@103
  1692
          } else {
alpar@103
  1693
            _lid = -1;
alpar@103
  1694
          }
alpar@103
  1695
        } else {
alpar@103
  1696
          _id = _huf->leftNode(id);
alpar@103
  1697
          _lid = -1;
alpar@103
  1698
        } 
alpar@103
  1699
      }
alpar@103
  1700
      
alpar@103
  1701
      /// \brief Increment operator
alpar@103
  1702
      ///
alpar@103
  1703
      /// It steps to the next item in the class.
alpar@103
  1704
      ItemIt& operator++() {
alpar@103
  1705
        _id = _huf->nextNode(_id);
alpar@103
  1706
        return *this;
alpar@103
  1707
      }
alpar@103
  1708
alpar@103
  1709
      /// \brief Conversion operator
alpar@103
  1710
      ///
alpar@103
  1711
      /// It converts the iterator to the current item.
alpar@103
  1712
      operator const Item&() const {
alpar@103
  1713
        return _huf->nodes[_id].item;
alpar@103
  1714
      }
alpar@103
  1715
      
alpar@103
  1716
      /// \brief Equality operator
alpar@103
  1717
      ///
alpar@103
  1718
      /// Equality operator
alpar@103
  1719
      bool operator==(const ItemIt& i) { 
alpar@103
  1720
        return i._id == _id;
alpar@103
  1721
      }
alpar@103
  1722
alpar@103
  1723
      /// \brief Inequality operator
alpar@103
  1724
      ///
alpar@103
  1725
      /// Inequality operator
alpar@103
  1726
      bool operator!=(const ItemIt& i) { 
alpar@103
  1727
        return i._id != _id;
alpar@103
  1728
      }
alpar@103
  1729
alpar@103
  1730
      /// \brief Equality operator
alpar@103
  1731
      ///
alpar@103
  1732
      /// Equality operator
alpar@103
  1733
      bool operator==(Invalid) { 
alpar@103
  1734
        return _id == _lid;
alpar@103
  1735
      }
alpar@103
  1736
alpar@103
  1737
      /// \brief Inequality operator
alpar@103
  1738
      ///
alpar@103
  1739
      /// Inequality operator
alpar@103
  1740
      bool operator!=(Invalid) { 
alpar@103
  1741
        return _id != _lid;
alpar@103
  1742
      }
alpar@103
  1743
      
alpar@103
  1744
    };
alpar@103
  1745
alpar@103
  1746
    /// \brief Class iterator
alpar@103
  1747
    ///
alpar@103
  1748
    /// The iterator stores 
alpar@103
  1749
    class ClassIt {
alpar@103
  1750
    private:
alpar@103
  1751
alpar@103
  1752
      const HeapUnionFind* _huf;
alpar@103
  1753
      int _id;
alpar@103
  1754
alpar@103
  1755
    public:
alpar@103
  1756
alpar@103
  1757
      ClassIt(const HeapUnionFind& huf) 
alpar@103
  1758
        : _huf(&huf), _id(huf.first_class) {}
alpar@103
  1759
alpar@103
  1760
      ClassIt(const HeapUnionFind& huf, int cls) 
alpar@103
  1761
        : _huf(&huf), _id(huf.classes[cls].left) {}
alpar@103
  1762
alpar@103
  1763
      ClassIt(Invalid) : _huf(0), _id(-1) {}
alpar@103
  1764
      
alpar@103
  1765
      const ClassIt& operator++() {
alpar@103
  1766
        _id = _huf->classes[_id].next;
alpar@103
  1767
	return *this;
alpar@103
  1768
      }
alpar@103
  1769
alpar@103
  1770
      /// \brief Equality operator
alpar@103
  1771
      ///
alpar@103
  1772
      /// Equality operator
alpar@103
  1773
      bool operator==(const ClassIt& i) { 
alpar@103
  1774
        return i._id == _id;
alpar@103
  1775
      }
alpar@103
  1776
alpar@103
  1777
      /// \brief Inequality operator
alpar@103
  1778
      ///
alpar@103
  1779
      /// Inequality operator
alpar@103
  1780
      bool operator!=(const ClassIt& i) { 
alpar@103
  1781
        return i._id != _id;
alpar@103
  1782
      }      
alpar@103
  1783
      
alpar@103
  1784
      operator int() const {
alpar@103
  1785
	return _id;
alpar@103
  1786
      }
alpar@103
  1787
            
alpar@103
  1788
    };
alpar@103
  1789
alpar@103
  1790
  };
alpar@103
  1791
alpar@103
  1792
  //! @}
alpar@103
  1793
alpar@103
  1794
} //namespace lemon
alpar@103
  1795
alpar@103
  1796
#endif //LEMON_UNION_FIND_H