lemon/min_mean_cycle.h
changeset 2555 a84e52e99f57
parent 2553 bfced05fa852
child 2583 7216b6a52ab9
     1.1 --- a/lemon/min_mean_cycle.h	Mon Jan 07 17:07:40 2008 +0000
     1.2 +++ b/lemon/min_mean_cycle.h	Sun Jan 13 10:26:55 2008 +0000
     1.3 @@ -22,224 +22,242 @@
     1.4  /// \ingroup shortest_path
     1.5  ///
     1.6  /// \file
     1.7 -/// \brief Karp's algorithm for finding a minimum mean (directed) cycle.
     1.8 +/// \brief Howard's algorithm for finding a minimum mean directed cycle.
     1.9  
    1.10  #include <vector>
    1.11  #include <lemon/graph_utils.h>
    1.12  #include <lemon/topology.h>
    1.13  #include <lemon/path.h>
    1.14 +#include <lemon/tolerance.h>
    1.15  
    1.16  namespace lemon {
    1.17  
    1.18    /// \addtogroup shortest_path
    1.19    /// @{
    1.20  
    1.21 -  /// \brief Implementation of Karp's algorithm for finding a
    1.22 -  /// minimum mean (directed) cycle.
    1.23 +  /// \brief Implementation of Howard's algorithm for finding a minimum
    1.24 +  /// mean directed cycle.
    1.25    ///
    1.26 -  /// The \ref lemon::MinMeanCycle "MinMeanCycle" implements Karp's
    1.27 -  /// algorithm for finding a minimum mean (directed) cycle.
    1.28 +  /// \ref MinMeanCycle implements Howard's algorithm for finding a
    1.29 +  /// minimum mean directed cycle.
    1.30    ///
    1.31    /// \param Graph The directed graph type the algorithm runs on.
    1.32    /// \param LengthMap The type of the length (cost) map.
    1.33    ///
    1.34 +  /// \warning \c LengthMap::Value must be convertible to \c double.
    1.35 +  ///
    1.36    /// \author Peter Kovacs
    1.37  
    1.38 -  template <typename Graph,
    1.39 -    typename LengthMap = typename Graph::template EdgeMap<int> >
    1.40 +  template < typename Graph,
    1.41 +             typename LengthMap = typename Graph::template EdgeMap<int> >
    1.42    class MinMeanCycle
    1.43    {
    1.44 -    typedef typename Graph::Node Node;
    1.45 -    typedef typename Graph::NodeIt NodeIt;
    1.46 -    typedef typename Graph::Edge Edge;
    1.47 -    typedef typename Graph::EdgeIt EdgeIt;
    1.48 -    typedef typename Graph::OutEdgeIt OutEdgeIt;
    1.49 +    GRAPH_TYPEDEFS(typename Graph);
    1.50  
    1.51      typedef typename LengthMap::Value Length;
    1.52 -
    1.53 -    typedef typename Graph::template NodeMap<int> IntNodeMap;
    1.54 -    typedef typename Graph::template NodeMap<Edge> PredNodeMap;
    1.55      typedef Path<Graph> Path;
    1.56 -    typedef std::vector<Node> NodeVector;
    1.57  
    1.58    protected:
    1.59  
    1.60 -    /// \brief Data structure for path data.
    1.61 -    struct PathData
    1.62 -    {
    1.63 -      bool found;
    1.64 -      Length dist;
    1.65 -      Edge pred;
    1.66 -      PathData(bool _found = false, Length _dist = 0) :
    1.67 -	found(_found), dist(_dist), pred(INVALID) {}
    1.68 -      PathData(bool _found, Length _dist, Edge _pred) :
    1.69 -	found(_found), dist(_dist), pred(_pred) {}
    1.70 -    };
    1.71 +    typename Graph::template NodeMap<bool> _reached;
    1.72 +    typename Graph::template NodeMap<double> _dist;
    1.73 +    typename Graph::template NodeMap<Edge> _policy;
    1.74  
    1.75 -  private:
    1.76 +    /// The directed graph the algorithm runs on.
    1.77 +    const Graph &_graph;
    1.78 +    /// The length of the edges.
    1.79 +    const LengthMap &_length;
    1.80  
    1.81 -    typedef typename Graph::template NodeMap<std::vector<PathData> >
    1.82 -      PathDataNodeMap;
    1.83 +    /// The total length of the found cycle.
    1.84 +    Length _cycle_length;
    1.85 +    /// The number of edges on the found cycle.
    1.86 +    int _cycle_size;
    1.87 +    /// The found cycle.
    1.88 +    Path *_cycle_path;
    1.89  
    1.90 -  protected:
    1.91 +    bool _local_path;
    1.92 +    bool _cycle_found;
    1.93 +    Node _cycle_node;
    1.94  
    1.95 -    /// \brief Node map for storing path data.
    1.96 -    ///
    1.97 -    /// Node map for storing path data of all nodes in the current
    1.98 -    /// component. dmap[v][k] is the length of a shortest directed walk
    1.99 -    /// to node v from the starting node containing exactly k edges.
   1.100 -    PathDataNodeMap dmap;
   1.101 +    typename Graph::template NodeMap<int> _component;
   1.102 +    int _component_num;
   1.103  
   1.104 -    /// \brief The directed graph the algorithm runs on.
   1.105 -    const Graph &graph;
   1.106 -    /// \brief The length of the edges.
   1.107 -    const LengthMap &length;
   1.108 +    std::vector<Node> _nodes;
   1.109 +    std::vector<Edge> _edges;
   1.110 +    Tolerance<double> _tolerance;
   1.111  
   1.112 -    /// \brief The total length of the found cycle.
   1.113 -    Length cycle_length;
   1.114 -    /// \brief The number of edges in the found cycle.
   1.115 -    int cycle_size;
   1.116 -    /// \brief A node for obtaining a minimum mean cycle.
   1.117 -    Node cycle_node;
   1.118 -
   1.119 -    /// \brief The found cycle.
   1.120 -    Path *cycle_path;
   1.121 -    /// \brief The algorithm uses local \ref lemon::Path "Path"
   1.122 -    /// structure to store the found cycle.
   1.123 -    bool local_path;
   1.124 -
   1.125 -    /// \brief Node map for identifying strongly connected components.
   1.126 -    IntNodeMap comp;
   1.127 -    /// \brief The number of strongly connected components.
   1.128 -    int comp_num;
   1.129 -    /// \brief Counter for identifying the current component.
   1.130 -    int comp_cnt;
   1.131 -    /// \brief Nodes of the current component.
   1.132 -    NodeVector nodes;
   1.133 -    /// \brief The processed nodes in the last round.
   1.134 -    NodeVector process;
   1.135 -
   1.136 -  public :
   1.137 +  public:
   1.138  
   1.139      /// \brief The constructor of the class.
   1.140      ///
   1.141      /// The constructor of the class.
   1.142      ///
   1.143 -    /// \param _graph The directed graph the algorithm runs on.
   1.144 -    /// \param _length The length (cost) of the edges.
   1.145 -    MinMeanCycle( const Graph &_graph,
   1.146 -		  const LengthMap &_length ) :
   1.147 -      graph(_graph), length(_length), dmap(_graph), comp(_graph),
   1.148 -      cycle_length(0), cycle_size(-1), cycle_node(INVALID),
   1.149 -      cycle_path(NULL), local_path(false)
   1.150 +    /// \param graph The directed graph the algorithm runs on.
   1.151 +    /// \param length The length (cost) of the edges.
   1.152 +    MinMeanCycle( const Graph &graph,
   1.153 +                  const LengthMap &length ) :
   1.154 +      _graph(graph), _length(length), _dist(graph), _reached(graph),
   1.155 +      _policy(graph), _component(graph), _cycle_length(0),
   1.156 +      _cycle_size(-1), _cycle_path(NULL), _local_path(false)
   1.157      { }
   1.158  
   1.159 -    /// \brief The destructor of the class.
   1.160 +    /// The destructor of the class.
   1.161      ~MinMeanCycle() {
   1.162 -      if (local_path) delete cycle_path;
   1.163 +      if (_local_path) delete _cycle_path;
   1.164      }
   1.165  
   1.166    protected:
   1.167  
   1.168 -    /// \brief Initializes the internal data structures for the current
   1.169 -    /// component.
   1.170 -    void initCurrent() {
   1.171 -      nodes.clear();
   1.172 +    // Initializes the internal data structures for the current strongly
   1.173 +    // connected component and creating the policy graph.
   1.174 +    // The policy graph can be represented by the _policy map because
   1.175 +    // the out degree of every node is 1.
   1.176 +    bool initCurrentComponent(int comp) {
   1.177        // Finding the nodes of the current component
   1.178 -      for (NodeIt v(graph); v != INVALID; ++v) {
   1.179 -	if (comp[v] == comp_cnt) nodes.push_back(v);
   1.180 +      _nodes.clear();
   1.181 +      for (NodeIt n(_graph); n != INVALID; ++n) {
   1.182 +        if (_component[n] == comp) _nodes.push_back(n);
   1.183        }
   1.184 -      // Creating vectors for all nodes
   1.185 -      int n = nodes.size();
   1.186 -      for (int i = 0; i < n; ++i) {
   1.187 -	dmap[nodes[i]].resize(n + 1);
   1.188 +      if (_nodes.size() <= 1) return false;
   1.189 +      // Finding the edges of the current component
   1.190 +      _edges.clear();
   1.191 +      for (EdgeIt e(_graph); e != INVALID; ++e) {
   1.192 +        if ( _component[_graph.source(e)] == comp &&
   1.193 +             _component[_graph.target(e)] == comp )
   1.194 +          _edges.push_back(e);
   1.195 +      }
   1.196 +      // Initializing _reached, _dist, _policy maps
   1.197 +      for (int i = 0; i < _nodes.size(); ++i) {
   1.198 +        _reached[_nodes[i]] = false;
   1.199 +        _policy[_nodes[i]] = INVALID;
   1.200 +      }
   1.201 +      Node u; Edge e;
   1.202 +      for (int j = 0; j < _edges.size(); ++j) {
   1.203 +        e = _edges[j];
   1.204 +        u = _graph.source(e);
   1.205 +        if (!_reached[u] || _length[e] < _dist[u]) {
   1.206 +          _dist[u] = _length[e];
   1.207 +          _policy[u] = e;
   1.208 +          _reached[u] = true;
   1.209 +        }
   1.210 +      }
   1.211 +      return true;
   1.212 +    }
   1.213 +
   1.214 +    // Finds all cycles in the policy graph.
   1.215 +    // Sets _cycle_found to true if a cycle is found and sets
   1.216 +    // _cycle_length, _cycle_size, _cycle_node to represent the minimum
   1.217 +    // mean cycle in the policy graph.
   1.218 +    bool findPolicyCycles(int comp) {
   1.219 +      typename Graph::template NodeMap<int> level(_graph, -1);
   1.220 +      _cycle_found = false;
   1.221 +      Length clength;
   1.222 +      int csize;
   1.223 +      int path_cnt = 0;
   1.224 +      Node u, v;
   1.225 +      // Searching for cycles
   1.226 +      for (int i = 0; i < _nodes.size(); ++i) {
   1.227 +        if (level[_nodes[i]] < 0) {
   1.228 +          u = _nodes[i];
   1.229 +          level[u] = path_cnt;
   1.230 +          while (level[u = _graph.target(_policy[u])] < 0)
   1.231 +            level[u] = path_cnt;
   1.232 +          if (level[u] == path_cnt) {
   1.233 +            // A cycle is found
   1.234 +            clength = _length[_policy[u]];
   1.235 +            csize = 1;
   1.236 +            for (v = u; (v = _graph.target(_policy[v])) != u; ) {
   1.237 +              clength += _length[_policy[v]];
   1.238 +              ++csize;
   1.239 +            }
   1.240 +            if ( !_cycle_found ||
   1.241 +                 clength * _cycle_size < _cycle_length * csize ) {
   1.242 +              _cycle_found = true;
   1.243 +              _cycle_length = clength;
   1.244 +              _cycle_size = csize;
   1.245 +              _cycle_node = u;
   1.246 +            }
   1.247 +          }
   1.248 +          ++path_cnt;
   1.249 +        }
   1.250 +      }
   1.251 +      return _cycle_found;
   1.252 +    }
   1.253 +
   1.254 +    // Contracts the policy graph to be connected by cutting all cycles
   1.255 +    // except for the main cycle (i.e. the minimum mean cycle).
   1.256 +    void contractPolicyGraph(int comp) {
   1.257 +      // Finding the component of the main cycle using
   1.258 +      // reverse BFS search
   1.259 +      typename Graph::template NodeMap<int> found(_graph, false);
   1.260 +      std::deque<Node> queue;
   1.261 +      queue.push_back(_cycle_node);
   1.262 +      found[_cycle_node] = true;
   1.263 +      Node u, v;
   1.264 +      while (!queue.empty()) {
   1.265 +        v = queue.front(); queue.pop_front();
   1.266 +        for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
   1.267 +          u = _graph.source(e);
   1.268 +          if (_component[u] == comp && !found[u] && _policy[u] == e) {
   1.269 +            found[u] = true;
   1.270 +            queue.push_back(u);
   1.271 +          }
   1.272 +        }
   1.273 +      }
   1.274 +      // Connecting all other nodes to this component using
   1.275 +      // reverse BFS search
   1.276 +      queue.clear();
   1.277 +      for (int i = 0; i < _nodes.size(); ++i)
   1.278 +        if (found[_nodes[i]]) queue.push_back(_nodes[i]);
   1.279 +      int found_cnt = queue.size();
   1.280 +      while (found_cnt < _nodes.size() && !queue.empty()) {
   1.281 +        v = queue.front(); queue.pop_front();
   1.282 +        for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
   1.283 +          u = _graph.source(e);
   1.284 +          if (_component[u] == comp && !found[u]) {
   1.285 +            found[u] = true;
   1.286 +            ++found_cnt;
   1.287 +            _policy[u] = e;
   1.288 +            queue.push_back(u);
   1.289 +          }
   1.290 +        }
   1.291        }
   1.292      }
   1.293  
   1.294 -    /// \brief Processes all rounds of computing required path data for
   1.295 -    /// the current component.
   1.296 -    void processRounds() {
   1.297 -      dmap[nodes[0]][0] = PathData(true, 0);
   1.298 -      process.clear();
   1.299 -      // Processing the first round
   1.300 -      for (OutEdgeIt e(graph, nodes[0]); e != INVALID; ++e) {
   1.301 -	Node v = graph.target(e);
   1.302 -	if (comp[v] != comp_cnt || v == nodes[0]) continue;
   1.303 -	dmap[v][1] = PathData(true, length[e], e);
   1.304 -	process.push_back(v);
   1.305 +    // Computes node distances in the policy graph and updates the
   1.306 +    // policy graph if the node distances can be improved.
   1.307 +    bool computeNodeDistances(int comp) {
   1.308 +      // Computing node distances using reverse BFS search
   1.309 +      double cycle_mean = (double)_cycle_length / _cycle_size;
   1.310 +      typename Graph::template NodeMap<int> found(_graph, false);
   1.311 +      std::deque<Node> queue;
   1.312 +      queue.push_back(_cycle_node);
   1.313 +      found[_cycle_node] = true;
   1.314 +      Node u, v;
   1.315 +      while (!queue.empty()) {
   1.316 +        v = queue.front(); queue.pop_front();
   1.317 +        for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
   1.318 +          u = _graph.source(e);
   1.319 +          if (_component[u] == comp && !found[u] && _policy[u] == e) {
   1.320 +            found[u] = true;
   1.321 +            _dist[u] = _dist[v] + _length[e] - cycle_mean;
   1.322 +            queue.push_back(u);
   1.323 +          }
   1.324 +        }
   1.325        }
   1.326 -      // Processing other rounds
   1.327 -      int n = nodes.size(), k;
   1.328 -      for (k = 2; k <= n && process.size() < n; ++k)
   1.329 -	processNextBuildRound(k);
   1.330 -      for ( ; k <= n; ++k)
   1.331 -	processNextFullRound(k);
   1.332 -    }
   1.333 -
   1.334 -    /// \brief Processes one round of computing required path data and
   1.335 -    /// rebuilds \ref process vector.
   1.336 -    void processNextBuildRound(int k) {
   1.337 -      NodeVector next;
   1.338 -      for (int i = 0; i < process.size(); ++i) {
   1.339 -	for (OutEdgeIt e(graph, process[i]); e != INVALID; ++e) {
   1.340 -	  Node v = graph.target(e);
   1.341 -	  if (comp[v] != comp_cnt) continue;
   1.342 -	  if (!dmap[v][k].found) {
   1.343 -	    next.push_back(v);
   1.344 -	    dmap[v][k] = PathData(true, dmap[process[i]][k-1].dist + length[e], e);
   1.345 -	  }
   1.346 -	  else if (dmap[process[i]][k-1].dist + length[e] < dmap[v][k].dist) {
   1.347 -	    dmap[v][k] = PathData(true, dmap[process[i]][k-1].dist + length[e], e);
   1.348 -	  }
   1.349 -	}
   1.350 +      // Improving node distances
   1.351 +      bool improved = false;
   1.352 +      for (int j = 0; j < _edges.size(); ++j) {
   1.353 +        Edge e = _edges[j];
   1.354 +        u = _graph.source(e); v = _graph.target(e);
   1.355 +        double delta = _dist[v] + _length[e] - cycle_mean;
   1.356 +        if (_tolerance.less(delta, _dist[u])) {
   1.357 +          improved = true;
   1.358 +          _dist[u] = delta;
   1.359 +          _policy[u] = e;
   1.360 +        }
   1.361        }
   1.362 -      process.swap(next);
   1.363 -    }
   1.364 -
   1.365 -    /// \brief Processes one round of computing required path data
   1.366 -    /// using \ref nodes vector instead of \ref process vector.
   1.367 -    void processNextFullRound(int k) {
   1.368 -      for (int i = 0; i < nodes.size(); ++i) {
   1.369 -	for (OutEdgeIt e(graph, nodes[i]); e != INVALID; ++e) {
   1.370 -	  Node v = graph.target(e);
   1.371 -	  if (comp[v] != comp_cnt) continue;
   1.372 -	  if ( !dmap[v][k].found ||
   1.373 -	       dmap[nodes[i]][k-1].dist + length[e] < dmap[v][k].dist ) {
   1.374 -	    dmap[v][k] = PathData(true, dmap[nodes[i]][k-1].dist + length[e], e);
   1.375 -	  }
   1.376 -	}
   1.377 -      }
   1.378 -    }
   1.379 -
   1.380 -    /// \brief Finds the minimum cycle mean value in the current
   1.381 -    /// component.
   1.382 -    bool findCurrentMin(Length &min_length, int &min_size, Node &min_node) {
   1.383 -      bool found_min = false;
   1.384 -      for (int i = 0; i < nodes.size(); ++i) {
   1.385 -	int n = nodes.size();
   1.386 -	if (!dmap[nodes[i]][n].found) continue;
   1.387 -	Length len;
   1.388 -	int size;
   1.389 -	bool found_one = false;
   1.390 -	for (int k = 0; k < n; ++k) {
   1.391 -	  if (!dmap[nodes[i]][k].found) continue;
   1.392 -	  Length _len = dmap[nodes[i]][n].dist - dmap[nodes[i]][k].dist;
   1.393 -	  int _size = n - k;
   1.394 -	  if (!found_one || len * _size < _len * size) {
   1.395 -	    found_one = true;
   1.396 -	    len = _len;
   1.397 -	    size = _size;
   1.398 -	  }
   1.399 -	}
   1.400 -	if ( found_one &&
   1.401 -	     (!found_min || len * min_size < min_length * size) ) {
   1.402 -	  found_min = true;
   1.403 -	  min_length = len;
   1.404 -	  min_size = size;
   1.405 -	  min_node = nodes[i];
   1.406 -	}
   1.407 -      }
   1.408 -      return found_min;
   1.409 +      return improved;
   1.410      }
   1.411  
   1.412    public:
   1.413 @@ -248,19 +266,18 @@
   1.414      ///
   1.415      /// Runs the algorithm.
   1.416      ///
   1.417 -    /// \return \c true if a cycle exists in the graph.
   1.418 +    /// \return Returns \c true if a directed cycle exists in the graph.
   1.419      ///
   1.420 -    /// \note Apart from the return value, <tt>m.run()</tt> is just a
   1.421 +    /// \note Apart from the return value, <tt>mmc.run()</tt> is just a
   1.422      /// shortcut of the following code.
   1.423      /// \code
   1.424 -    ///   m.init();
   1.425 -    ///   m.findMinMean();
   1.426 -    ///   m.findCycle();
   1.427 +    ///   mmc.init();
   1.428 +    ///   mmc.findMinMean();
   1.429 +    ///   mmc.findCycle();
   1.430      /// \endcode
   1.431      bool run() {
   1.432        init();
   1.433 -      findMinMean();
   1.434 -      return findCycle();
   1.435 +      return findMinMean() && findCycle();
   1.436      }
   1.437  
   1.438      /// \brief Initializes the internal data structures.
   1.439 @@ -269,75 +286,13 @@
   1.440      ///
   1.441      /// \sa reset()
   1.442      void init() {
   1.443 -      comp_num = stronglyConnectedComponents(graph, comp);
   1.444 -      if (!cycle_path) {
   1.445 -	local_path = true;
   1.446 -	cycle_path = new Path;
   1.447 +      _tolerance.epsilon(1e-8);
   1.448 +      if (!_cycle_path) {
   1.449 +        _local_path = true;
   1.450 +        _cycle_path = new Path;
   1.451        }
   1.452 -    }
   1.453 -
   1.454 -    /// \brief Finds the minimum cycle mean value in the graph.
   1.455 -    ///
   1.456 -    /// Computes all the required path data and finds the minimum cycle
   1.457 -    /// mean value in the graph.
   1.458 -    ///
   1.459 -    /// \return \c true if a cycle exists in the graph.
   1.460 -    ///
   1.461 -    /// \pre \ref init() must be called before using this function.
   1.462 -    bool findMinMean() {
   1.463 -      cycle_node = INVALID;
   1.464 -      for (comp_cnt = 0; comp_cnt < comp_num; ++comp_cnt) {
   1.465 -	initCurrent();
   1.466 -	processRounds();
   1.467 -
   1.468 -	Length min_length;
   1.469 -	int min_size;
   1.470 -	Node min_node;
   1.471 -	bool found_min = findCurrentMin(min_length, min_size, min_node);
   1.472 -
   1.473 -	if ( found_min && (cycle_node == INVALID ||
   1.474 -	     min_length * cycle_size < cycle_length * min_size) ) {
   1.475 -	  cycle_length = min_length;
   1.476 -	  cycle_size = min_size;
   1.477 -	  cycle_node = min_node;
   1.478 -	}
   1.479 -      }
   1.480 -      return (cycle_node != INVALID);
   1.481 -    }
   1.482 -
   1.483 -    /// \brief Finds a critical (minimum mean) cycle.
   1.484 -    ///
   1.485 -    /// Finds a critical (minimum mean) cycle using the path data
   1.486 -    /// stored in \ref dmap.
   1.487 -    ///
   1.488 -    /// \return \c true if a cycle exists in the graph.
   1.489 -    ///
   1.490 -    /// \pre \ref init() and \ref findMinMean() must be called before
   1.491 -    /// using this function.
   1.492 -    bool findCycle() {
   1.493 -      if (cycle_node == INVALID) return false;
   1.494 -      cycle_length = 0;
   1.495 -      cycle_size = 0;
   1.496 -      IntNodeMap reached(graph, -1);
   1.497 -      int r = reached[cycle_node] = dmap[cycle_node].size() - 1;
   1.498 -      Node u = graph.source(dmap[cycle_node][r].pred);
   1.499 -      while (reached[u] < 0) {
   1.500 -	reached[u] = --r;
   1.501 -	u = graph.source(dmap[u][r].pred);
   1.502 -      }
   1.503 -      r = reached[u];
   1.504 -      Edge e = dmap[u][r].pred;
   1.505 -      cycle_path->addFront(e);
   1.506 -      cycle_length = cycle_length + length[e];
   1.507 -      ++cycle_size;
   1.508 -      Node v;
   1.509 -      while ((v = graph.source(e)) != u) {
   1.510 -	e = dmap[v][--r].pred;
   1.511 -	cycle_path->addFront(e);
   1.512 -	cycle_length = cycle_length + length[e];
   1.513 -	++cycle_size;
   1.514 -      }
   1.515 -      return true;
   1.516 +      _cycle_found = false;
   1.517 +      _component_num = stronglyConnectedComponents(_graph, _component);
   1.518      }
   1.519  
   1.520      /// \brief Resets the internal data structures.
   1.521 @@ -348,33 +303,68 @@
   1.522      ///
   1.523      /// \sa init()
   1.524      void reset() {
   1.525 -      for (NodeIt u(graph); u != INVALID; ++u)
   1.526 -	dmap[u].clear();
   1.527 -      cycle_node = INVALID;
   1.528 -      if (cycle_path) cycle_path->clear();
   1.529 -      comp_num = stronglyConnectedComponents(graph, comp);
   1.530 +      if (_cycle_path) _cycle_path->clear();
   1.531 +      _cycle_found = false;
   1.532 +      _component_num = stronglyConnectedComponents(_graph, _component);
   1.533 +    }
   1.534 +
   1.535 +    /// \brief Finds the minimum cycle mean length in the graph.
   1.536 +    ///
   1.537 +    /// Computes all the required data and finds the minimum cycle mean
   1.538 +    /// length in the graph.
   1.539 +    ///
   1.540 +    /// \return Returns \c true if a directed cycle exists in the graph.
   1.541 +    ///
   1.542 +    /// \pre \ref init() must be called before using this function.
   1.543 +    bool findMinMean() {
   1.544 +      // Finding the minimum mean cycle in the components
   1.545 +      for (int comp = 0; comp < _component_num; ++comp) {
   1.546 +        if (!initCurrentComponent(comp)) continue;
   1.547 +        while (true) {
   1.548 +          if (!findPolicyCycles(comp)) return false;
   1.549 +          contractPolicyGraph(comp);
   1.550 +          if (!computeNodeDistances(comp)) return true;
   1.551 +        }
   1.552 +      }
   1.553 +    }
   1.554 +
   1.555 +    /// \brief Finds a critical (minimum mean) directed cycle.
   1.556 +    ///
   1.557 +    /// Finds a critical (minimum mean) directed cycle using the data
   1.558 +    /// computed in the \ref findMinMean() function.
   1.559 +    ///
   1.560 +    /// \return Returns \c true if a directed cycle exists in the graph.
   1.561 +    ///
   1.562 +    /// \pre \ref init() and \ref findMinMean() must be called before
   1.563 +    /// using this function.
   1.564 +    bool findCycle() {
   1.565 +      if (!_cycle_found) return false;
   1.566 +      _cycle_path->addBack(_policy[_cycle_node]);
   1.567 +      for ( Node v = _cycle_node;
   1.568 +            (v = _graph.target(_policy[v])) != _cycle_node; ) {
   1.569 +        _cycle_path->addBack(_policy[v]);
   1.570 +      }
   1.571 +      return true;
   1.572      }
   1.573  
   1.574      /// \brief Returns the total length of the found cycle.
   1.575      ///
   1.576      /// Returns the total length of the found cycle.
   1.577      ///
   1.578 -    /// \pre \ref run() or \ref findCycle() must be called before
   1.579 -    /// using this function. If only \ref findMinMean() is called,
   1.580 -    /// the return value is not valid.
   1.581 +    /// \pre \ref run() or \ref findMinMean() must be called before
   1.582 +    /// using this function.
   1.583      Length cycleLength() const {
   1.584 -      return cycle_length;
   1.585 +      return _cycle_length;
   1.586      }
   1.587  
   1.588 -    /// \brief Returns the number of edges in the found cycle.
   1.589 +    /// \brief Returns the number of edges on the found cycle.
   1.590      ///
   1.591 -    /// Returns the number of edges in the found cycle.
   1.592 +    /// Returns the number of edges on the found cycle.
   1.593      ///
   1.594 -    /// \pre \ref run() or \ref findCycle() must be called before
   1.595 -    /// using this function. If only \ref findMinMean() is called,
   1.596 -    /// the return value is not valid.
   1.597 +    /// \pre \ref run() or \ref findMinMean() must be called before
   1.598 +    /// using this function.
   1.599      int cycleEdgeNum() const {
   1.600 -      return cycle_size;
   1.601 +      return _cycle_size;
   1.602      }
   1.603  
   1.604      /// \brief Returns the mean length of the found cycle.
   1.605 @@ -384,52 +374,53 @@
   1.606      /// \pre \ref run() or \ref findMinMean() must be called before
   1.607      /// using this function.
   1.608      ///
   1.609 -    /// \warning LengthMap::Value must be convertible to double.
   1.610 -    ///
   1.611 -    /// \note <tt>m.minMean()</tt> is just a shortcut of the following
   1.612 -    /// code.
   1.613 +    /// \note <tt>mmc.cycleMean()</tt> is just a shortcut of the
   1.614 +    /// following code.
   1.615      /// \code
   1.616 -    ///   return m.cycleLength() / double(m.cycleEdgeNum());
   1.617 +    ///   return double(mmc.cycleLength()) / mmc.cycleEdgeNum();
   1.618      /// \endcode
   1.619 -    /// However if only \ref findMinMean() is called before using this
   1.620 -    /// function, <tt>m.cycleLength()</tt> and <tt>m.cycleEdgeNum()</tt>
   1.621 -    /// are not valid alone, only their ratio is valid.
   1.622 -    double minMean() const {
   1.623 -      return cycle_length / (double)cycle_size;
   1.624 +    double cycleMean() const {
   1.625 +      return (double)_cycle_length / _cycle_size;
   1.626      }
   1.627  
   1.628 -    /// \brief Returns a const reference to the \ref lemon::Path "Path"
   1.629 -    /// structure of the found cycle.
   1.630 +    /// \brief Returns a const reference to the \ref Path "path"
   1.631 +    /// structure storing the found cycle.
   1.632      ///
   1.633 -    /// Returns a const reference to the \ref lemon::Path "Path"
   1.634 -    /// structure of the found cycle.
   1.635 +    /// Returns a const reference to the \ref Path "path"
   1.636 +    /// structure storing the found cycle.
   1.637      ///
   1.638 -    /// \pre \ref run() must be called before using this function.
   1.639 +    /// \pre \ref run() or \ref findCycle() must be called before using
   1.640 +    /// this function.
   1.641      ///
   1.642 -    /// \sa \ref cyclePath()
   1.643 +    /// \sa cyclePath()
   1.644      const Path& cycle() const {
   1.645 -      return *cycle_path;
   1.646 +      return *_cycle_path;
   1.647      }
   1.648  
   1.649 -    /// \brief Sets the \ref lemon::Path "Path" structure storing the
   1.650 +    /// \brief Sets the \ref Path "path" structure for storing the found
   1.651 +    /// cycle.
   1.652 +    ///
   1.653 +    /// Sets an external \ref Path "path" structure for storing the
   1.654      /// found cycle.
   1.655      ///
   1.656 -    /// Sets the \ref lemon::Path "Path" structure storing the found
   1.657 -    /// cycle. If you don't use this function before calling
   1.658 -    /// \ref run(), it will allocate one. The destuctor deallocates
   1.659 -    /// this automatically allocated map, of course.
   1.660 +    /// If you don't call this function before calling \ref run() or
   1.661 +    /// \ref init(), it will allocate a local \ref Path "path"
   1.662 +    /// structure.
   1.663 +    /// The destuctor deallocates this automatically allocated map,
   1.664 +    /// of course.
   1.665      ///
   1.666 -    /// \note The algorithm calls only the \ref lemon::Path::addFront()
   1.667 -    /// "addFront()" method of the given \ref lemon::Path "Path"
   1.668 -    /// structure.
   1.669 +    /// \note The algorithm calls only the \ref lemon::Path::addBack()
   1.670 +    /// "addBack()" function of the given \ref Path "path" structure.
   1.671      ///
   1.672 -    /// \return \c (*this)
   1.673 +    /// \return <tt>(*this)</tt>
   1.674 +    ///
   1.675 +    /// \sa cycle()
   1.676      MinMeanCycle& cyclePath(Path &path) {
   1.677 -      if (local_path) {
   1.678 -	delete cycle_path;
   1.679 -	local_path = false;
   1.680 +      if (_local_path) {
   1.681 +        delete _cycle_path;
   1.682 +        _local_path = false;
   1.683        }
   1.684 -      cycle_path = &path;
   1.685 +      _cycle_path = &path;
   1.686        return *this;
   1.687      }
   1.688