lemon/fib_heap.h
author Peter Kovacs <kpeter@inf.elte.hu>
Tue, 15 Mar 2011 19:32:21 +0100
changeset 936 ddd3c0d3d9bf
parent 710 f1fe0ddad6f7
permissions -rw-r--r--
Implement the scaling Price Refinement heuristic in CostScaling (#417)
instead of Early Termination.

These two heuristics are similar, but the newer one is faster
and not only makes it possible to skip some epsilon phases, but
it can improve the performance of the other phases, as well.
deba@681
     1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
deba@681
     2
 *
deba@681
     3
 * This file is a part of LEMON, a generic C++ optimization library.
deba@681
     4
 *
deba@681
     5
 * Copyright (C) 2003-2009
deba@681
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
deba@681
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
deba@681
     8
 *
deba@681
     9
 * Permission to use, modify and distribute this software is granted
deba@681
    10
 * provided that this copyright notice appears in all copies. For
deba@681
    11
 * precise terms see the accompanying LICENSE file.
deba@681
    12
 *
deba@681
    13
 * This software is provided "AS IS" with no warranty of any kind,
deba@681
    14
 * express or implied, and with no claim as to its suitability for any
deba@681
    15
 * purpose.
deba@681
    16
 *
deba@681
    17
 */
deba@681
    18
deba@681
    19
#ifndef LEMON_FIB_HEAP_H
deba@681
    20
#define LEMON_FIB_HEAP_H
deba@681
    21
deba@681
    22
///\file
kpeter@710
    23
///\ingroup heaps
kpeter@709
    24
///\brief Fibonacci heap implementation.
deba@681
    25
deba@681
    26
#include <vector>
kpeter@709
    27
#include <utility>
deba@681
    28
#include <functional>
deba@681
    29
#include <lemon/math.h>
deba@681
    30
deba@681
    31
namespace lemon {
deba@681
    32
kpeter@710
    33
  /// \ingroup heaps
deba@681
    34
  ///
kpeter@709
    35
  /// \brief Fibonacci heap data structure.
deba@681
    36
  ///
kpeter@709
    37
  /// This class implements the \e Fibonacci \e heap data structure.
kpeter@709
    38
  /// It fully conforms to the \ref concepts::Heap "heap concept".
deba@681
    39
  ///
kpeter@709
    40
  /// The methods \ref increase() and \ref erase() are not efficient in a
kpeter@709
    41
  /// Fibonacci heap. In case of many calls of these operations, it is
kpeter@709
    42
  /// better to use other heap structure, e.g. \ref BinHeap "binary heap".
deba@681
    43
  ///
kpeter@709
    44
  /// \tparam PR Type of the priorities of the items.
kpeter@709
    45
  /// \tparam IM A read-writable item map with \c int values, used
kpeter@709
    46
  /// internally to handle the cross references.
kpeter@709
    47
  /// \tparam CMP A functor class for comparing the priorities.
kpeter@709
    48
  /// The default is \c std::less<PR>.
deba@681
    49
#ifdef DOXYGEN
kpeter@709
    50
  template <typename PR, typename IM, typename CMP>
deba@681
    51
#else
kpeter@709
    52
  template <typename PR, typename IM, typename CMP = std::less<PR> >
deba@681
    53
#endif
deba@681
    54
  class FibHeap {
deba@681
    55
  public:
kpeter@709
    56
kpeter@709
    57
    /// Type of the item-int map.
deba@683
    58
    typedef IM ItemIntMap;
kpeter@709
    59
    /// Type of the priorities.
kpeter@709
    60
    typedef PR Prio;
kpeter@709
    61
    /// Type of the items stored in the heap.
deba@681
    62
    typedef typename ItemIntMap::Key Item;
kpeter@709
    63
    /// Type of the item-priority pairs.
deba@681
    64
    typedef std::pair<Item,Prio> Pair;
kpeter@709
    65
    /// Functor type for comparing the priorities.
deba@683
    66
    typedef CMP Compare;
deba@681
    67
deba@681
    68
  private:
deba@683
    69
    class Store;
deba@681
    70
deba@683
    71
    std::vector<Store> _data;
deba@683
    72
    int _minimum;
deba@683
    73
    ItemIntMap &_iim;
deba@683
    74
    Compare _comp;
deba@683
    75
    int _num;
deba@681
    76
deba@681
    77
  public:
deba@683
    78
kpeter@709
    79
    /// \brief Type to represent the states of the items.
deba@683
    80
    ///
kpeter@709
    81
    /// Each item has a state associated to it. It can be "in heap",
kpeter@709
    82
    /// "pre-heap" or "post-heap". The latter two are indifferent from the
deba@683
    83
    /// heap's point of view, but may be useful to the user.
deba@683
    84
    ///
deba@683
    85
    /// The item-int map must be initialized in such way that it assigns
deba@683
    86
    /// \c PRE_HEAP (<tt>-1</tt>) to any element to be put in the heap.
deba@681
    87
    enum State {
deba@683
    88
      IN_HEAP = 0,    ///< = 0.
deba@683
    89
      PRE_HEAP = -1,  ///< = -1.
deba@683
    90
      POST_HEAP = -2  ///< = -2.
deba@681
    91
    };
deba@681
    92
kpeter@709
    93
    /// \brief Constructor.
deba@681
    94
    ///
kpeter@709
    95
    /// Constructor.
kpeter@709
    96
    /// \param map A map that assigns \c int values to the items.
kpeter@709
    97
    /// It is used internally to handle the cross references.
kpeter@709
    98
    /// The assigned value must be \c PRE_HEAP (<tt>-1</tt>) for each item.
deba@683
    99
    explicit FibHeap(ItemIntMap &map)
deba@683
   100
      : _minimum(0), _iim(map), _num() {}
deba@681
   101
kpeter@709
   102
    /// \brief Constructor.
deba@681
   103
    ///
kpeter@709
   104
    /// Constructor.
kpeter@709
   105
    /// \param map A map that assigns \c int values to the items.
kpeter@709
   106
    /// It is used internally to handle the cross references.
kpeter@709
   107
    /// The assigned value must be \c PRE_HEAP (<tt>-1</tt>) for each item.
kpeter@709
   108
    /// \param comp The function object used for comparing the priorities.
deba@683
   109
    FibHeap(ItemIntMap &map, const Compare &comp)
deba@683
   110
      : _minimum(0), _iim(map), _comp(comp), _num() {}
deba@681
   111
deba@681
   112
    /// \brief The number of items stored in the heap.
deba@681
   113
    ///
kpeter@709
   114
    /// This function returns the number of items stored in the heap.
deba@683
   115
    int size() const { return _num; }
deba@681
   116
kpeter@709
   117
    /// \brief Check if the heap is empty.
deba@681
   118
    ///
kpeter@709
   119
    /// This function returns \c true if the heap is empty.
deba@683
   120
    bool empty() const { return _num==0; }
deba@681
   121
kpeter@709
   122
    /// \brief Make the heap empty.
deba@681
   123
    ///
kpeter@709
   124
    /// This functon makes the heap empty.
kpeter@709
   125
    /// It does not change the cross reference map. If you want to reuse
kpeter@709
   126
    /// a heap that is not surely empty, you should first clear it and
kpeter@709
   127
    /// then you should set the cross reference map to \c PRE_HEAP
kpeter@709
   128
    /// for each item.
deba@681
   129
    void clear() {
deba@683
   130
      _data.clear(); _minimum = 0; _num = 0;
deba@681
   131
    }
deba@681
   132
kpeter@709
   133
    /// \brief Insert an item into the heap with the given priority.
deba@681
   134
    ///
kpeter@709
   135
    /// This function inserts the given item into the heap with the
kpeter@709
   136
    /// given priority.
kpeter@709
   137
    /// \param item The item to insert.
kpeter@709
   138
    /// \param prio The priority of the item.
kpeter@709
   139
    /// \pre \e item must not be stored in the heap.
kpeter@709
   140
    void push (const Item& item, const Prio& prio) {
deba@683
   141
      int i=_iim[item];
deba@681
   142
      if ( i < 0 ) {
deba@683
   143
        int s=_data.size();
deba@683
   144
        _iim.set( item, s );
deba@683
   145
        Store st;
deba@681
   146
        st.name=item;
deba@683
   147
        _data.push_back(st);
deba@681
   148
        i=s;
deba@681
   149
      } else {
deba@683
   150
        _data[i].parent=_data[i].child=-1;
deba@683
   151
        _data[i].degree=0;
deba@683
   152
        _data[i].in=true;
deba@683
   153
        _data[i].marked=false;
deba@681
   154
      }
deba@681
   155
deba@683
   156
      if ( _num ) {
deba@683
   157
        _data[_data[_minimum].right_neighbor].left_neighbor=i;
deba@683
   158
        _data[i].right_neighbor=_data[_minimum].right_neighbor;
deba@683
   159
        _data[_minimum].right_neighbor=i;
deba@683
   160
        _data[i].left_neighbor=_minimum;
kpeter@709
   161
        if ( _comp( prio, _data[_minimum].prio) ) _minimum=i;
deba@681
   162
      } else {
deba@683
   163
        _data[i].right_neighbor=_data[i].left_neighbor=i;
deba@683
   164
        _minimum=i;
deba@681
   165
      }
kpeter@709
   166
      _data[i].prio=prio;
deba@683
   167
      ++_num;
deba@681
   168
    }
deba@681
   169
kpeter@709
   170
    /// \brief Return the item having minimum priority.
deba@681
   171
    ///
kpeter@709
   172
    /// This function returns the item having minimum priority.
kpeter@709
   173
    /// \pre The heap must be non-empty.
deba@683
   174
    Item top() const { return _data[_minimum].name; }
deba@681
   175
kpeter@709
   176
    /// \brief The minimum priority.
deba@681
   177
    ///
kpeter@709
   178
    /// This function returns the minimum priority.
kpeter@709
   179
    /// \pre The heap must be non-empty.
kpeter@709
   180
    Prio prio() const { return _data[_minimum].prio; }
deba@681
   181
kpeter@709
   182
    /// \brief Remove the item having minimum priority.
deba@681
   183
    ///
kpeter@709
   184
    /// This function removes the item having minimum priority.
deba@681
   185
    /// \pre The heap must be non-empty.
deba@681
   186
    void pop() {
deba@681
   187
      /*The first case is that there are only one root.*/
deba@683
   188
      if ( _data[_minimum].left_neighbor==_minimum ) {
deba@683
   189
        _data[_minimum].in=false;
deba@683
   190
        if ( _data[_minimum].degree!=0 ) {
kpeter@711
   191
          makeRoot(_data[_minimum].child);
deba@683
   192
          _minimum=_data[_minimum].child;
deba@681
   193
          balance();
deba@681
   194
        }
deba@681
   195
      } else {
deba@683
   196
        int right=_data[_minimum].right_neighbor;
deba@683
   197
        unlace(_minimum);
deba@683
   198
        _data[_minimum].in=false;
deba@683
   199
        if ( _data[_minimum].degree > 0 ) {
deba@683
   200
          int left=_data[_minimum].left_neighbor;
deba@683
   201
          int child=_data[_minimum].child;
deba@683
   202
          int last_child=_data[child].left_neighbor;
deba@681
   203
kpeter@711
   204
          makeRoot(child);
deba@681
   205
deba@683
   206
          _data[left].right_neighbor=child;
deba@683
   207
          _data[child].left_neighbor=left;
deba@683
   208
          _data[right].left_neighbor=last_child;
deba@683
   209
          _data[last_child].right_neighbor=right;
deba@681
   210
        }
deba@683
   211
        _minimum=right;
deba@681
   212
        balance();
deba@681
   213
      } // the case where there are more roots
deba@683
   214
      --_num;
deba@681
   215
    }
deba@681
   216
kpeter@709
   217
    /// \brief Remove the given item from the heap.
deba@681
   218
    ///
kpeter@709
   219
    /// This function removes the given item from the heap if it is
kpeter@709
   220
    /// already stored.
kpeter@709
   221
    /// \param item The item to delete.
kpeter@709
   222
    /// \pre \e item must be in the heap.
deba@681
   223
    void erase (const Item& item) {
deba@683
   224
      int i=_iim[item];
deba@681
   225
deba@683
   226
      if ( i >= 0 && _data[i].in ) {
deba@683
   227
        if ( _data[i].parent!=-1 ) {
deba@683
   228
          int p=_data[i].parent;
deba@681
   229
          cut(i,p);
deba@681
   230
          cascade(p);
deba@681
   231
        }
deba@683
   232
        _minimum=i;     //As if its prio would be -infinity
deba@681
   233
        pop();
deba@681
   234
      }
deba@681
   235
    }
deba@681
   236
kpeter@709
   237
    /// \brief The priority of the given item.
deba@681
   238
    ///
kpeter@709
   239
    /// This function returns the priority of the given item.
kpeter@709
   240
    /// \param item The item.
kpeter@709
   241
    /// \pre \e item must be in the heap.
kpeter@709
   242
    Prio operator[](const Item& item) const {
kpeter@709
   243
      return _data[_iim[item]].prio;
kpeter@709
   244
    }
kpeter@709
   245
kpeter@709
   246
    /// \brief Set the priority of an item or insert it, if it is
kpeter@709
   247
    /// not stored in the heap.
kpeter@709
   248
    ///
kpeter@709
   249
    /// This method sets the priority of the given item if it is
kpeter@709
   250
    /// already stored in the heap. Otherwise it inserts the given
kpeter@709
   251
    /// item into the heap with the given priority.
kpeter@709
   252
    /// \param item The item.
kpeter@709
   253
    /// \param prio The priority.
kpeter@709
   254
    void set (const Item& item, const Prio& prio) {
deba@683
   255
      int i=_iim[item];
kpeter@709
   256
      if ( i >= 0 && _data[i].in ) {
kpeter@709
   257
        if ( _comp(prio, _data[i].prio) ) decrease(item, prio);
kpeter@709
   258
        if ( _comp(_data[i].prio, prio) ) increase(item, prio);
kpeter@709
   259
      } else push(item, prio);
kpeter@709
   260
    }
kpeter@709
   261
kpeter@709
   262
    /// \brief Decrease the priority of an item to the given value.
kpeter@709
   263
    ///
kpeter@709
   264
    /// This function decreases the priority of an item to the given value.
kpeter@709
   265
    /// \param item The item.
kpeter@709
   266
    /// \param prio The priority.
kpeter@709
   267
    /// \pre \e item must be stored in the heap with priority at least \e prio.
kpeter@709
   268
    void decrease (const Item& item, const Prio& prio) {
kpeter@709
   269
      int i=_iim[item];
kpeter@709
   270
      _data[i].prio=prio;
deba@683
   271
      int p=_data[i].parent;
deba@681
   272
kpeter@709
   273
      if ( p!=-1 && _comp(prio, _data[p].prio) ) {
deba@681
   274
        cut(i,p);
deba@681
   275
        cascade(p);
deba@681
   276
      }
kpeter@709
   277
      if ( _comp(prio, _data[_minimum].prio) ) _minimum=i;
deba@681
   278
    }
deba@681
   279
kpeter@709
   280
    /// \brief Increase the priority of an item to the given value.
deba@681
   281
    ///
kpeter@709
   282
    /// This function increases the priority of an item to the given value.
kpeter@709
   283
    /// \param item The item.
kpeter@709
   284
    /// \param prio The priority.
kpeter@709
   285
    /// \pre \e item must be stored in the heap with priority at most \e prio.
kpeter@709
   286
    void increase (const Item& item, const Prio& prio) {
deba@681
   287
      erase(item);
kpeter@709
   288
      push(item, prio);
deba@681
   289
    }
deba@681
   290
kpeter@709
   291
    /// \brief Return the state of an item.
deba@681
   292
    ///
kpeter@709
   293
    /// This method returns \c PRE_HEAP if the given item has never
kpeter@709
   294
    /// been in the heap, \c IN_HEAP if it is in the heap at the moment,
kpeter@709
   295
    /// and \c POST_HEAP otherwise.
kpeter@709
   296
    /// In the latter case it is possible that the item will get back
kpeter@709
   297
    /// to the heap again.
kpeter@709
   298
    /// \param item The item.
deba@681
   299
    State state(const Item &item) const {
deba@683
   300
      int i=_iim[item];
deba@681
   301
      if( i>=0 ) {
deba@683
   302
        if ( _data[i].in ) i=0;
deba@681
   303
        else i=-2;
deba@681
   304
      }
deba@681
   305
      return State(i);
deba@681
   306
    }
deba@681
   307
kpeter@709
   308
    /// \brief Set the state of an item in the heap.
deba@681
   309
    ///
kpeter@709
   310
    /// This function sets the state of the given item in the heap.
kpeter@709
   311
    /// It can be used to manually clear the heap when it is important
kpeter@709
   312
    /// to achive better time complexity.
deba@681
   313
    /// \param i The item.
deba@681
   314
    /// \param st The state. It should not be \c IN_HEAP.
deba@681
   315
    void state(const Item& i, State st) {
deba@681
   316
      switch (st) {
deba@681
   317
      case POST_HEAP:
deba@681
   318
      case PRE_HEAP:
deba@681
   319
        if (state(i) == IN_HEAP) {
deba@681
   320
          erase(i);
deba@681
   321
        }
deba@683
   322
        _iim[i] = st;
deba@681
   323
        break;
deba@681
   324
      case IN_HEAP:
deba@681
   325
        break;
deba@681
   326
      }
deba@681
   327
    }
deba@681
   328
deba@681
   329
  private:
deba@681
   330
deba@681
   331
    void balance() {
deba@681
   332
deba@683
   333
      int maxdeg=int( std::floor( 2.08*log(double(_data.size()))))+1;
deba@681
   334
deba@681
   335
      std::vector<int> A(maxdeg,-1);
deba@681
   336
deba@681
   337
      /*
deba@681
   338
       *Recall that now minimum does not point to the minimum prio element.
deba@681
   339
       *We set minimum to this during balance().
deba@681
   340
       */
deba@683
   341
      int anchor=_data[_minimum].left_neighbor;
deba@683
   342
      int next=_minimum;
deba@681
   343
      bool end=false;
deba@681
   344
deba@681
   345
      do {
deba@681
   346
        int active=next;
deba@681
   347
        if ( anchor==active ) end=true;
deba@683
   348
        int d=_data[active].degree;
deba@683
   349
        next=_data[active].right_neighbor;
deba@681
   350
deba@681
   351
        while (A[d]!=-1) {
deba@683
   352
          if( _comp(_data[active].prio, _data[A[d]].prio) ) {
deba@681
   353
            fuse(active,A[d]);
deba@681
   354
          } else {
deba@681
   355
            fuse(A[d],active);
deba@681
   356
            active=A[d];
deba@681
   357
          }
deba@681
   358
          A[d]=-1;
deba@681
   359
          ++d;
deba@681
   360
        }
deba@681
   361
        A[d]=active;
deba@681
   362
      } while ( !end );
deba@681
   363
deba@681
   364
deba@683
   365
      while ( _data[_minimum].parent >=0 )
deba@683
   366
        _minimum=_data[_minimum].parent;
deba@683
   367
      int s=_minimum;
deba@683
   368
      int m=_minimum;
deba@681
   369
      do {
deba@683
   370
        if ( _comp(_data[s].prio, _data[_minimum].prio) ) _minimum=s;
deba@683
   371
        s=_data[s].right_neighbor;
deba@681
   372
      } while ( s != m );
deba@681
   373
    }
deba@681
   374
kpeter@711
   375
    void makeRoot(int c) {
deba@681
   376
      int s=c;
deba@681
   377
      do {
deba@683
   378
        _data[s].parent=-1;
deba@683
   379
        s=_data[s].right_neighbor;
deba@681
   380
      } while ( s != c );
deba@681
   381
    }
deba@681
   382
deba@681
   383
    void cut(int a, int b) {
deba@681
   384
      /*
deba@681
   385
       *Replacing a from the children of b.
deba@681
   386
       */
deba@683
   387
      --_data[b].degree;
deba@681
   388
deba@683
   389
      if ( _data[b].degree !=0 ) {
deba@683
   390
        int child=_data[b].child;
deba@681
   391
        if ( child==a )
deba@683
   392
          _data[b].child=_data[child].right_neighbor;
deba@681
   393
        unlace(a);
deba@681
   394
      }
deba@681
   395
deba@681
   396
deba@681
   397
      /*Lacing a to the roots.*/
deba@683
   398
      int right=_data[_minimum].right_neighbor;
deba@683
   399
      _data[_minimum].right_neighbor=a;
deba@683
   400
      _data[a].left_neighbor=_minimum;
deba@683
   401
      _data[a].right_neighbor=right;
deba@683
   402
      _data[right].left_neighbor=a;
deba@681
   403
deba@683
   404
      _data[a].parent=-1;
deba@683
   405
      _data[a].marked=false;
deba@681
   406
    }
deba@681
   407
deba@681
   408
    void cascade(int a) {
deba@683
   409
      if ( _data[a].parent!=-1 ) {
deba@683
   410
        int p=_data[a].parent;
deba@681
   411
deba@683
   412
        if ( _data[a].marked==false ) _data[a].marked=true;
deba@681
   413
        else {
deba@681
   414
          cut(a,p);
deba@681
   415
          cascade(p);
deba@681
   416
        }
deba@681
   417
      }
deba@681
   418
    }
deba@681
   419
deba@681
   420
    void fuse(int a, int b) {
deba@681
   421
      unlace(b);
deba@681
   422
deba@681
   423
      /*Lacing b under a.*/
deba@683
   424
      _data[b].parent=a;
deba@681
   425
deba@683
   426
      if (_data[a].degree==0) {
deba@683
   427
        _data[b].left_neighbor=b;
deba@683
   428
        _data[b].right_neighbor=b;
deba@683
   429
        _data[a].child=b;
deba@681
   430
      } else {
deba@683
   431
        int child=_data[a].child;
deba@683
   432
        int last_child=_data[child].left_neighbor;
deba@683
   433
        _data[child].left_neighbor=b;
deba@683
   434
        _data[b].right_neighbor=child;
deba@683
   435
        _data[last_child].right_neighbor=b;
deba@683
   436
        _data[b].left_neighbor=last_child;
deba@681
   437
      }
deba@681
   438
deba@683
   439
      ++_data[a].degree;
deba@681
   440
deba@683
   441
      _data[b].marked=false;
deba@681
   442
    }
deba@681
   443
deba@681
   444
    /*
deba@681
   445
     *It is invoked only if a has siblings.
deba@681
   446
     */
deba@681
   447
    void unlace(int a) {
deba@683
   448
      int leftn=_data[a].left_neighbor;
deba@683
   449
      int rightn=_data[a].right_neighbor;
deba@683
   450
      _data[leftn].right_neighbor=rightn;
deba@683
   451
      _data[rightn].left_neighbor=leftn;
deba@681
   452
    }
deba@681
   453
deba@681
   454
deba@683
   455
    class Store {
deba@681
   456
      friend class FibHeap;
deba@681
   457
deba@681
   458
      Item name;
deba@681
   459
      int parent;
deba@681
   460
      int left_neighbor;
deba@681
   461
      int right_neighbor;
deba@681
   462
      int child;
deba@681
   463
      int degree;
deba@681
   464
      bool marked;
deba@681
   465
      bool in;
deba@681
   466
      Prio prio;
deba@681
   467
deba@683
   468
      Store() : parent(-1), child(-1), degree(), marked(false), in(true) {}
deba@681
   469
    };
deba@681
   470
  };
deba@681
   471
deba@681
   472
} //namespace lemon
deba@681
   473
deba@681
   474
#endif //LEMON_FIB_HEAP_H
deba@681
   475