lemon/min_mean_cycle.h
author deba
Sat, 17 Nov 2007 20:58:11 +0000
changeset 2514 57143c09dc20
parent 2437 02c7076bf894
child 2517 d9cfac072869
permissions -rw-r--r--
Redesign the maximum flow algorithms

Redesigned interface
Preflow changed to use elevator
Edmonds-Karp does not use the ResGraphAdaptor
Goldberg-Tarjan algorithm (Preflow with Dynamic Trees)
Dinitz-Sleator-Tarjan (Blocking flow with Dynamic Tree)
alpar@2409
     1
/* -*- C++ -*-
alpar@2409
     2
 *
alpar@2409
     3
 * This file is a part of LEMON, a generic C++ optimization library
alpar@2409
     4
 *
alpar@2409
     5
 * Copyright (C) 2003-2007
alpar@2409
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
alpar@2409
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
alpar@2409
     8
 *
alpar@2409
     9
 * Permission to use, modify and distribute this software is granted
alpar@2409
    10
 * provided that this copyright notice appears in all copies. For
alpar@2409
    11
 * precise terms see the accompanying LICENSE file.
alpar@2409
    12
 *
alpar@2409
    13
 * This software is provided "AS IS" with no warranty of any kind,
alpar@2409
    14
 * express or implied, and with no claim as to its suitability for any
alpar@2409
    15
 * purpose.
alpar@2409
    16
 *
alpar@2409
    17
 */
alpar@2409
    18
alpar@2409
    19
#ifndef LEMON_MIN_MEAN_CYCLE_H
alpar@2409
    20
#define LEMON_MIN_MEAN_CYCLE_H
alpar@2409
    21
alpar@2409
    22
/// \ingroup min_cost_flow
alpar@2409
    23
///
deba@2437
    24
/// \file
kpeter@2509
    25
/// \brief Karp's algorithm for finding a minimum mean (directed) cycle.
alpar@2409
    26
kpeter@2509
    27
#include <vector>
alpar@2409
    28
#include <lemon/graph_utils.h>
alpar@2409
    29
#include <lemon/topology.h>
alpar@2409
    30
#include <lemon/path.h>
alpar@2409
    31
alpar@2409
    32
