lemon/hao_orlin.h
author Balazs Dezso <deba@inf.elte.hu>
Sun, 14 Nov 2010 16:35:31 +0100
changeset 1018 2e959a5a0c2d
parent 877 141f9c0db4a3
child 1092 dceba191c00d
permissions -rw-r--r--
Add bipartite graph concepts (#69)
deba@409
     1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
deba@409
     2
 *
deba@409
     3
 * This file is a part of LEMON, a generic C++ optimization library.
deba@409
     4
 *
alpar@877
     5
 * Copyright (C) 2003-2010
deba@409
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
deba@409
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
deba@409
     8
 *
deba@409
     9
 * Permission to use, modify and distribute this software is granted
deba@409
    10
 * provided that this copyright notice appears in all copies. For
deba@409
    11
 * precise terms see the accompanying LICENSE file.
deba@409
    12
 *
deba@409
    13
 * This software is provided "AS IS" with no warranty of any kind,
deba@409
    14
 * express or implied, and with no claim as to its suitability for any
deba@409
    15
 * purpose.
deba@409
    16
 *
deba@409
    17
 */
deba@409
    18
deba@409
    19
#ifndef LEMON_HAO_ORLIN_H
deba@409
    20
#define LEMON_HAO_ORLIN_H
deba@409
    21
deba@409
    22
#include <vector>
deba@409
    23
#include <list>
deba@409
    24
#include <limits>
deba@409
    25
deba@409
    26
#include <lemon/maps.h>
deba@409
    27
#include <lemon/core.h>
deba@409
    28
#include <lemon/tolerance.h>
deba@409
    29
deba@409
    30
/// \file
deba@409
    31
/// \ingroup min_cut
deba@409
    32
/// \brief Implementation of the Hao-Orlin algorithm.
deba@409
    33
///
alpar@877
    34
/// Implementation of the Hao-Orlin algorithm for finding a minimum cut
kpeter@596
    35
/// in a digraph.
deba@409
    36
deba@409
    37
namespace lemon {
deba@409
    38
deba@409
    39
  /// \ingroup min_cut
deba@409
    40
  ///
kpeter@596
    41
  /// \brief Hao-Orlin algorithm for finding a minimum cut in a digraph.
deba@409
    42
  ///
kpeter@596
    43
  /// This class implements the Hao-Orlin algorithm for finding a minimum
alpar@877
    44
  /// value cut in a directed graph \f$D=(V,A)\f$.
kpeter@596
    45
  /// It takes a fixed node \f$ source \in V \f$ and
deba@409
    46
  /// consists of two phases: in the first phase it determines a
deba@409
    47
  /// minimum cut with \f$ source \f$ on the source-side (i.e. a set
kpeter@596
    48
  /// \f$ X\subsetneq V \f$ with \f$ source \in X \f$ and minimal outgoing
kpeter@596
    49
  /// capacity) and in the second phase it determines a minimum cut
deba@409
    50
  /// with \f$ source \f$ on the sink-side (i.e. a set
kpeter@596
    51
  /// \f$ X\subsetneq V \f$ with \f$ source \notin X \f$ and minimal outgoing
kpeter@596
    52
  /// capacity). Obviously, the smaller of these two cuts will be a
deba@409
    53
  /// minimum cut of \f$ D \f$. The algorithm is a modified
kpeter@596
    54
  /// preflow push-relabel algorithm. Our implementation calculates
deba@409
    55
  /// the minimum cut in \f$ O(n^2\sqrt{m}) \f$ time (we use the
alpar@915
    56
  /// highest-label rule), or in \f$O(nm)\f$ for unit capacities. A notable
alpar@915
    57
  /// use of this algorithm is testing network reliability.
deba@409
    58
  ///
kpeter@596
    59
  /// For an undirected graph you can run just the first phase of the
kpeter@596
    60
  /// algorithm or you can use the algorithm of Nagamochi and Ibaraki,
alpar@877
    61
  /// which solves the undirected problem in \f$ O(nm + n^2 \log n) \f$
kpeter@596
    62
  /// time. It is implemented in the NagamochiIbaraki algorithm class.
kpeter@596
    63
  ///
kpeter@596
    64
  /// \tparam GR The type of the digraph the algorithm runs on.
kpeter@596
    65
  /// \tparam CAP The type of the arc map containing the capacities,
kpeter@596
    66
  /// which can be any numreric type. The default map type is
kpeter@596
    67
  /// \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
kpeter@596
    68
  /// \tparam TOL Tolerance class for handling inexact computations. The
kpeter@559
    69
  /// default tolerance type is \ref Tolerance "Tolerance<CAP::Value>".
deba@409
    70
#ifdef DOXYGEN
kpeter@559
    71
  template <typename GR, typename CAP, typename TOL>
deba@409
    72
#else
kpeter@559
    73
  template <typename GR,
kpeter@559
    74
            typename CAP = typename GR::template ArcMap<int>,
kpeter@559
    75
            typename TOL = Tolerance<typename CAP::Value> >
deba@409
    76
#endif
deba@409
    77
  class HaoOrlin {
kpeter@596
    78
  public:
alpar@877
    79
kpeter@596
    80
    /// The digraph type of the algorithm
kpeter@596
    81
    typedef GR Digraph;
kpeter@596
    82
    /// The capacity map type of the algorithm
kpeter@596
    83
    typedef CAP CapacityMap;
kpeter@596
    84
    /// The tolerance type of the algorithm
kpeter@596
    85
    typedef TOL Tolerance;
kpeter@596
    86
deba@409
    87
  private:
deba@409
    88
deba@409
    89
    typedef typename CapacityMap::Value Value;
deba@409
    90
kpeter@596
    91
    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
deba@409
    92
deba@409
    93
    const Digraph& _graph;
deba@409
    94
    const CapacityMap* _capacity;
deba@409
    95
deba@409
    96
    typedef typename Digraph::template ArcMap<Value> FlowMap;
deba@409
    97
    FlowMap* _flow;
deba@409
    98
deba@409
    99
    Node _source;
deba@409
   100
deba@409
   101
    int _node_num;
deba@409
   102
deba@409
   103
    // Bucketing structure
deba@409
   104
    std::vector<Node> _first, _last;
deba@409
   105
    typename Digraph::template NodeMap<Node>* _next;
deba@409
   106
    typename Digraph::template NodeMap<Node>* _prev;
deba@409
   107
    typename Digraph::template NodeMap<bool>* _active;
deba@409
   108
    typename Digraph::template NodeMap<int>* _bucket;
deba@409
   109
deba@409
   110
    std::vector<bool> _dormant;
deba@409
   111
deba@409
   112
    std::list<std::list<int> > _sets;
deba@409
   113
    std::list<int>::iterator _highest;
deba@409
   114
deba@409
   115
    typedef typename Digraph::template NodeMap<Value> ExcessMap;
deba@409
   116
    ExcessMap* _excess;
deba@409
   117
deba@409
   118
    typedef typename Digraph::template NodeMap<bool> SourceSetMap;
deba@409
   119
    SourceSetMap* _source_set;
deba@409
   120
deba@409
   121
    Value _min_cut;
deba@409
   122
deba@409
   123
    typedef typename Digraph::template NodeMap<bool> MinCutMap;
deba@409
   124
    MinCutMap* _min_cut_map;
deba@409
   125
deba@409
   126
    Tolerance _tolerance;
deba@409
   127
deba@409
   128
  public:
deba@409
   129
deba@409
   130
    /// \brief Constructor
deba@409
   131
    ///
deba@409
   132
    /// Constructor of the algorithm class.
deba@409
   133
    HaoOrlin(const Digraph& graph, const CapacityMap& capacity,
deba@409
   134
             const Tolerance& tolerance = Tolerance()) :
deba@409
   135
      _graph(graph), _capacity(&capacity), _flow(0), _source(),
deba@409
   136
      _node_num(), _first(), _last(), _next(0), _prev(0),
deba@409
   137
      _active(0), _bucket(0), _dormant(), _sets(), _highest(),
deba@409
   138
      _excess(0), _source_set(0), _min_cut(), _min_cut_map(0),
deba@409
   139
      _tolerance(tolerance) {}
deba@409
   140
deba@409
   141
    ~HaoOrlin() {
deba@409
   142
      if (_min_cut_map) {
deba@409
   143
        delete _min_cut_map;
deba@409
   144
      }
deba@409
   145
      if (_source_set) {
deba@409
   146
        delete _source_set;
deba@409
   147
      }
deba@409
   148
      if (_excess) {
deba@409
   149
        delete _excess;
deba@409
   150
      }
deba@409
   151
      if (_next) {
deba@409
   152
        delete _next;
deba@409
   153
      }
deba@409
   154
      if (_prev) {
deba@409
   155
        delete _prev;
deba@409
   156
      }
deba@409
   157
      if (_active) {
deba@409
   158
        delete _active;
deba@409
   159
      }
deba@409
   160
      if (_bucket) {
deba@409
   161
        delete _bucket;
deba@409
   162
      }
deba@409
   163
      if (_flow) {
deba@409
   164
        delete _flow;
deba@409
   165
      }
deba@409
   166
    }
deba@409
   167
kpeter@860
   168
    /// \brief Set the tolerance used by the algorithm.
kpeter@860
   169
    ///
kpeter@860
   170
    /// This function sets the tolerance object used by the algorithm.
kpeter@860
   171
    /// \return <tt>(*this)</tt>
kpeter@860
   172
    HaoOrlin& tolerance(const Tolerance& tolerance) {
kpeter@860
   173
      _tolerance = tolerance;
kpeter@860
   174
      return *this;
kpeter@860
   175
    }
kpeter@860
   176
kpeter@860
   177
    /// \brief Returns a const reference to the tolerance.
kpeter@860
   178
    ///
kpeter@860
   179
    /// This function returns a const reference to the tolerance object
kpeter@860
   180
    /// used by the algorithm.
kpeter@860
   181
    const Tolerance& tolerance() const {
kpeter@860
   182
      return _tolerance;
kpeter@860
   183
    }
kpeter@860
   184
deba@409
   185
  private:
deba@409
   186
deba@409
   187
    void activate(const Node& i) {
kpeter@581
   188
      (*_active)[i] = true;
deba@409
   189
deba@409
   190
      int bucket = (*_bucket)[i];
deba@409
   191
deba@409
   192
      if ((*_prev)[i] == INVALID || (*_active)[(*_prev)[i]]) return;
deba@409
   193
      //unlace
kpeter@581
   194
      (*_next)[(*_prev)[i]] = (*_next)[i];
deba@409
   195
      if ((*_next)[i] != INVALID) {
kpeter@581
   196
        (*_prev)[(*_next)[i]] = (*_prev)[i];
deba@409
   197
      } else {
deba@409
   198
        _last[bucket] = (*_prev)[i];
deba@409
   199
      }
deba@409
   200
      //lace
kpeter@581
   201
      (*_next)[i] = _first[bucket];
kpeter@581
   202
      (*_prev)[_first[bucket]] = i;
kpeter@581
   203
      (*_prev)[i] = INVALID;
deba@409
   204
      _first[bucket] = i;
deba@409
   205
    }
deba@409
   206
deba@409
   207
    void deactivate(const Node& i) {
kpeter@581
   208
      (*_active)[i] = false;
deba@409
   209
      int bucket = (*_bucket)[i];
deba@409
   210
deba@409
   211
      if ((*_next)[i] == INVALID || !(*_active)[(*_next)[i]]) return;
deba@409
   212
deba@409
   213
      //unlace
kpeter@581
   214
      (*_prev)[(*_next)[i]] = (*_prev)[i];
deba@409
   215
      if ((*_prev)[i] != INVALID) {
kpeter@581
   216
        (*_next)[(*_prev)[i]] = (*_next)[i];
deba@409
   217
      } else {
deba@409
   218
        _first[bucket] = (*_next)[i];
deba@409
   219
      }
deba@409
   220
      //lace
kpeter@581
   221
      (*_prev)[i] = _last[bucket];
kpeter@581
   222
      (*_next)[_last[bucket]] = i;
kpeter@581
   223
      (*_next)[i] = INVALID;
deba@409
   224
      _last[bucket] = i;
deba@409
   225
    }
deba@409
   226
deba@409
   227
    void addItem(const Node& i, int bucket) {
deba@409
   228
      (*_bucket)[i] = bucket;
deba@409
   229
      if (_last[bucket] != INVALID) {
kpeter@581
   230
        (*_prev)[i] = _last[bucket];
kpeter@581
   231
        (*_next)[_last[bucket]] = i;
kpeter@581
   232
        (*_next)[i] = INVALID;
deba@409
   233
        _last[bucket] = i;
deba@409
   234
      } else {
kpeter@581
   235
        (*_prev)[i] = INVALID;
deba@409
   236
        _first[bucket] = i;
kpeter@581
   237
        (*_next)[i] = INVALID;
deba@409
   238
        _last[bucket] = i;
deba@409
   239
      }
deba@409
   240
    }
deba@409
   241
deba@409
   242
    void findMinCutOut() {
deba@409
   243
deba@409
   244
      for (NodeIt n(_graph); n != INVALID; ++n) {
kpeter@581
   245
        (*_excess)[n] = 0;
deba@597
   246
        (*_source_set)[n] = false;
deba@409
   247
      }
deba@409
   248
deba@409
   249
      for (ArcIt a(_graph); a != INVALID; ++a) {
kpeter@581
   250
        (*_flow)[a] = 0;
deba@409
   251
      }
deba@409
   252
deba@411
   253
      int bucket_num = 0;
deba@411
   254
      std::vector<Node> queue(_node_num);
deba@411
   255
      int qfirst = 0, qlast = 0, qsep = 0;
deba@409
   256
deba@409
   257
      {
deba@409
   258
        typename Digraph::template NodeMap<bool> reached(_graph, false);
deba@409
   259
kpeter@581
   260
        reached[_source] = true;
deba@409
   261
        bool first_set = true;
deba@409
   262
deba@409
   263
        for (NodeIt t(_graph); t != INVALID; ++t) {
deba@409
   264
          if (reached[t]) continue;
deba@409
   265
          _sets.push_front(std::list<int>());
alpar@440
   266
deba@411
   267
          queue[qlast++] = t;
kpeter@581
   268
          reached[t] = true;
deba@409
   269
deba@411
   270
          while (qfirst != qlast) {
deba@411
   271
            if (qsep == qfirst) {
deba@411
   272
              ++bucket_num;
deba@411
   273
              _sets.front().push_front(bucket_num);
deba@411
   274
              _dormant[bucket_num] = !first_set;
deba@411
   275
              _first[bucket_num] = _last[bucket_num] = INVALID;
deba@411
   276
              qsep = qlast;
deba@411
   277
            }
deba@409
   278
deba@411
   279
            Node n = queue[qfirst++];
deba@411
   280
            addItem(n, bucket_num);
deba@411
   281
deba@411
   282
            for (InArcIt a(_graph, n); a != INVALID; ++a) {
deba@411
   283
              Node u = _graph.source(a);
deba@411
   284
              if (!reached[u] && _tolerance.positive((*_capacity)[a])) {
kpeter@581
   285
                reached[u] = true;
deba@411
   286
                queue[qlast++] = u;
deba@409
   287
              }
deba@409
   288
            }
deba@409
   289
          }
deba@409
   290
          first_set = false;
deba@409
   291
        }
deba@409
   292
deba@411
   293
        ++bucket_num;
kpeter@581
   294
        (*_bucket)[_source] = 0;
deba@409
   295
        _dormant[0] = true;
deba@409
   296
      }
kpeter@581
   297
      (*_source_set)[_source] = true;
deba@409
   298
deba@409
   299
      Node target = _last[_sets.back().back()];
deba@409
   300
      {
deba@409
   301
        for (OutArcIt a(_graph, _source); a != INVALID; ++a) {
deba@409
   302
          if (_tolerance.positive((*_capacity)[a])) {
deba@409
   303
            Node u = _graph.target(a);
kpeter@581
   304
            (*_flow)[a] = (*_capacity)[a];
kpeter@581
   305
            (*_excess)[u] += (*_capacity)[a];
deba@409
   306
            if (!(*_active)[u] && u != _source) {
deba@409
   307
              activate(u);
deba@409
   308
            }
deba@409
   309
          }
deba@409
   310
        }
deba@409
   311
deba@409
   312
        if ((*_active)[target]) {
deba@409
   313
          deactivate(target);
deba@409
   314
        }
deba@409
   315
deba@409
   316
        _highest = _sets.back().begin();
deba@409
   317
        while (_highest != _sets.back().end() &&
deba@409
   318
               !(*_active)[_first[*_highest]]) {
deba@409
   319
          ++_highest;
deba@409
   320
        }
deba@409
   321
      }
deba@409
   322
deba@409
   323
      while (true) {
deba@409
   324
        while (_highest != _sets.back().end()) {
deba@409
   325
          Node n = _first[*_highest];
deba@409
   326
          Value excess = (*_excess)[n];
deba@409
   327
          int next_bucket = _node_num;
deba@409
   328
deba@409
   329
          int under_bucket;
deba@409
   330
          if (++std::list<int>::iterator(_highest) == _sets.back().end()) {
deba@409
   331
            under_bucket = -1;
deba@409
   332
          } else {
deba@409
   333
            under_bucket = *(++std::list<int>::iterator(_highest));
deba@409
   334
          }
deba@409
   335
deba@409
   336
          for (OutArcIt a(_graph, n); a != INVALID; ++a) {
deba@409
   337
            Node v = _graph.target(a);
deba@409
   338
            if (_dormant[(*_bucket)[v]]) continue;
deba@409
   339
            Value rem = (*_capacity)[a] - (*_flow)[a];
deba@409
   340
            if (!_tolerance.positive(rem)) continue;
deba@409
   341
            if ((*_bucket)[v] == under_bucket) {
deba@409
   342
              if (!(*_active)[v] && v != target) {
deba@409
   343
                activate(v);
deba@409
   344
              }
deba@409
   345
              if (!_tolerance.less(rem, excess)) {
kpeter@581
   346
                (*_flow)[a] += excess;
kpeter@581
   347
                (*_excess)[v] += excess;
deba@409
   348
                excess = 0;
deba@409
   349
                goto no_more_push;
deba@409
   350
              } else {
deba@409
   351
                excess -= rem;
kpeter@581
   352
                (*_excess)[v] += rem;
kpeter@581
   353
                (*_flow)[a] = (*_capacity)[a];
deba@409
   354
              }
deba@409
   355
            } else if (next_bucket > (*_bucket)[v]) {
deba@409
   356
              next_bucket = (*_bucket)[v];
deba@409
   357
            }
deba@409
   358
          }
deba@409
   359
deba@409
   360
          for (InArcIt a(_graph, n); a != INVALID; ++a) {
deba@409
   361
            Node v = _graph.source(a);
deba@409
   362
            if (_dormant[(*_bucket)[v]]) continue;
deba@409
   363
            Value rem = (*_flow)[a];
deba@409
   364
            if (!_tolerance.positive(rem)) continue;
deba@409
   365
            if ((*_bucket)[v] == under_bucket) {
deba@409
   366
              if (!(*_active)[v] && v != target) {
deba@409
   367
                activate(v);
deba@409
   368
              }
deba@409
   369
              if (!_tolerance.less(rem, excess)) {
kpeter@581
   370
                (*_flow)[a] -= excess;
kpeter@581
   371
                (*_excess)[v] += excess;
deba@409
   372
                excess = 0;
deba@409
   373
                goto no_more_push;
deba@409
   374
              } else {
deba@409
   375
                excess -= rem;
kpeter@581
   376
                (*_excess)[v] += rem;
kpeter@581
   377
                (*_flow)[a] = 0;
deba@409
   378
              }
deba@409
   379
            } else if (next_bucket > (*_bucket)[v]) {
deba@409
   380
              next_bucket = (*_bucket)[v];
deba@409
   381
            }
deba@409
   382
          }
deba@409
   383
deba@409
   384
        no_more_push:
deba@409
   385
kpeter@581
   386
          (*_excess)[n] = excess;
deba@409
   387
deba@409
   388
          if (excess != 0) {
deba@409
   389
            if ((*_next)[n] == INVALID) {
deba@409
   390
              typename std::list<std::list<int> >::iterator new_set =
deba@409
   391
                _sets.insert(--_sets.end(), std::list<int>());
deba@409
   392
              new_set->splice(new_set->end(), _sets.back(),
deba@409
   393
                              _sets.back().begin(), ++_highest);
deba@409
   394
              for (std::list<int>::iterator it = new_set->begin();
deba@409
   395
                   it != new_set->end(); ++it) {
deba@409
   396
                _dormant[*it] = true;
deba@409
   397
              }
deba@409
   398
              while (_highest != _sets.back().end() &&
deba@409
   399
                     !(*_active)[_first[*_highest]]) {
deba@409
   400
                ++_highest;
deba@409
   401
              }
deba@409
   402
            } else if (next_bucket == _node_num) {
deba@409
   403
              _first[(*_bucket)[n]] = (*_next)[n];
kpeter@581
   404
              (*_prev)[(*_next)[n]] = INVALID;
deba@409
   405
deba@409
   406
              std::list<std::list<int> >::iterator new_set =
deba@409
   407
                _sets.insert(--_sets.end(), std::list<int>());
deba@409
   408
deba@409
   409
              new_set->push_front(bucket_num);
kpeter@581
   410
              (*_bucket)[n] = bucket_num;
deba@409
   411
              _first[bucket_num] = _last[bucket_num] = n;
kpeter@581
   412
              (*_next)[n] = INVALID;
kpeter@581
   413
              (*_prev)[n] = INVALID;
deba@409
   414
              _dormant[bucket_num] = true;
deba@409
   415
              ++bucket_num;
deba@409
   416
deba@409
   417
              while (_highest != _sets.back().end() &&
deba@409
   418
                     !(*_active)[_first[*_highest]]) {
deba@409
   419
                ++_highest;
deba@409
   420
              }
deba@409
   421
            } else {
deba@409
   422
              _first[*_highest] = (*_next)[n];
kpeter@581
   423
              (*_prev)[(*_next)[n]] = INVALID;
deba@409
   424
deba@409
   425
              while (next_bucket != *_highest) {
deba@409
   426
                --_highest;
deba@409
   427
              }
deba@409
   428
deba@409
   429
              if (_highest == _sets.back().begin()) {
deba@409
   430
                _sets.back().push_front(bucket_num);
deba@409
   431
                _dormant[bucket_num] = false;
deba@409
   432
                _first[bucket_num] = _last[bucket_num] = INVALID;
deba@409
   433
                ++bucket_num;
deba@409
   434
              }
deba@409
   435
              --_highest;
deba@409
   436
kpeter@581
   437
              (*_bucket)[n] = *_highest;
kpeter@581
   438
              (*_next)[n] = _first[*_highest];
deba@409
   439
              if (_first[*_highest] != INVALID) {
kpeter@581
   440
                (*_prev)[_first[*_highest]] = n;
deba@409
   441
              } else {
deba@409
   442
                _last[*_highest] = n;
deba@409
   443
              }
deba@409
   444
              _first[*_highest] = n;
deba@409
   445
            }
deba@409
   446
          } else {
deba@409
   447
deba@409
   448
            deactivate(n);
deba@409
   449
            if (!(*_active)[_first[*_highest]]) {
deba@409
   450
              ++_highest;
deba@409
   451
              if (_highest != _sets.back().end() &&
deba@409
   452
                  !(*_active)[_first[*_highest]]) {
deba@409
   453
                _highest = _sets.back().end();
deba@409
   454
              }
deba@409
   455
            }
deba@409
   456
          }
deba@409
   457
        }
deba@409
   458
deba@409
   459
        if ((*_excess)[target] < _min_cut) {
deba@409
   460
          _min_cut = (*_excess)[target];
deba@409
   461
          for (NodeIt i(_graph); i != INVALID; ++i) {
kpeter@581
   462
            (*_min_cut_map)[i] = true;
deba@409
   463
          }
deba@409
   464
          for (std::list<int>::iterator it = _sets.back().begin();
deba@409
   465
               it != _sets.back().end(); ++it) {
deba@409
   466
            Node n = _first[*it];
deba@409
   467
            while (n != INVALID) {
kpeter@581
   468
              (*_min_cut_map)[n] = false;
deba@409
   469
              n = (*_next)[n];
deba@409
   470
            }
deba@409
   471
          }
deba@409
   472
        }
deba@409
   473
deba@409
   474
        {
deba@409
   475
          Node new_target;
deba@409
   476
          if ((*_prev)[target] != INVALID || (*_next)[target] != INVALID) {
deba@409
   477
            if ((*_next)[target] == INVALID) {
deba@409
   478
              _last[(*_bucket)[target]] = (*_prev)[target];
deba@409
   479
              new_target = (*_prev)[target];
deba@409
   480
            } else {
kpeter@581
   481
              (*_prev)[(*_next)[target]] = (*_prev)[target];
deba@409
   482
              new_target = (*_next)[target];
deba@409
   483
            }
deba@409
   484
            if ((*_prev)[target] == INVALID) {
deba@409
   485
              _first[(*_bucket)[target]] = (*_next)[target];
deba@409
   486
            } else {
kpeter@581
   487
              (*_next)[(*_prev)[target]] = (*_next)[target];
deba@409
   488
            }
deba@409
   489
          } else {
deba@409
   490
            _sets.back().pop_back();
deba@409
   491
            if (_sets.back().empty()) {
deba@409
   492
              _sets.pop_back();
deba@409
   493
              if (_sets.empty())
deba@409
   494
                break;
deba@409
   495
              for (std::list<int>::iterator it = _sets.back().begin();
deba@409
   496
                   it != _sets.back().end(); ++it) {
deba@409
   497
                _dormant[*it] = false;
deba@409
   498
              }
deba@409
   499
            }
deba@409
   500
            new_target = _last[_sets.back().back()];
deba@409
   501
          }
deba@409
   502
kpeter@581
   503
          (*_bucket)[target] = 0;
deba@409
   504
kpeter@581
   505
          (*_source_set)[target] = true;
deba@409
   506
          for (OutArcIt a(_graph, target); a != INVALID; ++a) {
deba@409
   507
            Value rem = (*_capacity)[a] - (*_flow)[a];
deba@409
   508
            if (!_tolerance.positive(rem)) continue;
deba@409
   509
            Node v = _graph.target(a);
deba@409
   510
            if (!(*_active)[v] && !(*_source_set)[v]) {
deba@409
   511
              activate(v);
deba@409
   512
            }
kpeter@581
   513
            (*_excess)[v] += rem;
kpeter@581
   514
            (*_flow)[a] = (*_capacity)[a];
deba@409
   515
          }
deba@409
   516
deba@409
   517
          for (InArcIt a(_graph, target); a != INVALID; ++a) {
deba@409
   518
            Value rem = (*_flow)[a];
deba@409
   519
            if (!_tolerance.positive(rem)) continue;
deba@409
   520
            Node v = _graph.source(a);
deba@409
   521
            if (!(*_active)[v] && !(*_source_set)[v]) {
deba@409
   522
              activate(v);
deba@409
   523
            }
kpeter@581
   524
            (*_excess)[v] += rem;
kpeter@581
   525
            (*_flow)[a] = 0;
deba@409
   526
          }
deba@409
   527
deba@409
   528
          target = new_target;
deba@409
   529
          if ((*_active)[target]) {
deba@409
   530
            deactivate(target);
deba@409
   531
          }
deba@409
   532
deba@409
   533
          _highest = _sets.back().begin();
deba@409
   534
          while (_highest != _sets.back().end() &&
deba@409
   535
                 !(*_active)[_first[*_highest]]) {
deba@409
   536
            ++_highest;
deba@409
   537
          }
deba@409
   538
        }
deba@409
   539
      }
deba@409
   540
    }
deba@409
   541
deba@409
   542
    void findMinCutIn() {
deba@409
   543
deba@409
   544
      for (NodeIt n(_graph); n != INVALID; ++n) {
kpeter@581
   545
        (*_excess)[n] = 0;
deba@597
   546
        (*_source_set)[n] = false;
deba@409
   547
      }
deba@409
   548
deba@409
   549
      for (ArcIt a(_graph); a != INVALID; ++a) {
kpeter@581
   550
        (*_flow)[a] = 0;
deba@409
   551
      }
deba@409
   552
deba@411
   553
      int bucket_num = 0;
deba@411
   554
      std::vector<Node> queue(_node_num);
deba@411
   555
      int qfirst = 0, qlast = 0, qsep = 0;
deba@409
   556
deba@409
   557
      {
deba@409
   558
        typename Digraph::template NodeMap<bool> reached(_graph, false);
deba@409
   559
kpeter@581
   560
        reached[_source] = true;
deba@409
   561
deba@409
   562
        bool first_set = true;
deba@409
   563
deba@409
   564
        for (NodeIt t(_graph); t != INVALID; ++t) {
deba@409
   565
          if (reached[t]) continue;
deba@409
   566
          _sets.push_front(std::list<int>());
alpar@440
   567
deba@411
   568
          queue[qlast++] = t;
kpeter@581
   569
          reached[t] = true;
deba@409
   570
deba@411
   571
          while (qfirst != qlast) {
deba@411
   572
            if (qsep == qfirst) {
deba@411
   573
              ++bucket_num;
deba@411
   574
              _sets.front().push_front(bucket_num);
deba@411
   575
              _dormant[bucket_num] = !first_set;
deba@411
   576
              _first[bucket_num] = _last[bucket_num] = INVALID;
deba@411
   577
              qsep = qlast;
deba@411
   578
            }
deba@409
   579
deba@411
   580
            Node n = queue[qfirst++];
deba@411
   581
            addItem(n, bucket_num);
deba@411
   582
deba@411
   583
            for (OutArcIt a(_graph, n); a != INVALID; ++a) {
deba@411
   584
              Node u = _graph.target(a);
deba@411
   585
              if (!reached[u] && _tolerance.positive((*_capacity)[a])) {
kpeter@581
   586
                reached[u] = true;
deba@411
   587
                queue[qlast++] = u;
deba@409
   588
              }
deba@409
   589
            }
deba@409
   590
          }
deba@409
   591
          first_set = false;
deba@409
   592
        }
deba@409
   593
deba@411
   594
        ++bucket_num;
kpeter@581
   595
        (*_bucket)[_source] = 0;
deba@409
   596
        _dormant[0] = true;
deba@409
   597
      }
kpeter@581
   598
      (*_source_set)[_source] = true;
deba@409
   599
deba@409
   600
      Node target = _last[_sets.back().back()];
deba@409
   601
      {
deba@409
   602
        for (InArcIt a(_graph, _source); a != INVALID; ++a) {
deba@409
   603
          if (_tolerance.positive((*_capacity)[a])) {
deba@409
   604
            Node u = _graph.source(a);
kpeter@581
   605
            (*_flow)[a] = (*_capacity)[a];
kpeter@581
   606
            (*_excess)[u] += (*_capacity)[a];
deba@409
   607
            if (!(*_active)[u] && u != _source) {
deba@409
   608
              activate(u);
deba@409
   609
            }
deba@409
   610
          }
deba@409
   611
        }
deba@409
   612
        if ((*_active)[target]) {
deba@409
   613
          deactivate(target);
deba@409
   614
        }
deba@409
   615
deba@409
   616
        _highest = _sets.back().begin();
deba@409
   617
        while (_highest != _sets.back().end() &&
deba@409
   618
               !(*_active)[_first[*_highest]]) {
deba@409
   619
          ++_highest;
deba@409
   620
        }
deba@409
   621
      }
deba@409
   622
deba@409
   623
deba@409
   624
      while (true) {
deba@409
   625
        while (_highest != _sets.back().end()) {
deba@409
   626
          Node n = _first[*_highest];
deba@409
   627
          Value excess = (*_excess)[n];
deba@409
   628
          int next_bucket = _node_num;
deba@409
   629
deba@409
   630
          int under_bucket;
deba@409
   631
          if (++std::list<int>::iterator(_highest) == _sets.back().end()) {
deba@409
   632
            under_bucket = -1;
deba@409
   633
          } else {
deba@409
   634
            under_bucket = *(++std::list<int>::iterator(_highest));
deba@409
   635
          }
deba@409
   636
deba@409
   637
          for (InArcIt a(_graph, n); a != INVALID; ++a) {
deba@409
   638
            Node v = _graph.source(a);
deba@409
   639
            if (_dormant[(*_bucket)[v]]) continue;
deba@409
   640
            Value rem = (*_capacity)[a] - (*_flow)[a];
deba@409
   641
            if (!_tolerance.positive(rem)) continue;
deba@409
   642
            if ((*_bucket)[v] == under_bucket) {
deba@409
   643
              if (!(*_active)[v] && v != target) {
deba@409
   644
                activate(v);
deba@409
   645
              }
deba@409
   646
              if (!_tolerance.less(rem, excess)) {
kpeter@581
   647
                (*_flow)[a] += excess;
kpeter@581
   648
                (*_excess)[v] += excess;
deba@409
   649
                excess = 0;
deba@409
   650
                goto no_more_push;
deba@409
   651
              } else {
deba@409
   652
                excess -= rem;
kpeter@581
   653
                (*_excess)[v] += rem;
kpeter@581
   654
                (*_flow)[a] = (*_capacity)[a];
deba@409
   655
              }
deba@409
   656
            } else if (next_bucket > (*_bucket)[v]) {
deba@409
   657
              next_bucket = (*_bucket)[v];
deba@409
   658
            }
deba@409
   659
          }
deba@409
   660
deba@409
   661
          for (OutArcIt a(_graph, n); a != INVALID; ++a) {
deba@409
   662
            Node v = _graph.target(a);
deba@409
   663
            if (_dormant[(*_bucket)[v]]) continue;
deba@409
   664
            Value rem = (*_flow)[a];
deba@409
   665
            if (!_tolerance.positive(rem)) continue;
deba@409
   666
            if ((*_bucket)[v] == under_bucket) {
deba@409
   667
              if (!(*_active)[v] && v != target) {
deba@409
   668
                activate(v);
deba@409
   669
              }
deba@409
   670
              if (!_tolerance.less(rem, excess)) {
kpeter@581
   671
                (*_flow)[a] -= excess;
kpeter@581
   672
                (*_excess)[v] += excess;
deba@409
   673
                excess = 0;
deba@409
   674
                goto no_more_push;
deba@409
   675
              } else {
deba@409
   676
                excess -= rem;
kpeter@581
   677
                (*_excess)[v] += rem;
kpeter@581
   678
                (*_flow)[a] = 0;
deba@409
   679
              }
deba@409
   680
            } else if (next_bucket > (*_bucket)[v]) {
deba@409
   681
              next_bucket = (*_bucket)[v];
deba@409
   682
            }
deba@409
   683
          }
deba@409
   684
deba@409
   685
        no_more_push:
deba@409
   686
kpeter@581
   687
          (*_excess)[n] = excess;
deba@409
   688
deba@409
   689
          if (excess != 0) {
deba@409
   690
            if ((*_next)[n] == INVALID) {
deba@409
   691
              typename std::list<std::list<int> >::iterator new_set =
deba@409
   692
                _sets.insert(--_sets.end(), std::list<int>());
deba@409
   693
              new_set->splice(new_set->end(), _sets.back(),
deba@409
   694
                              _sets.back().begin(), ++_highest);
deba@409
   695
              for (std::list<int>::iterator it = new_set->begin();
deba@409
   696
                   it != new_set->end(); ++it) {
deba@409
   697
                _dormant[*it] = true;
deba@409
   698
              }
deba@409
   699
              while (_highest != _sets.back().end() &&
deba@409
   700
                     !(*_active)[_first[*_highest]]) {
deba@409
   701
                ++_highest;
deba@409
   702
              }
deba@409
   703
            } else if (next_bucket == _node_num) {
deba@409
   704
              _first[(*_bucket)[n]] = (*_next)[n];
kpeter@581
   705
              (*_prev)[(*_next)[n]] = INVALID;
deba@409
   706
deba@409
   707
              std::list<std::list<int> >::iterator new_set =
deba@409
   708
                _sets.insert(--_sets.end(), std::list<int>());
deba@409
   709
deba@409
   710
              new_set->push_front(bucket_num);
kpeter@581
   711
              (*_bucket)[n] = bucket_num;
deba@409
   712
              _first[bucket_num] = _last[bucket_num] = n;
kpeter@581
   713
              (*_next)[n] = INVALID;
kpeter@581
   714
              (*_prev)[n] = INVALID;
deba@409
   715
              _dormant[bucket_num] = true;
deba@409
   716
              ++bucket_num;
deba@409
   717
deba@409
   718
              while (_highest != _sets.back().end() &&
deba@409
   719
                     !(*_active)[_first[*_highest]]) {
deba@409
   720
                ++_highest;
deba@409
   721
              }
deba@409
   722
            } else {
deba@409
   723
              _first[*_highest] = (*_next)[n];
kpeter@581
   724
              (*_prev)[(*_next)[n]] = INVALID;
deba@409
   725
deba@409
   726
              while (next_bucket != *_highest) {
deba@409
   727
                --_highest;
deba@409
   728
              }
deba@409
   729
              if (_highest == _sets.back().begin()) {
deba@409
   730
                _sets.back().push_front(bucket_num);
deba@409
   731
                _dormant[bucket_num] = false;
deba@409
   732
                _first[bucket_num] = _last[bucket_num] = INVALID;
deba@409
   733
                ++bucket_num;
deba@409
   734
              }
deba@409
   735
              --_highest;
deba@409
   736
kpeter@581
   737
              (*_bucket)[n] = *_highest;
kpeter@581
   738
              (*_next)[n] = _first[*_highest];
deba@409
   739
              if (_first[*_highest] != INVALID) {
kpeter@581
   740
                (*_prev)[_first[*_highest]] = n;
deba@409
   741
              } else {
deba@409
   742
                _last[*_highest] = n;
deba@409
   743
              }
deba@409
   744
              _first[*_highest] = n;
deba@409
   745
            }
deba@409
   746
          } else {
deba@409
   747
deba@409
   748
            deactivate(n);
deba@409
   749
            if (!(*_active)[_first[*_highest]]) {
deba@409
   750
              ++_highest;
deba@409
   751
              if (_highest != _sets.back().end() &&
deba@409
   752
                  !(*_active)[_first[*_highest]]) {
deba@409
   753
                _highest = _sets.back().end();
deba@409
   754
              }
deba@409
   755
            }
deba@409
   756
          }
deba@409
   757
        }
deba@409
   758
deba@409
   759
        if ((*_excess)[target] < _min_cut) {
deba@409
   760
          _min_cut = (*_excess)[target];
deba@409
   761
          for (NodeIt i(_graph); i != INVALID; ++i) {
kpeter@581
   762
            (*_min_cut_map)[i] = false;
deba@409
   763
          }
deba@409
   764
          for (std::list<int>::iterator it = _sets.back().begin();
deba@409
   765
               it != _sets.back().end(); ++it) {
deba@409
   766
            Node n = _first[*it];
deba@409
   767
            while (n != INVALID) {
kpeter@581
   768
              (*_min_cut_map)[n] = true;
deba@409
   769
              n = (*_next)[n];
deba@409
   770
            }
deba@409
   771
          }
deba@409
   772
        }
deba@409
   773
deba@409
   774
        {
deba@409
   775
          Node new_target;
deba@409
   776
          if ((*_prev)[target] != INVALID || (*_next)[target] != INVALID) {
deba@409
   777
            if ((*_next)[target] == INVALID) {
deba@409
   778
              _last[(*_bucket)[target]] = (*_prev)[target];
deba@409
   779
              new_target = (*_prev)[target];
deba@409
   780
            } else {
kpeter@581
   781
              (*_prev)[(*_next)[target]] = (*_prev)[target];
deba@409
   782
              new_target = (*_next)[target];
deba@409
   783
            }
deba@409
   784
            if ((*_prev)[target] == INVALID) {
deba@409
   785
              _first[(*_bucket)[target]] = (*_next)[target];
deba@409
   786
            } else {
kpeter@581
   787
              (*_next)[(*_prev)[target]] = (*_next)[target];
deba@409
   788
            }
deba@409
   789
          } else {
deba@409
   790
            _sets.back().pop_back();
deba@409
   791
            if (_sets.back().empty()) {
deba@409
   792
              _sets.pop_back();
deba@409
   793
              if (_sets.empty())
deba@409
   794
                break;
deba@409
   795
              for (std::list<int>::iterator it = _sets.back().begin();
deba@409
   796
                   it != _sets.back().end(); ++it) {
deba@409
   797
                _dormant[*it] = false;
deba@409
   798
              }
deba@409
   799
            }
deba@409
   800
            new_target = _last[_sets.back().back()];
deba@409
   801
          }
deba@409
   802
kpeter@581
   803
          (*_bucket)[target] = 0;
deba@409
   804
kpeter@581
   805
          (*_source_set)[target] = true;
deba@409
   806
          for (InArcIt a(_graph, target); a != INVALID; ++a) {
deba@409
   807
            Value rem = (*_capacity)[a] - (*_flow)[a];
deba@409
   808
            if (!_tolerance.positive(rem)) continue;
deba@409
   809
            Node v = _graph.source(a);
deba@409
   810
            if (!(*_active)[v] && !(*_source_set)[v]) {
deba@409
   811
              activate(v);
deba@409
   812
            }
kpeter@581
   813
            (*_excess)[v] += rem;
kpeter@581
   814
            (*_flow)[a] = (*_capacity)[a];
deba@409
   815
          }
deba@409
   816
deba@409
   817
          for (OutArcIt a(_graph, target); a != INVALID; ++a) {
deba@409
   818
            Value rem = (*_flow)[a];
deba@409
   819
            if (!_tolerance.positive(rem)) continue;
deba@409
   820
            Node v = _graph.target(a);
deba@409
   821
            if (!(*_active)[v] && !(*_source_set)[v]) {
deba@409
   822
              activate(v);
deba@409
   823
            }
kpeter@581
   824
            (*_excess)[v] += rem;
kpeter@581
   825
            (*_flow)[a] = 0;
deba@409
   826
          }
deba@409
   827
deba@409
   828
          target = new_target;
deba@409
   829
          if ((*_active)[target]) {
deba@409
   830
            deactivate(target);
deba@409
   831
          }
deba@409
   832
deba@409
   833
          _highest = _sets.back().begin();
deba@409
   834
          while (_highest != _sets.back().end() &&
deba@409
   835
                 !(*_active)[_first[*_highest]]) {
deba@409
   836
            ++_highest;
deba@409
   837
          }
deba@409
   838
        }
deba@409
   839
      }
deba@409
   840
    }
deba@409
   841
deba@409
   842
  public:
deba@409
   843
kpeter@596
   844
    /// \name Execution Control
deba@409
   845
    /// The simplest way to execute the algorithm is to use
kpeter@559
   846
    /// one of the member functions called \ref run().
deba@409
   847
    /// \n
kpeter@596
   848
    /// If you need better control on the execution,
kpeter@596
   849
    /// you have to call one of the \ref init() functions first, then
kpeter@596
   850
    /// \ref calculateOut() and/or \ref calculateIn().
deba@409
   851
deba@409
   852
    /// @{
deba@409
   853
kpeter@596
   854
    /// \brief Initialize the internal data structures.
deba@409
   855
    ///
kpeter@596
   856
    /// This function initializes the internal data structures. It creates
kpeter@596
   857
    /// the maps and some bucket structures for the algorithm.
kpeter@596
   858
    /// The first node is used as the source node for the push-relabel
kpeter@596
   859
    /// algorithm.
deba@409
   860
    void init() {
deba@409
   861
      init(NodeIt(_graph));
deba@409
   862
    }
deba@409
   863
kpeter@596
   864
    /// \brief Initialize the internal data structures.
deba@409
   865
    ///
kpeter@596
   866
    /// This function initializes the internal data structures. It creates
alpar@877
   867
    /// the maps and some bucket structures for the algorithm.
kpeter@596
   868
    /// The given node is used as the source node for the push-relabel
kpeter@596
   869
    /// algorithm.
deba@409
   870
    void init(const Node& source) {
deba@409
   871
      _source = source;
deba@409
   872
deba@409
   873
      _node_num = countNodes(_graph);
deba@409
   874
deba@411
   875
      _first.resize(_node_num);
deba@411
   876
      _last.resize(_node_num);
deba@409
   877
deba@411
   878
      _dormant.resize(_node_num);
deba@409
   879
deba@409
   880
      if (!_flow) {
deba@409
   881
        _flow = new FlowMap(_graph);
deba@409
   882
      }
deba@409
   883
      if (!_next) {
deba@409
   884
        _next = new typename Digraph::template NodeMap<Node>(_graph);
deba@409
   885
      }
deba@409
   886
      if (!_prev) {
deba@409
   887
        _prev = new typename Digraph::template NodeMap<Node>(_graph);
deba@409
   888
      }
deba@409
   889
      if (!_active) {
deba@409
   890
        _active = new typename Digraph::template NodeMap<bool>(_graph);
deba@409
   891
      }
deba@409
   892
      if (!_bucket) {
deba@409
   893
        _bucket = new typename Digraph::template NodeMap<int>(_graph);
deba@409
   894
      }
deba@409
   895
      if (!_excess) {
deba@409
   896
        _excess = new ExcessMap(_graph);
deba@409
   897
      }
deba@409
   898
      if (!_source_set) {
deba@409
   899
        _source_set = new SourceSetMap(_graph);
deba@409
   900
      }
deba@409
   901
      if (!_min_cut_map) {
deba@409
   902
        _min_cut_map = new MinCutMap(_graph);
deba@409
   903
      }
deba@409
   904
deba@409
   905
      _min_cut = std::numeric_limits<Value>::max();
deba@409
   906
    }
deba@409
   907
deba@409
   908
kpeter@596
   909
    /// \brief Calculate a minimum cut with \f$ source \f$ on the
deba@409
   910
    /// source-side.
deba@409
   911
    ///
kpeter@596
   912
    /// This function calculates a minimum cut with \f$ source \f$ on the
alpar@412
   913
    /// source-side (i.e. a set \f$ X\subsetneq V \f$ with
kpeter@596
   914
    /// \f$ source \in X \f$ and minimal outgoing capacity).
alpar@915
   915
    /// It updates the stored cut if (and only if) the newly found one
alpar@915
   916
    /// is better.
kpeter@596
   917
    ///
kpeter@596
   918
    /// \pre \ref init() must be called before using this function.
deba@409
   919
    void calculateOut() {
deba@409
   920
      findMinCutOut();
deba@409
   921
    }
deba@409
   922
kpeter@596
   923
    /// \brief Calculate a minimum cut with \f$ source \f$ on the
kpeter@596
   924
    /// sink-side.
deba@409
   925
    ///
kpeter@596
   926
    /// This function calculates a minimum cut with \f$ source \f$ on the
kpeter@596
   927
    /// sink-side (i.e. a set \f$ X\subsetneq V \f$ with
kpeter@596
   928
    /// \f$ source \notin X \f$ and minimal outgoing capacity).
alpar@915
   929
    /// It updates the stored cut if (and only if) the newly found one
alpar@915
   930
    /// is better.
kpeter@596
   931
    ///
kpeter@596
   932
    /// \pre \ref init() must be called before using this function.
deba@409
   933
    void calculateIn() {
deba@409
   934
      findMinCutIn();
deba@409
   935
    }
deba@409
   936
deba@409
   937
kpeter@596
   938
    /// \brief Run the algorithm.
deba@409
   939
    ///
alpar@915
   940
    /// This function runs the algorithm. It chooses source node,
alpar@915
   941
    /// then calls \ref init(), \ref calculateOut()
deba@409
   942
    /// and \ref calculateIn().
deba@409
   943
    void run() {
deba@409
   944
      init();
deba@409
   945
      calculateOut();
deba@409
   946
      calculateIn();
deba@409
   947
    }
deba@409
   948
kpeter@596
   949
    /// \brief Run the algorithm.
deba@409
   950
    ///
alpar@915
   951
    /// This function runs the algorithm. It calls \ref init(),
alpar@915
   952
    /// \ref calculateOut() and \ref calculateIn() with the given
alpar@915
   953
    /// source node.
deba@409
   954
    void run(const Node& s) {
deba@409
   955
      init(s);
deba@409
   956
      calculateOut();
deba@409
   957
      calculateIn();
deba@409
   958
    }
deba@409
   959
deba@409
   960
    /// @}
deba@409
   961
deba@409
   962
    /// \name Query Functions
deba@409
   963
    /// The result of the %HaoOrlin algorithm
kpeter@596
   964
    /// can be obtained using these functions.\n
alpar@877
   965
    /// \ref run(), \ref calculateOut() or \ref calculateIn()
kpeter@596
   966
    /// should be called before using them.
deba@409
   967
deba@409
   968
    /// @{
deba@409
   969
kpeter@596
   970
    /// \brief Return the value of the minimum cut.
deba@409
   971
    ///
alpar@915
   972
    /// This function returns the value of the best cut found by the
alpar@915
   973
    /// previously called \ref run(), \ref calculateOut() or \ref
alpar@915
   974
    /// calculateIn().
kpeter@596
   975
    ///
alpar@877
   976
    /// \pre \ref run(), \ref calculateOut() or \ref calculateIn()
kpeter@596
   977
    /// must be called before using this function.
deba@409
   978
    Value minCutValue() const {
deba@409
   979
      return _min_cut;
deba@409
   980
    }
deba@409
   981
deba@409
   982
kpeter@596
   983
    /// \brief Return a minimum cut.
deba@409
   984
    ///
alpar@915
   985
    /// This function gives the best cut found by the
alpar@915
   986
    /// previously called \ref run(), \ref calculateOut() or \ref
alpar@915
   987
    /// calculateIn().
alpar@915
   988
    ///
alpar@915
   989
    /// It sets \c cutMap to the characteristic vector of the found
alpar@915
   990
    /// minimum value cut - a non-empty set \f$ X\subsetneq V \f$
alpar@915
   991
    /// of minimum outgoing capacity (i.e. \c cutMap will be \c true exactly
kpeter@596
   992
    /// for the nodes of \f$ X \f$).
kpeter@596
   993
    ///
kpeter@596
   994
    /// \param cutMap A \ref concepts::WriteMap "writable" node map with
kpeter@596
   995
    /// \c bool (or convertible) value type.
kpeter@596
   996
    ///
kpeter@596
   997
    /// \return The value of the minimum cut.
kpeter@596
   998
    ///
alpar@877
   999
    /// \pre \ref run(), \ref calculateOut() or \ref calculateIn()
kpeter@596
  1000
    /// must be called before using this function.
kpeter@596
  1001
    template <typename CutMap>
kpeter@596
  1002
    Value minCutMap(CutMap& cutMap) const {
deba@409
  1003
      for (NodeIt it(_graph); it != INVALID; ++it) {
kpeter@596
  1004
        cutMap.set(it, (*_min_cut_map)[it]);
deba@409
  1005
      }
deba@409
  1006
      return _min_cut;
deba@409
  1007
    }
deba@409
  1008
deba@409
  1009
    /// @}
deba@409
  1010
deba@409
  1011
  }; //class HaoOrlin
deba@409
  1012
deba@409
  1013
} //namespace lemon
deba@409
  1014
deba@409
  1015
#endif //LEMON_HAO_ORLIN_H