namespace lemon {
alpar@2409
    33
alpar@2409
    34
  /// \addtogroup min_cost_flow
alpar@2409
    35
  /// @{
alpar@2409
    36
deba@2413
    37
  /// \brief Implementation of Karp's algorithm for finding a
deba@2413
    38
  /// minimum mean (directed) cycle.
alpar@2409
    39
  ///
deba@2413
    40
  /// The \ref lemon::MinMeanCycle "MinMeanCycle" implements Karp's
deba@2413
    41
  /// algorithm for finding a minimum mean (directed) cycle.
alpar@2409
    42
  ///
alpar@2409
    43
  /// \param Graph The directed graph type the algorithm runs on.
alpar@2409
    44
  /// \param LengthMap The type of the length (cost) map.
alpar@2409
    45
  ///
alpar@2409
    46
  /// \author Peter Kovacs
alpar@2409
    47
alpar@2409
    48
  template <typename Graph,
alpar@2409
    49
    typename LengthMap = typename Graph::template EdgeMap<int> >
alpar@2409
    50
  class MinMeanCycle
alpar@2409
    51
  {
alpar@2409
    52
    typedef typename Graph::Node Node;
alpar@2409
    53
    typedef typename Graph::NodeIt NodeIt;
alpar@2409
    54
    typedef typename Graph::Edge Edge;
alpar@2409
    55
    typedef typename Graph::EdgeIt EdgeIt;
alpar@2409
    56
    typedef typename Graph::OutEdgeIt OutEdgeIt;
alpar@2409
    57
alpar@2409
    58
    typedef typename LengthMap::Value Length;
deba@2413
    59
alpar@2409
    60
    typedef typename Graph::template NodeMap<int> IntNodeMap;
alpar@2409
    61
    typedef typename Graph::template NodeMap<Edge> PredNodeMap;
alpar@2409
    62
    typedef Path<Graph> Path;
alpar@2409
    63
    typedef std::vector<Node> NodeVector;
alpar@2409
    64
alpar@2409
    65
  protected:
deba@2413
    66
kpeter@2509
    67
    /// \brief Data structure for path data.
deba@2413
    68
    struct PathData
deba@2413
    69
    {
alpar@2409
    70
      bool found;
alpar@2409
    71
      Length dist;
alpar@2409
    72
      Edge pred;
deba@2413
    73
      PathData(bool _found = false, Length _dist = 0) :
alpar@2409
    74
	found(_found), dist(_dist), pred(INVALID) {}
deba@2413
    75
      PathData(bool _found, Length _dist, Edge _pred) :
alpar@2409
    76
	found(_found), dist(_dist), pred(_pred) {}
alpar@2409
    77
    };
deba@2413
    78
alpar@2409
    79
  private:
deba@2413
    80
deba@2413
    81
    typedef typename Graph::template NodeMap<std::vector<PathData> >
deba@2413
    82
      PathDataNodeMap;
deba@2413
    83
alpar@2409
    84
  protected:
alpar@2409
    85
alpar@2409
    86
    /// \brief Node map for storing path data.
deba@2437
    87
    ///
alpar@2409
    88
    /// Node map for storing path data of all nodes in the current
deba@2437
    89
    /// component. dmap[v][k] is the length of a shortest directed walk
alpar@2409
    90
    /// to node v from the starting node containing exactly k edges.
deba@2413
    91
    PathDataNodeMap dmap;
deba@2437
    92
alpar@2409
    93
    /// \brief The directed graph the algorithm runs on.
deba@2437
    94
    const Graph &graph;
alpar@2409
    95
    /// \brief The length of the edges.
deba@2437
    96
    const LengthMap &length;
deba@2437
    97
alpar@2409
    98
    /// \brief The total length of the found cycle.
alpar@2409
    99
    Length cycle_length;
alpar@2409
   100
    /// \brief The number of edges in the found cycle.
alpar@2409
   101
    int cycle_size;
deba@2437
   102
    /// \brief A node for obtaining a minimum mean cycle.
deba@2413
   103
    Node cycle_node;
deba@2413
   104
alpar@2409
   105
    /// \brief The found cycle.
alpar@2409
   106
    Path *cycle_path;
deba@2437
   107
    /// \brief The algorithm uses local \ref lemon::Path "Path"
deba@2413
   108
    /// structure to store the found cycle.
alpar@2409
   109
    bool local_path;
deba@2437
   110
alpar@2409
   111
    /// \brief Node map for identifying strongly connected components.
alpar@2409
   112
    IntNodeMap comp;
alpar@2409
   113
    /// \brief The number of strongly connected components.
alpar@2409
   114
    int comp_num;
deba@2413
   115
    /// \brief Counter for identifying the current component.
alpar@2409
   116
    int comp_cnt;
alpar@2409
   117
    /// \brief Nodes of the current component.
alpar@2409
   118
    NodeVector nodes;
alpar@2409
   119
    /// \brief The processed nodes in the last round.
alpar@2409
   120
    NodeVector process;
deba@2437
   121
alpar@2409
   122
  public :
alpar@2409
   123
alpar@2409
   124
    /// \brief The constructor of the class.
alpar@2409
   125
    ///
alpar@2409
   126
    /// The constructor of the class.
alpar@2409
   127
    ///
alpar@2409
   128
    /// \param _graph The directed graph the algorithm runs on.
alpar@2409
   129
    /// \param _length The length (cost) of the edges.
deba@2437
   130
    MinMeanCycle( const Graph &_graph,
deba@2437
   131
		  const LengthMap &_length ) :
deba@2413
   132
      graph(_graph), length(_length), dmap(_graph), comp(_graph),
deba@2413
   133
      cycle_length(0), cycle_size(-1), cycle_node(INVALID),
alpar@2409
   134
      cycle_path(NULL), local_path(false)
alpar@2409
   135
    { }
alpar@2409
   136
alpar@2409
   137
    /// \brief The destructor of the class.
deba@2437
   138
    ~MinMeanCycle() {
alpar@2409
   139
      if (local_path) delete cycle_path;
alpar@2409
   140
    }
alpar@2409
   141
alpar@2409
   142
  protected:
alpar@2409
   143
alpar@2409
   144
    /// \brief Initializes the internal data structures for the current
alpar@2409
   145
    /// component.
alpar@2409
   146
    void initCurrent() {
alpar@2409
   147
      nodes.clear();
alpar@2409
   148
      // Finding the nodes of the current component
alpar@2409
   149
      for (NodeIt v(graph); v != INVALID; ++v) {
alpar@2409
   150
	if (comp[v] == comp_cnt) nodes.push_back(v);
alpar@2409
   151
      }
alpar@2409
   152
      // Creating vectors for all nodes
alpar@2409
   153
      int n = nodes.size();
kpeter@2509
   154
      for (int i = 0; i < n; ++i) {
deba@2437
   155
	dmap[nodes[i]].resize(n + 1);
alpar@2409
   156
      }
alpar@2409
   157
    }
deba@2437
   158
deba@2437
   159
    /// \brief Processes all rounds of computing required path data for
deba@2413
   160
    /// the current component.
alpar@2409
   161
    void processRounds() {
alpar@2409
   162
      dmap[nodes[0]][0] = PathData(true, 0);
alpar@2409
   163
      process.clear();
deba@2413
   164
      // Processing the first round
deba@2413
   165
      for (OutEdgeIt e(graph, nodes[0]); e != INVALID; ++e) {
deba@2413
   166
	Node v = graph.target(e);
alpar@2409
   167
	if (comp[v] != comp_cnt || v == nodes[0]) continue;
deba@2413
   168
	dmap[v][1] = PathData(true, length[e], e);
alpar@2409
   169
	process.push_back(v);
alpar@2409
   170
      }
deba@2437
   171
      // Processing other rounds
alpar@2409
   172
      int n = nodes.size(), k;
deba@2437
   173
      for (k = 2; k <= n && process.size() < n; ++k)
alpar@2409
   174
	processNextBuildRound(k);
deba@2437
   175
      for ( ; k <= n; ++k)
alpar@2409
   176
	processNextFullRound(k);
alpar@2409
   177
    }
deba@2437
   178
alpar@2409
   179
    /// \brief Processes one round of computing required path data and
alpar@2409
   180
    /// rebuilds \ref process vector.
alpar@2409
   181
    void processNextBuildRound(int k) {
alpar@2409
   182
      NodeVector next;
deba@2437
   183
      for (int i = 0; i < process.size(); ++i) {
deba@2437
   184
	for (OutEdgeIt e(graph, process[i]); e != INVALID; ++e) {
deba@2413
   185
	  Node v = graph.target(e);
alpar@2409
   186
	  if (comp[v] != comp_cnt) continue;
deba@2413
   187
	  if (!dmap[v][k].found) {
deba@2413
   188
	    next.push_back(v);
deba@2437
   189
	    dmap[v][k] = PathData(true, dmap[process[i]][k-1].dist + length[e], e);
deba@2413
   190
	  }
deba@2437
   191
	  else if (dmap[process[i]][k-1].dist + length[e] < dmap[v][k].dist) {
deba@2437
   192
	    dmap[v][k] = PathData(true, dmap[process[i]][k-1].dist + length[e], e);
alpar@2409
   193
	  }
alpar@2409
   194
	}
alpar@2409
   195
      }
alpar@2409
   196
      process.swap(next);
alpar@2409
   197
    }
alpar@2409
   198
alpar@2409
   199
    /// \brief Processes one round of computing required path data
alpar@2409
   200
    /// using \ref nodes vector instead of \ref process vector.
alpar@2409
   201
    void processNextFullRound(int k) {
deba@2437
   202
      for (int i = 0; i < nodes.size(); ++i) {
deba@2437
   203
	for (OutEdgeIt e(graph, nodes[i]); e != INVALID; ++e) {
deba@2413
   204
	  Node v = graph.target(e);
alpar@2409
   205
	  if (comp[v] != comp_cnt) continue;
deba@2437
   206
	  if ( !dmap[v][k].found ||
deba@2437
   207
	       dmap[nodes[i]][k-1].dist + length[e] < dmap[v][k].dist ) {
deba@2437
   208
	    dmap[v][k] = PathData(true, dmap[nodes[i]][k-1].dist + length[e], e);
alpar@2409
   209
	  }
alpar@2409
   210
	}
alpar@2409
   211
      }
alpar@2409
   212
    }
deba@2437
   213
deba@2437
   214
    /// \brief Finds the minimum cycle mean value in the current
deba@2413
   215
    /// component.
deba@2413
   216
    bool findCurrentMin(Length &min_length, int &min_size, Node &min_node) {
alpar@2409
   217
      bool found_min = false;
deba@2437
   218
      for (int i = 0; i < nodes.size(); ++i) {
deba@2413
   219
	int n = nodes.size();
deba@2437
   220
	if (!dmap[nodes[i]][n].found) continue;
alpar@2409
   221
	Length len;
alpar@2409
   222
	int size;
alpar@2409
   223
	bool found_one = false;
alpar@2409
   224
	for (int k = 0; k < n; ++k) {
deba@2437
   225
	  if (!dmap[nodes[i]][k].found) continue;
deba@2437
   226
	  Length _len = dmap[nodes[i]][n].dist - dmap[nodes[i]][k].dist;
deba@2437
   227
	  int _size = n - k;
deba@2413
   228
	  if (!found_one || len * _size < _len * size) {
alpar@2409
   229
	    found_one = true;
alpar@2409
   230
	    len = _len;
alpar@2409
   231
	    size = _size;
alpar@2409
   232
	  }
alpar@2409
   233
	}
deba@2437
   234
	if ( found_one &&
deba@2413
   235
	     (!found_min || len * min_size < min_length * size) ) {
alpar@2409
   236
	  found_min = true;
alpar@2409
   237
	  min_length = len;
alpar@2409
   238
	  min_size = size;
deba@2437
   239
	  min_node = nodes[i];
alpar@2409
   240
	}
alpar@2409
   241
      }
alpar@2409
   242
      return found_min;
alpar@2409
   243
    }
deba@2437
   244
deba@2437
   245
  public:
deba@2437
   246
alpar@2409
   247
    /// \brief Runs the algorithm.
alpar@2409
   248
    ///
alpar@2409
   249
    /// Runs the algorithm.
alpar@2409
   250
    ///
alpar@2409
   251
    /// \return \c true if a cycle exists in the graph.
deba@2413
   252
    ///
deba@2437
   253
    /// \note Apart from the return value, m.run() is just a shortcut
deba@2413
   254
    /// of the following code.
deba@2413
   255
    /// \code
deba@2413
   256
    ///   m.init();
deba@2413
   257
    ///   m.findMinMean();
deba@2413
   258
    ///   m.findCycle();
deba@2413
   259
    /// \endcode
alpar@2409
   260
    bool run() {
alpar@2409
   261
      init();
deba@2413
   262
      findMinMean();
deba@2413
   263
      return findCycle();
deba@2413
   264
    }
deba@2437
   265
deba@2413
   266
    /// \brief Initializes the internal data structures.
deba@2413
   267
    void init() {
deba@2413
   268
      comp_num = stronglyConnectedComponents(graph, comp);
deba@2413
   269
      if (!cycle_path) {
deba@2413
   270
	local_path = true;
deba@2413
   271
	cycle_path = new Path;
deba@2413
   272
      }
deba@2413
   273
    }
deba@2413
   274
deba@2413
   275
    /// \brief Finds the minimum cycle mean value in the graph.
deba@2413
   276
    ///
deba@2413
   277
    /// Computes all the required path data and finds the minimum cycle
deba@2413
   278
    /// mean value in the graph.
deba@2413
   279
    ///
deba@2413
   280
    /// \return \c true if a cycle exists in the graph.
deba@2437
   281
    ///
deba@2413
   282
    /// \pre \ref init() must be called before using this function.
deba@2413
   283
    bool findMinMean() {
deba@2413
   284
      cycle_node = INVALID;
alpar@2409
   285
      for (comp_cnt = 0; comp_cnt < comp_num; ++comp_cnt) {
alpar@2409
   286
	initCurrent();
alpar@2409
   287
	processRounds();
deba@2437
   288
alpar@2409
   289
	Length min_length;
alpar@2409
   290
	int min_size;
alpar@2409
   291
	Node min_node;
deba@2413
   292
	bool found_min = findCurrentMin(min_length, min_size, min_node);
deba@2437
   293
deba@2437
   294
	if ( found_min && (cycle_node == INVALID ||
deba@2413
   295
	     min_length * cycle_size < cycle_length * min_size) ) {
alpar@2409
   296
	  cycle_length = min_length;
alpar@2409
   297
	  cycle_size = min_size;
alpar@2409
   298
	  cycle_node = min_node;
alpar@2409
   299
	}
alpar@2409
   300
      }
deba@2413
   301
      return (cycle_node != INVALID);
alpar@2409
   302
    }
deba@2437
   303
deba@2413
   304
    /// \brief Finds a critical (minimum mean) cycle.
alpar@2409
   305
    ///
deba@2413
   306
    /// Finds a critical (minimum mean) cycle using the path data
deba@2413
   307
    /// stored in \ref dmap.
deba@2413
   308
    ///
deba@2413
   309
    /// \return \c true if a cycle exists in the graph.
deba@2437
   310
    ///
deba@2437
   311
    /// \pre \ref init() and \ref findMinMean() must be called before
deba@2413
   312
    /// using this function.
deba@2413
   313
    bool findCycle() {
deba@2413
   314
      if (cycle_node == INVALID) return false;
deba@2413
   315
      cycle_length = 0;
deba@2413
   316
      cycle_size = 0;
deba@2413
   317
      IntNodeMap reached(graph, -1);
deba@2413
   318
      int r = reached[cycle_node] = dmap[cycle_node].size() - 1;
deba@2413
   319
      Node u = graph.source(dmap[cycle_node][r].pred);
deba@2413
   320
      while (reached[u] < 0) {
deba@2413
   321
	reached[u] = --r;
deba@2413
   322
	u = graph.source(dmap[u][r].pred);
deba@2413
   323
      }
deba@2413
   324
      r = reached[u];
deba@2413
   325
      Edge e = dmap[u][r].pred;
deba@2413
   326
      cycle_path->addFront(e);
deba@2413
   327
      cycle_length = cycle_length + length[e];
deba@2413
   328
      ++cycle_size;
deba@2413
   329
      Node v;
deba@2413
   330
      while ((v = graph.source(e)) != u) {
deba@2413
   331
	e = dmap[v][--r].pred;
deba@2413
   332
	cycle_path->addFront(e);
deba@2413
   333
	cycle_length = cycle_length + length[e];
deba@2413
   334
	++cycle_size;
deba@2413
   335
      }
deba@2413
   336
      return true;
deba@2413
   337
    }
deba@2413
   338
deba@2413
   339
    /// \brief Resets the internal data structures.
deba@2413
   340
    ///
deba@2437
   341
    /// Resets the internal data structures so that \ref findMinMean()
deba@2437
   342
    /// and \ref findCycle() can be called again (e.g. when the
deba@2413
   343
    /// underlaying graph has been modified).
deba@2413
   344
    void reset() {
deba@2413
   345
      for (NodeIt u(graph); u != INVALID; ++u)
deba@2413
   346
	dmap[u].clear();
deba@2413
   347
      cycle_node = INVALID;
deba@2413
   348
      if (cycle_path) cycle_path->clear();
deba@2413
   349
      comp_num = stronglyConnectedComponents(graph, comp);
deba@2413
   350
    }
deba@2437
   351
deba@2413
   352
    /// \brief Returns the total length of the found cycle.
deba@2413
   353
    ///
deba@2413
   354
    /// Returns the total length of the found cycle.
alpar@2409
   355
    ///
alpar@2409
   356
    /// \pre \ref run() must be called before using this function.
deba@2437
   357
    Length cycleLength() const {
alpar@2409
   358
      return cycle_length;
alpar@2409
   359
    }
deba@2437
   360
deba@2413
   361
    /// \brief Returns the number of edges in the found cycle.
alpar@2409
   362
    ///
deba@2413
   363
    /// Returns the number of edges in the found cycle.
alpar@2409
   364
    ///
alpar@2409
   365
    /// \pre \ref run() must be called before using this function.
deba@2437
   366
    int cycleEdgeNum() const {
alpar@2409
   367
      return cycle_size;
alpar@2409
   368
    }
deba@2437
   369
deba@2413
   370
    /// \brief Returns the mean length of the found cycle.
alpar@2409
   371
    ///
deba@2413
   372
    /// Returns the mean length of the found cycle.
alpar@2409
   373
    ///
alpar@2409
   374
    /// \pre \ref run() must be called before using this function.
alpar@2409
   375
    ///
alpar@2409
   376
    /// \warning LengthMap::Value must be convertible to double.
deba@2413
   377
    ///
deba@2413
   378
    /// \note m.minMean() is just a shortcut of the following code.
deba@2413
   379
    /// \code
deba@2413
   380
    ///   return m.cycleEdgeNum() / double(m.cycleLength());
deba@2413
   381
    /// \endcode
deba@2437
   382
    double minMean() const {
deba@2413
   383
      return cycle_length / (double)cycle_size;
alpar@2409
   384
    }
alpar@2409
   385
alpar@2409
   386
    /// \brief Returns a const reference to the \ref lemon::Path "Path"
deba@2413
   387
    /// structure of the found cycle.
alpar@2409
   388
    ///
alpar@2409
   389
    /// Returns a const reference to the \ref lemon::Path "Path"
deba@2413
   390
    /// structure of the found cycle.
alpar@2409
   391
    ///
alpar@2409
   392
    /// \pre \ref run() must be called before using this function.
alpar@2409
   393
    ///
alpar@2409
   394
    /// \sa \ref cyclePath()
deba@2437
   395
    const Path& cycle() const {
alpar@2409
   396
      return *cycle_path;
alpar@2409
   397
    }
deba@2437
   398
deba@2437
   399
    /// \brief Sets the \ref lemon::Path "Path" structure storing the
alpar@2409
   400
    /// found cycle.
deba@2437
   401
    ///
deba@2437
   402
    /// Sets the \ref lemon::Path "Path" structure storing the found
deba@2437
   403
    /// cycle. If you don't use this function before calling
deba@2437
   404
    /// \ref run(), it will allocate one. The destuctor deallocates
alpar@2409
   405
    /// this automatically allocated map, of course.
alpar@2409
   406
    ///
alpar@2409
   407
    /// \note The algorithm calls only the \ref lemon::Path::addFront()
deba@2437
   408
    /// "addFront()" method of the given \ref lemon::Path "Path"
alpar@2409
   409
    /// structure.
deba@2437
   410
    ///
alpar@2409
   411
    /// \return \c (*this)
deba@2437
   412
    MinMeanCycle& cyclePath(Path &path) {
alpar@2409
   413
      if (local_path) {
alpar@2409
   414
	delete cycle_path;
alpar@2409
   415
	local_path = false;
alpar@2409
   416
      }
alpar@2409
   417
      cycle_path = &path;
alpar@2409
   418
      return *this;
alpar@2409
   419
    }
alpar@2409
   420
alpar@2409
   421
  }; //class MinMeanCycle
alpar@2409
   422
alpar@2409
   423
  ///@}
alpar@2409
   424
alpar@2409
   425
} //namespace lemon
alpar@2409
   426
alpar@2409
   427
#endif //LEMON_MIN_MEAN_CYCLE_H