Rename min mean cycle classes and their members (#179)
authorPeter Kovacs <kpeter@inf.elte.hu>
Sat, 13 Mar 2010 22:01:38 +0100
changeset 864d3ea191c3412
parent 863 a93f1a27d831
child 865 d48d79b11f5b
Rename min mean cycle classes and their members (#179)
with respect to the possible introduction of min ratio
cycle algorithms in the future.

The renamed classes:
- Karp --> KarpMmc
- HartmannOrlin --> HartmannOrlinMmc
- Howard --> HowardMmc

The renamed members:
- cycleLength() --> cycleCost()
- cycleArcNum() --> cycleSize()
- findMinMean() --> findCycleMean()
- Value --> Cost
- LargeValue --> LargeCost
- SetLargeValue --> SetLargeCost
lemon/Makefile.am
lemon/cycle_canceling.h
lemon/hartmann_orlin.h
lemon/hartmann_orlin_mmc.h
lemon/howard.h
lemon/howard_mmc.h
lemon/karp.h
lemon/karp_mmc.h
test/min_mean_cycle_test.cc
     1.1 --- a/lemon/Makefile.am	Mon Mar 08 08:33:41 2010 +0100
     1.2 +++ b/lemon/Makefile.am	Sat Mar 13 22:01:38 2010 +0100
     1.3 @@ -89,10 +89,10 @@
     1.4  	lemon/gomory_hu.h \
     1.5  	lemon/graph_to_eps.h \
     1.6  	lemon/grid_graph.h \
     1.7 -	lemon/hartmann_orlin.h \
     1.8 -	lemon/howard.h \
     1.9 +	lemon/hartmann_orlin_mmc.h \
    1.10 +	lemon/howard_mmc.h \
    1.11  	lemon/hypercube_graph.h \
    1.12 -	lemon/karp.h \
    1.13 +	lemon/karp_mmc.h \
    1.14  	lemon/kruskal.h \
    1.15  	lemon/hao_orlin.h \
    1.16  	lemon/lgf_reader.h \
     2.1 --- a/lemon/cycle_canceling.h	Mon Mar 08 08:33:41 2010 +0100
     2.2 +++ b/lemon/cycle_canceling.h	Sat Mar 13 22:01:38 2010 +0100
     2.3 @@ -34,7 +34,7 @@
     2.4  #include <lemon/adaptors.h>
     2.5  #include <lemon/circulation.h>
     2.6  #include <lemon/bellman_ford.h>
     2.7 -#include <lemon/howard.h>
     2.8 +#include <lemon/howard_mmc.h>
     2.9  
    2.10  namespace lemon {
    2.11  
    2.12 @@ -924,14 +924,14 @@
    2.13      void startMinMeanCycleCanceling() {
    2.14        typedef SimplePath<StaticDigraph> SPath;
    2.15        typedef typename SPath::ArcIt SPathArcIt;
    2.16 -      typedef typename Howard<StaticDigraph, CostArcMap>
    2.17 +      typedef typename HowardMmc<StaticDigraph, CostArcMap>
    2.18          ::template SetPath<SPath>::Create MMC;
    2.19        
    2.20        SPath cycle;
    2.21        MMC mmc(_sgr, _cost_map);
    2.22        mmc.cycle(cycle);
    2.23        buildResidualNetwork();
    2.24 -      while (mmc.findMinMean() && mmc.cycleLength() < 0) {
    2.25 +      while (mmc.findCycleMean() && mmc.cycleCost() < 0) {
    2.26          // Find the cycle
    2.27          mmc.findCycle();
    2.28  
    2.29 @@ -1132,17 +1132,17 @@
    2.30              }
    2.31            }
    2.32          } else {
    2.33 -          typedef Howard<StaticDigraph, CostArcMap> MMC;
    2.34 +          typedef HowardMmc<StaticDigraph, CostArcMap> MMC;
    2.35            typedef typename BellmanFord<StaticDigraph, CostArcMap>
    2.36              ::template SetDistMap<CostNodeMap>::Create BF;
    2.37  
    2.38            // Set epsilon to the minimum cycle mean
    2.39            buildResidualNetwork();
    2.40            MMC mmc(_sgr, _cost_map);
    2.41 -          mmc.findMinMean();
    2.42 +          mmc.findCycleMean();
    2.43            epsilon = -mmc.cycleMean();
    2.44 -          Cost cycle_cost = mmc.cycleLength();
    2.45 -          int cycle_size = mmc.cycleArcNum();
    2.46 +          Cost cycle_cost = mmc.cycleCost();
    2.47 +          int cycle_size = mmc.cycleSize();
    2.48            
    2.49            // Compute feasible potentials for the current epsilon
    2.50            for (int i = 0; i != int(_cost_vec.size()); ++i) {
     3.1 --- a/lemon/hartmann_orlin.h	Mon Mar 08 08:33:41 2010 +0100
     3.2 +++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
     3.3 @@ -1,650 +0,0 @@
     3.4 -/* -*- C++ -*-
     3.5 - *
     3.6 - * This file is a part of LEMON, a generic C++ optimization library
     3.7 - *
     3.8 - * Copyright (C) 2003-2008
     3.9 - * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    3.10 - * (Egervary Research Group on Combinatorial Optimization, EGRES).
    3.11 - *
    3.12 - * Permission to use, modify and distribute this software is granted
    3.13 - * provided that this copyright notice appears in all copies. For
    3.14 - * precise terms see the accompanying LICENSE file.
    3.15 - *
    3.16 - * This software is provided "AS IS" with no warranty of any kind,
    3.17 - * express or implied, and with no claim as to its suitability for any
    3.18 - * purpose.
    3.19 - *
    3.20 - */
    3.21 -
    3.22 -#ifndef LEMON_HARTMANN_ORLIN_H
    3.23 -#define LEMON_HARTMANN_ORLIN_H
    3.24 -
    3.25 -/// \ingroup min_mean_cycle
    3.26 -///
    3.27 -/// \file
    3.28 -/// \brief Hartmann-Orlin's algorithm for finding a minimum mean cycle.
    3.29 -
    3.30 -#include <vector>
    3.31 -#include <limits>
    3.32 -#include <lemon/core.h>
    3.33 -#include <lemon/path.h>
    3.34 -#include <lemon/tolerance.h>
    3.35 -#include <lemon/connectivity.h>
    3.36 -
    3.37 -namespace lemon {
    3.38 -
    3.39 -  /// \brief Default traits class of HartmannOrlin algorithm.
    3.40 -  ///
    3.41 -  /// Default traits class of HartmannOrlin algorithm.
    3.42 -  /// \tparam GR The type of the digraph.
    3.43 -  /// \tparam LEN The type of the length map.
    3.44 -  /// It must conform to the \ref concepts::Rea_data "Rea_data" concept.
    3.45 -#ifdef DOXYGEN
    3.46 -  template <typename GR, typename LEN>
    3.47 -#else
    3.48 -  template <typename GR, typename LEN,
    3.49 -    bool integer = std::numeric_limits<typename LEN::Value>::is_integer>
    3.50 -#endif
    3.51 -  struct HartmannOrlinDefaultTraits
    3.52 -  {
    3.53 -    /// The type of the digraph
    3.54 -    typedef GR Digraph;
    3.55 -    /// The type of the length map
    3.56 -    typedef LEN LengthMap;
    3.57 -    /// The type of the arc lengths
    3.58 -    typedef typename LengthMap::Value Value;
    3.59 -
    3.60 -    /// \brief The large value type used for internal computations
    3.61 -    ///
    3.62 -    /// The large value type used for internal computations.
    3.63 -    /// It is \c long \c long if the \c Value type is integer,
    3.64 -    /// otherwise it is \c double.
    3.65 -    /// \c Value must be convertible to \c LargeValue.
    3.66 -    typedef double LargeValue;
    3.67 -
    3.68 -    /// The tolerance type used for internal computations
    3.69 -    typedef lemon::Tolerance<LargeValue> Tolerance;
    3.70 -
    3.71 -    /// \brief The path type of the found cycles
    3.72 -    ///
    3.73 -    /// The path type of the found cycles.
    3.74 -    /// It must conform to the \ref lemon::concepts::Path "Path" concept
    3.75 -    /// and it must have an \c addFront() function.
    3.76 -    typedef lemon::Path<Digraph> Path;
    3.77 -  };
    3.78 -
    3.79 -  // Default traits class for integer value types
    3.80 -  template <typename GR, typename LEN>
    3.81 -  struct HartmannOrlinDefaultTraits<GR, LEN, true>
    3.82 -  {
    3.83 -    typedef GR Digraph;
    3.84 -    typedef LEN LengthMap;
    3.85 -    typedef typename LengthMap::Value Value;
    3.86 -#ifdef LEMON_HAVE_LONG_LONG
    3.87 -    typedef long long LargeValue;
    3.88 -#else
    3.89 -    typedef long LargeValue;
    3.90 -#endif
    3.91 -    typedef lemon::Tolerance<LargeValue> Tolerance;
    3.92 -    typedef lemon::Path<Digraph> Path;
    3.93 -  };
    3.94 -
    3.95 -
    3.96 -  /// \addtogroup min_mean_cycle
    3.97 -  /// @{
    3.98 -
    3.99 -  /// \brief Implementation of the Hartmann-Orlin algorithm for finding
   3.100 -  /// a minimum mean cycle.
   3.101 -  ///
   3.102 -  /// This class implements the Hartmann-Orlin algorithm for finding
   3.103 -  /// a directed cycle of minimum mean length (cost) in a digraph
   3.104 -  /// \ref amo93networkflows, \ref dasdan98minmeancycle.
   3.105 -  /// It is an improved version of \ref Karp "Karp"'s original algorithm,
   3.106 -  /// it applies an efficient early termination scheme.
   3.107 -  /// It runs in time O(ne) and uses space O(n<sup>2</sup>+e).
   3.108 -  ///
   3.109 -  /// \tparam GR The type of the digraph the algorithm runs on.
   3.110 -  /// \tparam LEN The type of the length map. The default
   3.111 -  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
   3.112 -  /// \tparam TR The traits class that defines various types used by the
   3.113 -  /// algorithm. By default, it is \ref HartmannOrlinDefaultTraits
   3.114 -  /// "HartmannOrlinDefaultTraits<GR, LEN>".
   3.115 -  /// In most cases, this parameter should not be set directly,
   3.116 -  /// consider to use the named template parameters instead.
   3.117 -#ifdef DOXYGEN
   3.118 -  template <typename GR, typename LEN, typename TR>
   3.119 -#else
   3.120 -  template < typename GR,
   3.121 -             typename LEN = typename GR::template ArcMap<int>,
   3.122 -             typename TR = HartmannOrlinDefaultTraits<GR, LEN> >
   3.123 -#endif
   3.124 -  class HartmannOrlin
   3.125 -  {
   3.126 -  public:
   3.127 -
   3.128 -    /// The type of the digraph
   3.129 -    typedef typename TR::Digraph Digraph;
   3.130 -    /// The type of the length map
   3.131 -    typedef typename TR::LengthMap LengthMap;
   3.132 -    /// The type of the arc lengths
   3.133 -    typedef typename TR::Value Value;
   3.134 -
   3.135 -    /// \brief The large value type
   3.136 -    ///
   3.137 -    /// The large value type used for internal computations.
   3.138 -    /// By default, it is \c long \c long if the \c Value type is integer,
   3.139 -    /// otherwise it is \c double.
   3.140 -    typedef typename TR::LargeValue LargeValue;
   3.141 -
   3.142 -    /// The tolerance type
   3.143 -    typedef typename TR::Tolerance Tolerance;
   3.144 -
   3.145 -    /// \brief The path type of the found cycles
   3.146 -    ///
   3.147 -    /// The path type of the found cycles.
   3.148 -    /// Using the \ref HartmannOrlinDefaultTraits "default traits class",
   3.149 -    /// it is \ref lemon::Path "Path<Digraph>".
   3.150 -    typedef typename TR::Path Path;
   3.151 -
   3.152 -    /// The \ref HartmannOrlinDefaultTraits "traits class" of the algorithm
   3.153 -    typedef TR Traits;
   3.154 -
   3.155 -  private:
   3.156 -
   3.157 -    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   3.158 -
   3.159 -    // Data sturcture for path data
   3.160 -    struct PathData
   3.161 -    {
   3.162 -      LargeValue dist;
   3.163 -      Arc pred;
   3.164 -      PathData(LargeValue d, Arc p = INVALID) :
   3.165 -        dist(d), pred(p) {}
   3.166 -    };
   3.167 -
   3.168 -    typedef typename Digraph::template NodeMap<std::vector<PathData> >
   3.169 -      PathDataNodeMap;
   3.170 -
   3.171 -  private:
   3.172 -
   3.173 -    // The digraph the algorithm runs on
   3.174 -    const Digraph &_gr;
   3.175 -    // The length of the arcs
   3.176 -    const LengthMap &_length;
   3.177 -
   3.178 -    // Data for storing the strongly connected components
   3.179 -    int _comp_num;
   3.180 -    typename Digraph::template NodeMap<int> _comp;
   3.181 -    std::vector<std::vector<Node> > _comp_nodes;
   3.182 -    std::vector<Node>* _nodes;
   3.183 -    typename Digraph::template NodeMap<std::vector<Arc> > _out_arcs;
   3.184 -
   3.185 -    // Data for the found cycles
   3.186 -    bool _curr_found, _best_found;
   3.187 -    LargeValue _curr_length, _best_length;
   3.188 -    int _curr_size, _best_size;
   3.189 -    Node _curr_node, _best_node;
   3.190 -    int _curr_level, _best_level;
   3.191 -
   3.192 -    Path *_cycle_path;
   3.193 -    bool _local_path;
   3.194 -
   3.195 -    // Node map for storing path data
   3.196 -    PathDataNodeMap _data;
   3.197 -    // The processed nodes in the last round
   3.198 -    std::vector<Node> _process;
   3.199 -
   3.200 -    Tolerance _tolerance;
   3.201 -
   3.202 -    // Infinite constant
   3.203 -    const LargeValue INF;
   3.204 -
   3.205 -  public:
   3.206 -
   3.207 -    /// \name Named Template Parameters
   3.208 -    /// @{
   3.209 -
   3.210 -    template <typename T>
   3.211 -    struct SetLargeValueTraits : public Traits {
   3.212 -      typedef T LargeValue;
   3.213 -      typedef lemon::Tolerance<T> Tolerance;
   3.214 -    };
   3.215 -
   3.216 -    /// \brief \ref named-templ-param "Named parameter" for setting
   3.217 -    /// \c LargeValue type.
   3.218 -    ///
   3.219 -    /// \ref named-templ-param "Named parameter" for setting \c LargeValue
   3.220 -    /// type. It is used for internal computations in the algorithm.
   3.221 -    template <typename T>
   3.222 -    struct SetLargeValue
   3.223 -      : public HartmannOrlin<GR, LEN, SetLargeValueTraits<T> > {
   3.224 -      typedef HartmannOrlin<GR, LEN, SetLargeValueTraits<T> > Create;
   3.225 -    };
   3.226 -
   3.227 -    template <typename T>
   3.228 -    struct SetPathTraits : public Traits {
   3.229 -      typedef T Path;
   3.230 -    };
   3.231 -
   3.232 -    /// \brief \ref named-templ-param "Named parameter" for setting
   3.233 -    /// \c %Path type.
   3.234 -    ///
   3.235 -    /// \ref named-templ-param "Named parameter" for setting the \c %Path
   3.236 -    /// type of the found cycles.
   3.237 -    /// It must conform to the \ref lemon::concepts::Path "Path" concept
   3.238 -    /// and it must have an \c addFront() function.
   3.239 -    template <typename T>
   3.240 -    struct SetPath
   3.241 -      : public HartmannOrlin<GR, LEN, SetPathTraits<T> > {
   3.242 -      typedef HartmannOrlin<GR, LEN, SetPathTraits<T> > Create;
   3.243 -    };
   3.244 -
   3.245 -    /// @}
   3.246 -
   3.247 -  protected:
   3.248 -
   3.249 -    HartmannOrlin() {}
   3.250 -
   3.251 -  public:
   3.252 -
   3.253 -    /// \brief Constructor.
   3.254 -    ///
   3.255 -    /// The constructor of the class.
   3.256 -    ///
   3.257 -    /// \param digraph The digraph the algorithm runs on.
   3.258 -    /// \param length The lengths (costs) of the arcs.
   3.259 -    HartmannOrlin( const Digraph &digraph,
   3.260 -                   const LengthMap &length ) :
   3.261 -      _gr(digraph), _length(length), _comp(digraph), _out_arcs(digraph),
   3.262 -      _best_found(false), _best_length(0), _best_size(1),
   3.263 -      _cycle_path(NULL), _local_path(false), _data(digraph),
   3.264 -      INF(std::numeric_limits<LargeValue>::has_infinity ?
   3.265 -          std::numeric_limits<LargeValue>::infinity() :
   3.266 -          std::numeric_limits<LargeValue>::max())
   3.267 -    {}
   3.268 -
   3.269 -    /// Destructor.
   3.270 -    ~HartmannOrlin() {
   3.271 -      if (_local_path) delete _cycle_path;
   3.272 -    }
   3.273 -
   3.274 -    /// \brief Set the path structure for storing the found cycle.
   3.275 -    ///
   3.276 -    /// This function sets an external path structure for storing the
   3.277 -    /// found cycle.
   3.278 -    ///
   3.279 -    /// If you don't call this function before calling \ref run() or
   3.280 -    /// \ref findMinMean(), it will allocate a local \ref Path "path"
   3.281 -    /// structure. The destuctor deallocates this automatically
   3.282 -    /// allocated object, of course.
   3.283 -    ///
   3.284 -    /// \note The algorithm calls only the \ref lemon::Path::addFront()
   3.285 -    /// "addFront()" function of the given path structure.
   3.286 -    ///
   3.287 -    /// \return <tt>(*this)</tt>
   3.288 -    HartmannOrlin& cycle(Path &path) {
   3.289 -      if (_local_path) {
   3.290 -        delete _cycle_path;
   3.291 -        _local_path = false;
   3.292 -      }
   3.293 -      _cycle_path = &path;
   3.294 -      return *this;
   3.295 -    }
   3.296 -
   3.297 -    /// \brief Set the tolerance used by the algorithm.
   3.298 -    ///
   3.299 -    /// This function sets the tolerance object used by the algorithm.
   3.300 -    ///
   3.301 -    /// \return <tt>(*this)</tt>
   3.302 -    HartmannOrlin& tolerance(const Tolerance& tolerance) {
   3.303 -      _tolerance = tolerance;
   3.304 -      return *this;
   3.305 -    }
   3.306 -
   3.307 -    /// \brief Return a const reference to the tolerance.
   3.308 -    ///
   3.309 -    /// This function returns a const reference to the tolerance object
   3.310 -    /// used by the algorithm.
   3.311 -    const Tolerance& tolerance() const {
   3.312 -      return _tolerance;
   3.313 -    }
   3.314 -
   3.315 -    /// \name Execution control
   3.316 -    /// The simplest way to execute the algorithm is to call the \ref run()
   3.317 -    /// function.\n
   3.318 -    /// If you only need the minimum mean length, you may call
   3.319 -    /// \ref findMinMean().
   3.320 -
   3.321 -    /// @{
   3.322 -
   3.323 -    /// \brief Run the algorithm.
   3.324 -    ///
   3.325 -    /// This function runs the algorithm.
   3.326 -    /// It can be called more than once (e.g. if the underlying digraph
   3.327 -    /// and/or the arc lengths have been modified).
   3.328 -    ///
   3.329 -    /// \return \c true if a directed cycle exists in the digraph.
   3.330 -    ///
   3.331 -    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.
   3.332 -    /// \code
   3.333 -    ///   return mmc.findMinMean() && mmc.findCycle();
   3.334 -    /// \endcode
   3.335 -    bool run() {
   3.336 -      return findMinMean() && findCycle();
   3.337 -    }
   3.338 -
   3.339 -    /// \brief Find the minimum cycle mean.
   3.340 -    ///
   3.341 -    /// This function finds the minimum mean length of the directed
   3.342 -    /// cycles in the digraph.
   3.343 -    ///
   3.344 -    /// \return \c true if a directed cycle exists in the digraph.
   3.345 -    bool findMinMean() {
   3.346 -      // Initialization and find strongly connected components
   3.347 -      init();
   3.348 -      findComponents();
   3.349 -      
   3.350 -      // Find the minimum cycle mean in the components
   3.351 -      for (int comp = 0; comp < _comp_num; ++comp) {
   3.352 -        if (!initComponent(comp)) continue;
   3.353 -        processRounds();
   3.354 -        
   3.355 -        // Update the best cycle (global minimum mean cycle)
   3.356 -        if ( _curr_found && (!_best_found || 
   3.357 -             _curr_length * _best_size < _best_length * _curr_size) ) {
   3.358 -          _best_found = true;
   3.359 -          _best_length = _curr_length;
   3.360 -          _best_size = _curr_size;
   3.361 -          _best_node = _curr_node;
   3.362 -          _best_level = _curr_level;
   3.363 -        }
   3.364 -      }
   3.365 -      return _best_found;
   3.366 -    }
   3.367 -
   3.368 -    /// \brief Find a minimum mean directed cycle.
   3.369 -    ///
   3.370 -    /// This function finds a directed cycle of minimum mean length
   3.371 -    /// in the digraph using the data computed by findMinMean().
   3.372 -    ///
   3.373 -    /// \return \c true if a directed cycle exists in the digraph.
   3.374 -    ///
   3.375 -    /// \pre \ref findMinMean() must be called before using this function.
   3.376 -    bool findCycle() {
   3.377 -      if (!_best_found) return false;
   3.378 -      IntNodeMap reached(_gr, -1);
   3.379 -      int r = _best_level + 1;
   3.380 -      Node u = _best_node;
   3.381 -      while (reached[u] < 0) {
   3.382 -        reached[u] = --r;
   3.383 -        u = _gr.source(_data[u][r].pred);
   3.384 -      }
   3.385 -      r = reached[u];
   3.386 -      Arc e = _data[u][r].pred;
   3.387 -      _cycle_path->addFront(e);
   3.388 -      _best_length = _length[e];
   3.389 -      _best_size = 1;
   3.390 -      Node v;
   3.391 -      while ((v = _gr.source(e)) != u) {
   3.392 -        e = _data[v][--r].pred;
   3.393 -        _cycle_path->addFront(e);
   3.394 -        _best_length += _length[e];
   3.395 -        ++_best_size;
   3.396 -      }
   3.397 -      return true;
   3.398 -    }
   3.399 -
   3.400 -    /// @}
   3.401 -
   3.402 -    /// \name Query Functions
   3.403 -    /// The results of the algorithm can be obtained using these
   3.404 -    /// functions.\n
   3.405 -    /// The algorithm should be executed before using them.
   3.406 -
   3.407 -    /// @{
   3.408 -
   3.409 -    /// \brief Return the total length of the found cycle.
   3.410 -    ///
   3.411 -    /// This function returns the total length of the found cycle.
   3.412 -    ///
   3.413 -    /// \pre \ref run() or \ref findMinMean() must be called before
   3.414 -    /// using this function.
   3.415 -    Value cycleLength() const {
   3.416 -      return static_cast<Value>(_best_length);
   3.417 -    }
   3.418 -
   3.419 -    /// \brief Return the number of arcs on the found cycle.
   3.420 -    ///
   3.421 -    /// This function returns the number of arcs on the found cycle.
   3.422 -    ///
   3.423 -    /// \pre \ref run() or \ref findMinMean() must be called before
   3.424 -    /// using this function.
   3.425 -    int cycleArcNum() const {
   3.426 -      return _best_size;
   3.427 -    }
   3.428 -
   3.429 -    /// \brief Return the mean length of the found cycle.
   3.430 -    ///
   3.431 -    /// This function returns the mean length of the found cycle.
   3.432 -    ///
   3.433 -    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the
   3.434 -    /// following code.
   3.435 -    /// \code
   3.436 -    ///   return static_cast<double>(alg.cycleLength()) / alg.cycleArcNum();
   3.437 -    /// \endcode
   3.438 -    ///
   3.439 -    /// \pre \ref run() or \ref findMinMean() must be called before
   3.440 -    /// using this function.
   3.441 -    double cycleMean() const {
   3.442 -      return static_cast<double>(_best_length) / _best_size;
   3.443 -    }
   3.444 -
   3.445 -    /// \brief Return the found cycle.
   3.446 -    ///
   3.447 -    /// This function returns a const reference to the path structure
   3.448 -    /// storing the found cycle.
   3.449 -    ///
   3.450 -    /// \pre \ref run() or \ref findCycle() must be called before using
   3.451 -    /// this function.
   3.452 -    const Path& cycle() const {
   3.453 -      return *_cycle_path;
   3.454 -    }
   3.455 -
   3.456 -    ///@}
   3.457 -
   3.458 -  private:
   3.459 -
   3.460 -    // Initialization
   3.461 -    void init() {
   3.462 -      if (!_cycle_path) {
   3.463 -        _local_path = true;
   3.464 -        _cycle_path = new Path;
   3.465 -      }
   3.466 -      _cycle_path->clear();
   3.467 -      _best_found = false;
   3.468 -      _best_length = 0;
   3.469 -      _best_size = 1;
   3.470 -      _cycle_path->clear();
   3.471 -      for (NodeIt u(_gr); u != INVALID; ++u)
   3.472 -        _data[u].clear();
   3.473 -    }
   3.474 -
   3.475 -    // Find strongly connected components and initialize _comp_nodes
   3.476 -    // and _out_arcs
   3.477 -    void findComponents() {
   3.478 -      _comp_num = stronglyConnectedComponents(_gr, _comp);
   3.479 -      _comp_nodes.resize(_comp_num);
   3.480 -      if (_comp_num == 1) {
   3.481 -        _comp_nodes[0].clear();
   3.482 -        for (NodeIt n(_gr); n != INVALID; ++n) {
   3.483 -          _comp_nodes[0].push_back(n);
   3.484 -          _out_arcs[n].clear();
   3.485 -          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   3.486 -            _out_arcs[n].push_back(a);
   3.487 -          }
   3.488 -        }
   3.489 -      } else {
   3.490 -        for (int i = 0; i < _comp_num; ++i)
   3.491 -          _comp_nodes[i].clear();
   3.492 -        for (NodeIt n(_gr); n != INVALID; ++n) {
   3.493 -          int k = _comp[n];
   3.494 -          _comp_nodes[k].push_back(n);
   3.495 -          _out_arcs[n].clear();
   3.496 -          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   3.497 -            if (_comp[_gr.target(a)] == k) _out_arcs[n].push_back(a);
   3.498 -          }
   3.499 -        }
   3.500 -      }
   3.501 -    }
   3.502 -
   3.503 -    // Initialize path data for the current component
   3.504 -    bool initComponent(int comp) {
   3.505 -      _nodes = &(_comp_nodes[comp]);
   3.506 -      int n = _nodes->size();
   3.507 -      if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) {
   3.508 -        return false;
   3.509 -      }      
   3.510 -      for (int i = 0; i < n; ++i) {
   3.511 -        _data[(*_nodes)[i]].resize(n + 1, PathData(INF));
   3.512 -      }
   3.513 -      return true;
   3.514 -    }
   3.515 -
   3.516 -    // Process all rounds of computing path data for the current component.
   3.517 -    // _data[v][k] is the length of a shortest directed walk from the root
   3.518 -    // node to node v containing exactly k arcs.
   3.519 -    void processRounds() {
   3.520 -      Node start = (*_nodes)[0];
   3.521 -      _data[start][0] = PathData(0);
   3.522 -      _process.clear();
   3.523 -      _process.push_back(start);
   3.524 -
   3.525 -      int k, n = _nodes->size();
   3.526 -      int next_check = 4;
   3.527 -      bool terminate = false;
   3.528 -      for (k = 1; k <= n && int(_process.size()) < n && !terminate; ++k) {
   3.529 -        processNextBuildRound(k);
   3.530 -        if (k == next_check || k == n) {
   3.531 -          terminate = checkTermination(k);
   3.532 -          next_check = next_check * 3 / 2;
   3.533 -        }
   3.534 -      }
   3.535 -      for ( ; k <= n && !terminate; ++k) {
   3.536 -        processNextFullRound(k);
   3.537 -        if (k == next_check || k == n) {
   3.538 -          terminate = checkTermination(k);
   3.539 -          next_check = next_check * 3 / 2;
   3.540 -        }
   3.541 -      }
   3.542 -    }
   3.543 -
   3.544 -    // Process one round and rebuild _process
   3.545 -    void processNextBuildRound(int k) {
   3.546 -      std::vector<Node> next;
   3.547 -      Node u, v;
   3.548 -      Arc e;
   3.549 -      LargeValue d;
   3.550 -      for (int i = 0; i < int(_process.size()); ++i) {
   3.551 -        u = _process[i];
   3.552 -        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   3.553 -          e = _out_arcs[u][j];
   3.554 -          v = _gr.target(e);
   3.555 -          d = _data[u][k-1].dist + _length[e];
   3.556 -          if (_tolerance.less(d, _data[v][k].dist)) {
   3.557 -            if (_data[v][k].dist == INF) next.push_back(v);
   3.558 -            _data[v][k] = PathData(d, e);
   3.559 -          }
   3.560 -        }
   3.561 -      }
   3.562 -      _process.swap(next);
   3.563 -    }
   3.564 -
   3.565 -    // Process one round using _nodes instead of _process
   3.566 -    void processNextFullRound(int k) {
   3.567 -      Node u, v;
   3.568 -      Arc e;
   3.569 -      LargeValue d;
   3.570 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   3.571 -        u = (*_nodes)[i];
   3.572 -        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   3.573 -          e = _out_arcs[u][j];
   3.574 -          v = _gr.target(e);
   3.575 -          d = _data[u][k-1].dist + _length[e];
   3.576 -          if (_tolerance.less(d, _data[v][k].dist)) {
   3.577 -            _data[v][k] = PathData(d, e);
   3.578 -          }
   3.579 -        }
   3.580 -      }
   3.581 -    }
   3.582 -    
   3.583 -    // Check early termination
   3.584 -    bool checkTermination(int k) {
   3.585 -      typedef std::pair<int, int> Pair;
   3.586 -      typename GR::template NodeMap<Pair> level(_gr, Pair(-1, 0));
   3.587 -      typename GR::template NodeMap<LargeValue> pi(_gr);
   3.588 -      int n = _nodes->size();
   3.589 -      LargeValue length;
   3.590 -      int size;
   3.591 -      Node u;
   3.592 -      
   3.593 -      // Search for cycles that are already found
   3.594 -      _curr_found = false;
   3.595 -      for (int i = 0; i < n; ++i) {
   3.596 -        u = (*_nodes)[i];
   3.597 -        if (_data[u][k].dist == INF) continue;
   3.598 -        for (int j = k; j >= 0; --j) {
   3.599 -          if (level[u].first == i && level[u].second > 0) {
   3.600 -            // A cycle is found
   3.601 -            length = _data[u][level[u].second].dist - _data[u][j].dist;
   3.602 -            size = level[u].second - j;
   3.603 -            if (!_curr_found || length * _curr_size < _curr_length * size) {
   3.604 -              _curr_length = length;
   3.605 -              _curr_size = size;
   3.606 -              _curr_node = u;
   3.607 -              _curr_level = level[u].second;
   3.608 -              _curr_found = true;
   3.609 -            }
   3.610 -          }
   3.611 -          level[u] = Pair(i, j);
   3.612 -          if (j != 0) {
   3.613 -	    u = _gr.source(_data[u][j].pred);
   3.614 -	  }
   3.615 -        }
   3.616 -      }
   3.617 -
   3.618 -      // If at least one cycle is found, check the optimality condition
   3.619 -      LargeValue d;
   3.620 -      if (_curr_found && k < n) {
   3.621 -        // Find node potentials
   3.622 -        for (int i = 0; i < n; ++i) {
   3.623 -          u = (*_nodes)[i];
   3.624 -          pi[u] = INF;
   3.625 -          for (int j = 0; j <= k; ++j) {
   3.626 -            if (_data[u][j].dist < INF) {
   3.627 -              d = _data[u][j].dist * _curr_size - j * _curr_length;
   3.628 -              if (_tolerance.less(d, pi[u])) pi[u] = d;
   3.629 -            }
   3.630 -          }
   3.631 -        }
   3.632 -
   3.633 -        // Check the optimality condition for all arcs
   3.634 -        bool done = true;
   3.635 -        for (ArcIt a(_gr); a != INVALID; ++a) {
   3.636 -          if (_tolerance.less(_length[a] * _curr_size - _curr_length,
   3.637 -                              pi[_gr.target(a)] - pi[_gr.source(a)]) ) {
   3.638 -            done = false;
   3.639 -            break;
   3.640 -          }
   3.641 -        }
   3.642 -        return done;
   3.643 -      }
   3.644 -      return (k == n);
   3.645 -    }
   3.646 -
   3.647 -  }; //class HartmannOrlin
   3.648 -
   3.649 -  ///@}
   3.650 -
   3.651 -} //namespace lemon
   3.652 -
   3.653 -#endif //LEMON_HARTMANN_ORLIN_H
     4.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     4.2 +++ b/lemon/hartmann_orlin_mmc.h	Sat Mar 13 22:01:38 2010 +0100
     4.3 @@ -0,0 +1,650 @@
     4.4 +/* -*- C++ -*-
     4.5 + *
     4.6 + * This file is a part of LEMON, a generic C++ optimization library
     4.7 + *
     4.8 + * Copyright (C) 2003-2008
     4.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    4.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    4.11 + *
    4.12 + * Permission to use, modify and distribute this software is granted
    4.13 + * provided that this copyright notice appears in all copies. For
    4.14 + * precise terms see the accompanying LICENSE file.
    4.15 + *
    4.16 + * This software is provided "AS IS" with no warranty of any kind,
    4.17 + * express or implied, and with no claim as to its suitability for any
    4.18 + * purpose.
    4.19 + *
    4.20 + */
    4.21 +
    4.22 +#ifndef LEMON_HARTMANN_ORLIN_MMC_H
    4.23 +#define LEMON_HARTMANN_ORLIN_MMC_H
    4.24 +
    4.25 +/// \ingroup min_mean_cycle
    4.26 +///
    4.27 +/// \file
    4.28 +/// \brief Hartmann-Orlin's algorithm for finding a minimum mean cycle.
    4.29 +
    4.30 +#include <vector>
    4.31 +#include <limits>
    4.32 +#include <lemon/core.h>
    4.33 +#include <lemon/path.h>
    4.34 +#include <lemon/tolerance.h>
    4.35 +#include <lemon/connectivity.h>
    4.36 +
    4.37 +namespace lemon {
    4.38 +
    4.39 +  /// \brief Default traits class of HartmannOrlinMmc class.
    4.40 +  ///
    4.41 +  /// Default traits class of HartmannOrlinMmc class.
    4.42 +  /// \tparam GR The type of the digraph.
    4.43 +  /// \tparam CM The type of the cost map.
    4.44 +  /// It must conform to the \ref concepts::Rea_data "Rea_data" concept.
    4.45 +#ifdef DOXYGEN
    4.46 +  template <typename GR, typename CM>
    4.47 +#else
    4.48 +  template <typename GR, typename CM,
    4.49 +    bool integer = std::numeric_limits<typename CM::Value>::is_integer>
    4.50 +#endif
    4.51 +  struct HartmannOrlinMmcDefaultTraits
    4.52 +  {
    4.53 +    /// The type of the digraph
    4.54 +    typedef GR Digraph;
    4.55 +    /// The type of the cost map
    4.56 +    typedef CM CostMap;
    4.57 +    /// The type of the arc costs
    4.58 +    typedef typename CostMap::Value Cost;
    4.59 +
    4.60 +    /// \brief The large cost type used for internal computations
    4.61 +    ///
    4.62 +    /// The large cost type used for internal computations.
    4.63 +    /// It is \c long \c long if the \c Cost type is integer,
    4.64 +    /// otherwise it is \c double.
    4.65 +    /// \c Cost must be convertible to \c LargeCost.
    4.66 +    typedef double LargeCost;
    4.67 +
    4.68 +    /// The tolerance type used for internal computations
    4.69 +    typedef lemon::Tolerance<LargeCost> Tolerance;
    4.70 +
    4.71 +    /// \brief The path type of the found cycles
    4.72 +    ///
    4.73 +    /// The path type of the found cycles.
    4.74 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
    4.75 +    /// and it must have an \c addFront() function.
    4.76 +    typedef lemon::Path<Digraph> Path;
    4.77 +  };
    4.78 +
    4.79 +  // Default traits class for integer cost types
    4.80 +  template <typename GR, typename CM>
    4.81 +  struct HartmannOrlinMmcDefaultTraits<GR, CM, true>
    4.82 +  {
    4.83 +    typedef GR Digraph;
    4.84 +    typedef CM CostMap;
    4.85 +    typedef typename CostMap::Value Cost;
    4.86 +#ifdef LEMON_HAVE_LONG_LONG
    4.87 +    typedef long long LargeCost;
    4.88 +#else
    4.89 +    typedef long LargeCost;
    4.90 +#endif
    4.91 +    typedef lemon::Tolerance<LargeCost> Tolerance;
    4.92 +    typedef lemon::Path<Digraph> Path;
    4.93 +  };
    4.94 +
    4.95 +
    4.96 +  /// \addtogroup min_mean_cycle
    4.97 +  /// @{
    4.98 +
    4.99 +  /// \brief Implementation of the Hartmann-Orlin algorithm for finding
   4.100 +  /// a minimum mean cycle.
   4.101 +  ///
   4.102 +  /// This class implements the Hartmann-Orlin algorithm for finding
   4.103 +  /// a directed cycle of minimum mean cost in a digraph
   4.104 +  /// \ref amo93networkflows, \ref dasdan98minmeancycle.
   4.105 +  /// It is an improved version of \ref Karp "Karp"'s original algorithm,
   4.106 +  /// it applies an efficient early termination scheme.
   4.107 +  /// It runs in time O(ne) and uses space O(n<sup>2</sup>+e).
   4.108 +  ///
   4.109 +  /// \tparam GR The type of the digraph the algorithm runs on.
   4.110 +  /// \tparam CM The type of the cost map. The default
   4.111 +  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
   4.112 +  /// \tparam TR The traits class that defines various types used by the
   4.113 +  /// algorithm. By default, it is \ref HartmannOrlinMmcDefaultTraits
   4.114 +  /// "HartmannOrlinMmcDefaultTraits<GR, CM>".
   4.115 +  /// In most cases, this parameter should not be set directly,
   4.116 +  /// consider to use the named template parameters instead.
   4.117 +#ifdef DOXYGEN
   4.118 +  template <typename GR, typename CM, typename TR>
   4.119 +#else
   4.120 +  template < typename GR,
   4.121 +             typename CM = typename GR::template ArcMap<int>,
   4.122 +             typename TR = HartmannOrlinMmcDefaultTraits<GR, CM> >
   4.123 +#endif
   4.124 +  class HartmannOrlinMmc
   4.125 +  {
   4.126 +  public:
   4.127 +
   4.128 +    /// The type of the digraph
   4.129 +    typedef typename TR::Digraph Digraph;
   4.130 +    /// The type of the cost map
   4.131 +    typedef typename TR::CostMap CostMap;
   4.132 +    /// The type of the arc costs
   4.133 +    typedef typename TR::Cost Cost;
   4.134 +
   4.135 +    /// \brief The large cost type
   4.136 +    ///
   4.137 +    /// The large cost type used for internal computations.
   4.138 +    /// By default, it is \c long \c long if the \c Cost type is integer,
   4.139 +    /// otherwise it is \c double.
   4.140 +    typedef typename TR::LargeCost LargeCost;
   4.141 +
   4.142 +    /// The tolerance type
   4.143 +    typedef typename TR::Tolerance Tolerance;
   4.144 +
   4.145 +    /// \brief The path type of the found cycles
   4.146 +    ///
   4.147 +    /// The path type of the found cycles.
   4.148 +    /// Using the \ref HartmannOrlinMmcDefaultTraits "default traits class",
   4.149 +    /// it is \ref lemon::Path "Path<Digraph>".
   4.150 +    typedef typename TR::Path Path;
   4.151 +
   4.152 +    /// The \ref HartmannOrlinMmcDefaultTraits "traits class" of the algorithm
   4.153 +    typedef TR Traits;
   4.154 +
   4.155 +  private:
   4.156 +
   4.157 +    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   4.158 +
   4.159 +    // Data sturcture for path data
   4.160 +    struct PathData
   4.161 +    {
   4.162 +      LargeCost dist;
   4.163 +      Arc pred;
   4.164 +      PathData(LargeCost d, Arc p = INVALID) :
   4.165 +        dist(d), pred(p) {}
   4.166 +    };
   4.167 +
   4.168 +    typedef typename Digraph::template NodeMap<std::vector<PathData> >
   4.169 +      PathDataNodeMap;
   4.170 +
   4.171 +  private:
   4.172 +
   4.173 +    // The digraph the algorithm runs on
   4.174 +    const Digraph &_gr;
   4.175 +    // The cost of the arcs
   4.176 +    const CostMap &_cost;
   4.177 +
   4.178 +    // Data for storing the strongly connected components
   4.179 +    int _comp_num;
   4.180 +    typename Digraph::template NodeMap<int> _comp;
   4.181 +    std::vector<std::vector<Node> > _comp_nodes;
   4.182 +    std::vector<Node>* _nodes;
   4.183 +    typename Digraph::template NodeMap<std::vector<Arc> > _out_arcs;
   4.184 +
   4.185 +    // Data for the found cycles
   4.186 +    bool _curr_found, _best_found;
   4.187 +    LargeCost _curr_cost, _best_cost;
   4.188 +    int _curr_size, _best_size;
   4.189 +    Node _curr_node, _best_node;
   4.190 +    int _curr_level, _best_level;
   4.191 +
   4.192 +    Path *_cycle_path;
   4.193 +    bool _local_path;
   4.194 +
   4.195 +    // Node map for storing path data
   4.196 +    PathDataNodeMap _data;
   4.197 +    // The processed nodes in the last round
   4.198 +    std::vector<Node> _process;
   4.199 +
   4.200 +    Tolerance _tolerance;
   4.201 +
   4.202 +    // Infinite constant
   4.203 +    const LargeCost INF;
   4.204 +
   4.205 +  public:
   4.206 +
   4.207 +    /// \name Named Template Parameters
   4.208 +    /// @{
   4.209 +
   4.210 +    template <typename T>
   4.211 +    struct SetLargeCostTraits : public Traits {
   4.212 +      typedef T LargeCost;
   4.213 +      typedef lemon::Tolerance<T> Tolerance;
   4.214 +    };
   4.215 +
   4.216 +    /// \brief \ref named-templ-param "Named parameter" for setting
   4.217 +    /// \c LargeCost type.
   4.218 +    ///
   4.219 +    /// \ref named-templ-param "Named parameter" for setting \c LargeCost
   4.220 +    /// type. It is used for internal computations in the algorithm.
   4.221 +    template <typename T>
   4.222 +    struct SetLargeCost
   4.223 +      : public HartmannOrlinMmc<GR, CM, SetLargeCostTraits<T> > {
   4.224 +      typedef HartmannOrlinMmc<GR, CM, SetLargeCostTraits<T> > Create;
   4.225 +    };
   4.226 +
   4.227 +    template <typename T>
   4.228 +    struct SetPathTraits : public Traits {
   4.229 +      typedef T Path;
   4.230 +    };
   4.231 +
   4.232 +    /// \brief \ref named-templ-param "Named parameter" for setting
   4.233 +    /// \c %Path type.
   4.234 +    ///
   4.235 +    /// \ref named-templ-param "Named parameter" for setting the \c %Path
   4.236 +    /// type of the found cycles.
   4.237 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
   4.238 +    /// and it must have an \c addFront() function.
   4.239 +    template <typename T>
   4.240 +    struct SetPath
   4.241 +      : public HartmannOrlinMmc<GR, CM, SetPathTraits<T> > {
   4.242 +      typedef HartmannOrlinMmc<GR, CM, SetPathTraits<T> > Create;
   4.243 +    };
   4.244 +
   4.245 +    /// @}
   4.246 +
   4.247 +  protected:
   4.248 +
   4.249 +    HartmannOrlinMmc() {}
   4.250 +
   4.251 +  public:
   4.252 +
   4.253 +    /// \brief Constructor.
   4.254 +    ///
   4.255 +    /// The constructor of the class.
   4.256 +    ///
   4.257 +    /// \param digraph The digraph the algorithm runs on.
   4.258 +    /// \param cost The costs of the arcs.
   4.259 +    HartmannOrlinMmc( const Digraph &digraph,
   4.260 +                      const CostMap &cost ) :
   4.261 +      _gr(digraph), _cost(cost), _comp(digraph), _out_arcs(digraph),
   4.262 +      _best_found(false), _best_cost(0), _best_size(1),
   4.263 +      _cycle_path(NULL), _local_path(false), _data(digraph),
   4.264 +      INF(std::numeric_limits<LargeCost>::has_infinity ?
   4.265 +          std::numeric_limits<LargeCost>::infinity() :
   4.266 +          std::numeric_limits<LargeCost>::max())
   4.267 +    {}
   4.268 +
   4.269 +    /// Destructor.
   4.270 +    ~HartmannOrlinMmc() {
   4.271 +      if (_local_path) delete _cycle_path;
   4.272 +    }
   4.273 +
   4.274 +    /// \brief Set the path structure for storing the found cycle.
   4.275 +    ///
   4.276 +    /// This function sets an external path structure for storing the
   4.277 +    /// found cycle.
   4.278 +    ///
   4.279 +    /// If you don't call this function before calling \ref run() or
   4.280 +    /// \ref findCycleMean(), it will allocate a local \ref Path "path"
   4.281 +    /// structure. The destuctor deallocates this automatically
   4.282 +    /// allocated object, of course.
   4.283 +    ///
   4.284 +    /// \note The algorithm calls only the \ref lemon::Path::addFront()
   4.285 +    /// "addFront()" function of the given path structure.
   4.286 +    ///
   4.287 +    /// \return <tt>(*this)</tt>
   4.288 +    HartmannOrlinMmc& cycle(Path &path) {
   4.289 +      if (_local_path) {
   4.290 +        delete _cycle_path;
   4.291 +        _local_path = false;
   4.292 +      }
   4.293 +      _cycle_path = &path;
   4.294 +      return *this;
   4.295 +    }
   4.296 +
   4.297 +    /// \brief Set the tolerance used by the algorithm.
   4.298 +    ///
   4.299 +    /// This function sets the tolerance object used by the algorithm.
   4.300 +    ///
   4.301 +    /// \return <tt>(*this)</tt>
   4.302 +    HartmannOrlinMmc& tolerance(const Tolerance& tolerance) {
   4.303 +      _tolerance = tolerance;
   4.304 +      return *this;
   4.305 +    }
   4.306 +
   4.307 +    /// \brief Return a const reference to the tolerance.
   4.308 +    ///
   4.309 +    /// This function returns a const reference to the tolerance object
   4.310 +    /// used by the algorithm.
   4.311 +    const Tolerance& tolerance() const {
   4.312 +      return _tolerance;
   4.313 +    }
   4.314 +
   4.315 +    /// \name Execution control
   4.316 +    /// The simplest way to execute the algorithm is to call the \ref run()
   4.317 +    /// function.\n
   4.318 +    /// If you only need the minimum mean cost, you may call
   4.319 +    /// \ref findCycleMean().
   4.320 +
   4.321 +    /// @{
   4.322 +
   4.323 +    /// \brief Run the algorithm.
   4.324 +    ///
   4.325 +    /// This function runs the algorithm.
   4.326 +    /// It can be called more than once (e.g. if the underlying digraph
   4.327 +    /// and/or the arc costs have been modified).
   4.328 +    ///
   4.329 +    /// \return \c true if a directed cycle exists in the digraph.
   4.330 +    ///
   4.331 +    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.
   4.332 +    /// \code
   4.333 +    ///   return mmc.findCycleMean() && mmc.findCycle();
   4.334 +    /// \endcode
   4.335 +    bool run() {
   4.336 +      return findCycleMean() && findCycle();
   4.337 +    }
   4.338 +
   4.339 +    /// \brief Find the minimum cycle mean.
   4.340 +    ///
   4.341 +    /// This function finds the minimum mean cost of the directed
   4.342 +    /// cycles in the digraph.
   4.343 +    ///
   4.344 +    /// \return \c true if a directed cycle exists in the digraph.
   4.345 +    bool findCycleMean() {
   4.346 +      // Initialization and find strongly connected components
   4.347 +      init();
   4.348 +      findComponents();
   4.349 +      
   4.350 +      // Find the minimum cycle mean in the components
   4.351 +      for (int comp = 0; comp < _comp_num; ++comp) {
   4.352 +        if (!initComponent(comp)) continue;
   4.353 +        processRounds();
   4.354 +        
   4.355 +        // Update the best cycle (global minimum mean cycle)
   4.356 +        if ( _curr_found && (!_best_found || 
   4.357 +             _curr_cost * _best_size < _best_cost * _curr_size) ) {
   4.358 +          _best_found = true;
   4.359 +          _best_cost = _curr_cost;
   4.360 +          _best_size = _curr_size;
   4.361 +          _best_node = _curr_node;
   4.362 +          _best_level = _curr_level;
   4.363 +        }
   4.364 +      }
   4.365 +      return _best_found;
   4.366 +    }
   4.367 +
   4.368 +    /// \brief Find a minimum mean directed cycle.
   4.369 +    ///
   4.370 +    /// This function finds a directed cycle of minimum mean cost
   4.371 +    /// in the digraph using the data computed by findCycleMean().
   4.372 +    ///
   4.373 +    /// \return \c true if a directed cycle exists in the digraph.
   4.374 +    ///
   4.375 +    /// \pre \ref findCycleMean() must be called before using this function.
   4.376 +    bool findCycle() {
   4.377 +      if (!_best_found) return false;
   4.378 +      IntNodeMap reached(_gr, -1);
   4.379 +      int r = _best_level + 1;
   4.380 +      Node u = _best_node;
   4.381 +      while (reached[u] < 0) {
   4.382 +        reached[u] = --r;
   4.383 +        u = _gr.source(_data[u][r].pred);
   4.384 +      }
   4.385 +      r = reached[u];
   4.386 +      Arc e = _data[u][r].pred;
   4.387 +      _cycle_path->addFront(e);
   4.388 +      _best_cost = _cost[e];
   4.389 +      _best_size = 1;
   4.390 +      Node v;
   4.391 +      while ((v = _gr.source(e)) != u) {
   4.392 +        e = _data[v][--r].pred;
   4.393 +        _cycle_path->addFront(e);
   4.394 +        _best_cost += _cost[e];
   4.395 +        ++_best_size;
   4.396 +      }
   4.397 +      return true;
   4.398 +    }
   4.399 +
   4.400 +    /// @}
   4.401 +
   4.402 +    /// \name Query Functions
   4.403 +    /// The results of the algorithm can be obtained using these
   4.404 +    /// functions.\n
   4.405 +    /// The algorithm should be executed before using them.
   4.406 +
   4.407 +    /// @{
   4.408 +
   4.409 +    /// \brief Return the total cost of the found cycle.
   4.410 +    ///
   4.411 +    /// This function returns the total cost of the found cycle.
   4.412 +    ///
   4.413 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   4.414 +    /// using this function.
   4.415 +    Cost cycleCost() const {
   4.416 +      return static_cast<Cost>(_best_cost);
   4.417 +    }
   4.418 +
   4.419 +    /// \brief Return the number of arcs on the found cycle.
   4.420 +    ///
   4.421 +    /// This function returns the number of arcs on the found cycle.
   4.422 +    ///
   4.423 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   4.424 +    /// using this function.
   4.425 +    int cycleSize() const {
   4.426 +      return _best_size;
   4.427 +    }
   4.428 +
   4.429 +    /// \brief Return the mean cost of the found cycle.
   4.430 +    ///
   4.431 +    /// This function returns the mean cost of the found cycle.
   4.432 +    ///
   4.433 +    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the
   4.434 +    /// following code.
   4.435 +    /// \code
   4.436 +    ///   return static_cast<double>(alg.cycleCost()) / alg.cycleSize();
   4.437 +    /// \endcode
   4.438 +    ///
   4.439 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   4.440 +    /// using this function.
   4.441 +    double cycleMean() const {
   4.442 +      return static_cast<double>(_best_cost) / _best_size;
   4.443 +    }
   4.444 +
   4.445 +    /// \brief Return the found cycle.
   4.446 +    ///
   4.447 +    /// This function returns a const reference to the path structure
   4.448 +    /// storing the found cycle.
   4.449 +    ///
   4.450 +    /// \pre \ref run() or \ref findCycle() must be called before using
   4.451 +    /// this function.
   4.452 +    const Path& cycle() const {
   4.453 +      return *_cycle_path;
   4.454 +    }
   4.455 +
   4.456 +    ///@}
   4.457 +
   4.458 +  private:
   4.459 +
   4.460 +    // Initialization
   4.461 +    void init() {
   4.462 +      if (!_cycle_path) {
   4.463 +        _local_path = true;
   4.464 +        _cycle_path = new Path;
   4.465 +      }
   4.466 +      _cycle_path->clear();
   4.467 +      _best_found = false;
   4.468 +      _best_cost = 0;
   4.469 +      _best_size = 1;
   4.470 +      _cycle_path->clear();
   4.471 +      for (NodeIt u(_gr); u != INVALID; ++u)
   4.472 +        _data[u].clear();
   4.473 +    }
   4.474 +
   4.475 +    // Find strongly connected components and initialize _comp_nodes
   4.476 +    // and _out_arcs
   4.477 +    void findComponents() {
   4.478 +      _comp_num = stronglyConnectedComponents(_gr, _comp);
   4.479 +      _comp_nodes.resize(_comp_num);
   4.480 +      if (_comp_num == 1) {
   4.481 +        _comp_nodes[0].clear();
   4.482 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   4.483 +          _comp_nodes[0].push_back(n);
   4.484 +          _out_arcs[n].clear();
   4.485 +          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   4.486 +            _out_arcs[n].push_back(a);
   4.487 +          }
   4.488 +        }
   4.489 +      } else {
   4.490 +        for (int i = 0; i < _comp_num; ++i)
   4.491 +          _comp_nodes[i].clear();
   4.492 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   4.493 +          int k = _comp[n];
   4.494 +          _comp_nodes[k].push_back(n);
   4.495 +          _out_arcs[n].clear();
   4.496 +          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   4.497 +            if (_comp[_gr.target(a)] == k) _out_arcs[n].push_back(a);
   4.498 +          }
   4.499 +        }
   4.500 +      }
   4.501 +    }
   4.502 +
   4.503 +    // Initialize path data for the current component
   4.504 +    bool initComponent(int comp) {
   4.505 +      _nodes = &(_comp_nodes[comp]);
   4.506 +      int n = _nodes->size();
   4.507 +      if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) {
   4.508 +        return false;
   4.509 +      }      
   4.510 +      for (int i = 0; i < n; ++i) {
   4.511 +        _data[(*_nodes)[i]].resize(n + 1, PathData(INF));
   4.512 +      }
   4.513 +      return true;
   4.514 +    }
   4.515 +
   4.516 +    // Process all rounds of computing path data for the current component.
   4.517 +    // _data[v][k] is the cost of a shortest directed walk from the root
   4.518 +    // node to node v containing exactly k arcs.
   4.519 +    void processRounds() {
   4.520 +      Node start = (*_nodes)[0];
   4.521 +      _data[start][0] = PathData(0);
   4.522 +      _process.clear();
   4.523 +      _process.push_back(start);
   4.524 +
   4.525 +      int k, n = _nodes->size();
   4.526 +      int next_check = 4;
   4.527 +      bool terminate = false;
   4.528 +      for (k = 1; k <= n && int(_process.size()) < n && !terminate; ++k) {
   4.529 +        processNextBuildRound(k);
   4.530 +        if (k == next_check || k == n) {
   4.531 +          terminate = checkTermination(k);
   4.532 +          next_check = next_check * 3 / 2;
   4.533 +        }
   4.534 +      }
   4.535 +      for ( ; k <= n && !terminate; ++k) {
   4.536 +        processNextFullRound(k);
   4.537 +        if (k == next_check || k == n) {
   4.538 +          terminate = checkTermination(k);
   4.539 +          next_check = next_check * 3 / 2;
   4.540 +        }
   4.541 +      }
   4.542 +    }
   4.543 +
   4.544 +    // Process one round and rebuild _process
   4.545 +    void processNextBuildRound(int k) {
   4.546 +      std::vector<Node> next;
   4.547 +      Node u, v;
   4.548 +      Arc e;
   4.549 +      LargeCost d;
   4.550 +      for (int i = 0; i < int(_process.size()); ++i) {
   4.551 +        u = _process[i];
   4.552 +        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   4.553 +          e = _out_arcs[u][j];
   4.554 +          v = _gr.target(e);
   4.555 +          d = _data[u][k-1].dist + _cost[e];
   4.556 +          if (_tolerance.less(d, _data[v][k].dist)) {
   4.557 +            if (_data[v][k].dist == INF) next.push_back(v);
   4.558 +            _data[v][k] = PathData(d, e);
   4.559 +          }
   4.560 +        }
   4.561 +      }
   4.562 +      _process.swap(next);
   4.563 +    }
   4.564 +
   4.565 +    // Process one round using _nodes instead of _process
   4.566 +    void processNextFullRound(int k) {
   4.567 +      Node u, v;
   4.568 +      Arc e;
   4.569 +      LargeCost d;
   4.570 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   4.571 +        u = (*_nodes)[i];
   4.572 +        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   4.573 +          e = _out_arcs[u][j];
   4.574 +          v = _gr.target(e);
   4.575 +          d = _data[u][k-1].dist + _cost[e];
   4.576 +          if (_tolerance.less(d, _data[v][k].dist)) {
   4.577 +            _data[v][k] = PathData(d, e);
   4.578 +          }
   4.579 +        }
   4.580 +      }
   4.581 +    }
   4.582 +    
   4.583 +    // Check early termination
   4.584 +    bool checkTermination(int k) {
   4.585 +      typedef std::pair<int, int> Pair;
   4.586 +      typename GR::template NodeMap<Pair> level(_gr, Pair(-1, 0));
   4.587 +      typename GR::template NodeMap<LargeCost> pi(_gr);
   4.588 +      int n = _nodes->size();
   4.589 +      LargeCost cost;
   4.590 +      int size;
   4.591 +      Node u;
   4.592 +      
   4.593 +      // Search for cycles that are already found
   4.594 +      _curr_found = false;
   4.595 +      for (int i = 0; i < n; ++i) {
   4.596 +        u = (*_nodes)[i];
   4.597 +        if (_data[u][k].dist == INF) continue;
   4.598 +        for (int j = k; j >= 0; --j) {
   4.599 +          if (level[u].first == i && level[u].second > 0) {
   4.600 +            // A cycle is found
   4.601 +            cost = _data[u][level[u].second].dist - _data[u][j].dist;
   4.602 +            size = level[u].second - j;
   4.603 +            if (!_curr_found || cost * _curr_size < _curr_cost * size) {
   4.604 +              _curr_cost = cost;
   4.605 +              _curr_size = size;
   4.606 +              _curr_node = u;
   4.607 +              _curr_level = level[u].second;
   4.608 +              _curr_found = true;
   4.609 +            }
   4.610 +          }
   4.611 +          level[u] = Pair(i, j);
   4.612 +          if (j != 0) {
   4.613 +	    u = _gr.source(_data[u][j].pred);
   4.614 +	  }
   4.615 +        }
   4.616 +      }
   4.617 +
   4.618 +      // If at least one cycle is found, check the optimality condition
   4.619 +      LargeCost d;
   4.620 +      if (_curr_found && k < n) {
   4.621 +        // Find node potentials
   4.622 +        for (int i = 0; i < n; ++i) {
   4.623 +          u = (*_nodes)[i];
   4.624 +          pi[u] = INF;
   4.625 +          for (int j = 0; j <= k; ++j) {
   4.626 +            if (_data[u][j].dist < INF) {
   4.627 +              d = _data[u][j].dist * _curr_size - j * _curr_cost;
   4.628 +              if (_tolerance.less(d, pi[u])) pi[u] = d;
   4.629 +            }
   4.630 +          }
   4.631 +        }
   4.632 +
   4.633 +        // Check the optimality condition for all arcs
   4.634 +        bool done = true;
   4.635 +        for (ArcIt a(_gr); a != INVALID; ++a) {
   4.636 +          if (_tolerance.less(_cost[a] * _curr_size - _curr_cost,
   4.637 +                              pi[_gr.target(a)] - pi[_gr.source(a)]) ) {
   4.638 +            done = false;
   4.639 +            break;
   4.640 +          }
   4.641 +        }
   4.642 +        return done;
   4.643 +      }
   4.644 +      return (k == n);
   4.645 +    }
   4.646 +
   4.647 +  }; //class HartmannOrlinMmc
   4.648 +
   4.649 +  ///@}
   4.650 +
   4.651 +} //namespace lemon
   4.652 +
   4.653 +#endif //LEMON_HARTMANN_ORLIN_MMC_H
     5.1 --- a/lemon/howard.h	Mon Mar 08 08:33:41 2010 +0100
     5.2 +++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
     5.3 @@ -1,605 +0,0 @@
     5.4 -/* -*- C++ -*-
     5.5 - *
     5.6 - * This file is a part of LEMON, a generic C++ optimization library
     5.7 - *
     5.8 - * Copyright (C) 2003-2008
     5.9 - * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    5.10 - * (Egervary Research Group on Combinatorial Optimization, EGRES).
    5.11 - *
    5.12 - * Permission to use, modify and distribute this software is granted
    5.13 - * provided that this copyright notice appears in all copies. For
    5.14 - * precise terms see the accompanying LICENSE file.
    5.15 - *
    5.16 - * This software is provided "AS IS" with no warranty of any kind,
    5.17 - * express or implied, and with no claim as to its suitability for any
    5.18 - * purpose.
    5.19 - *
    5.20 - */
    5.21 -
    5.22 -#ifndef LEMON_HOWARD_H
    5.23 -#define LEMON_HOWARD_H
    5.24 -
    5.25 -/// \ingroup min_mean_cycle
    5.26 -///
    5.27 -/// \file
    5.28 -/// \brief Howard's algorithm for finding a minimum mean cycle.
    5.29 -
    5.30 -#include <vector>
    5.31 -#include <limits>
    5.32 -#include <lemon/core.h>
    5.33 -#include <lemon/path.h>
    5.34 -#include <lemon/tolerance.h>
    5.35 -#include <lemon/connectivity.h>
    5.36 -
    5.37 -namespace lemon {
    5.38 -
    5.39 -  /// \brief Default traits class of Howard class.
    5.40 -  ///
    5.41 -  /// Default traits class of Howard class.
    5.42 -  /// \tparam GR The type of the digraph.
    5.43 -  /// \tparam LEN The type of the length map.
    5.44 -  /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
    5.45 -#ifdef DOXYGEN
    5.46 -  template <typename GR, typename LEN>
    5.47 -#else
    5.48 -  template <typename GR, typename LEN,
    5.49 -    bool integer = std::numeric_limits<typename LEN::Value>::is_integer>
    5.50 -#endif
    5.51 -  struct HowardDefaultTraits
    5.52 -  {
    5.53 -    /// The type of the digraph
    5.54 -    typedef GR Digraph;
    5.55 -    /// The type of the length map
    5.56 -    typedef LEN LengthMap;
    5.57 -    /// The type of the arc lengths
    5.58 -    typedef typename LengthMap::Value Value;
    5.59 -
    5.60 -    /// \brief The large value type used for internal computations
    5.61 -    ///
    5.62 -    /// The large value type used for internal computations.
    5.63 -    /// It is \c long \c long if the \c Value type is integer,
    5.64 -    /// otherwise it is \c double.
    5.65 -    /// \c Value must be convertible to \c LargeValue.
    5.66 -    typedef double LargeValue;
    5.67 -
    5.68 -    /// The tolerance type used for internal computations
    5.69 -    typedef lemon::Tolerance<LargeValue> Tolerance;
    5.70 -
    5.71 -    /// \brief The path type of the found cycles
    5.72 -    ///
    5.73 -    /// The path type of the found cycles.
    5.74 -    /// It must conform to the \ref lemon::concepts::Path "Path" concept
    5.75 -    /// and it must have an \c addBack() function.
    5.76 -    typedef lemon::Path<Digraph> Path;
    5.77 -  };
    5.78 -
    5.79 -  // Default traits class for integer value types
    5.80 -  template <typename GR, typename LEN>
    5.81 -  struct HowardDefaultTraits<GR, LEN, true>
    5.82 -  {
    5.83 -    typedef GR Digraph;
    5.84 -    typedef LEN LengthMap;
    5.85 -    typedef typename LengthMap::Value Value;
    5.86 -#ifdef LEMON_HAVE_LONG_LONG
    5.87 -    typedef long long LargeValue;
    5.88 -#else
    5.89 -    typedef long LargeValue;
    5.90 -#endif
    5.91 -    typedef lemon::Tolerance<LargeValue> Tolerance;
    5.92 -    typedef lemon::Path<Digraph> Path;
    5.93 -  };
    5.94 -
    5.95 -
    5.96 -  /// \addtogroup min_mean_cycle
    5.97 -  /// @{
    5.98 -
    5.99 -  /// \brief Implementation of Howard's algorithm for finding a minimum
   5.100 -  /// mean cycle.
   5.101 -  ///
   5.102 -  /// This class implements Howard's policy iteration algorithm for finding
   5.103 -  /// a directed cycle of minimum mean length (cost) in a digraph
   5.104 -  /// \ref amo93networkflows, \ref dasdan98minmeancycle.
   5.105 -  /// This class provides the most efficient algorithm for the
   5.106 -  /// minimum mean cycle problem, though the best known theoretical
   5.107 -  /// bound on its running time is exponential.
   5.108 -  ///
   5.109 -  /// \tparam GR The type of the digraph the algorithm runs on.
   5.110 -  /// \tparam LEN The type of the length map. The default
   5.111 -  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
   5.112 -  /// \tparam TR The traits class that defines various types used by the
   5.113 -  /// algorithm. By default, it is \ref HowardDefaultTraits
   5.114 -  /// "HowardDefaultTraits<GR, LEN>".
   5.115 -  /// In most cases, this parameter should not be set directly,
   5.116 -  /// consider to use the named template parameters instead.
   5.117 -#ifdef DOXYGEN
   5.118 -  template <typename GR, typename LEN, typename TR>
   5.119 -#else
   5.120 -  template < typename GR,
   5.121 -             typename LEN = typename GR::template ArcMap<int>,
   5.122 -             typename TR = HowardDefaultTraits<GR, LEN> >
   5.123 -#endif
   5.124 -  class Howard
   5.125 -  {
   5.126 -  public:
   5.127 -  
   5.128 -    /// The type of the digraph
   5.129 -    typedef typename TR::Digraph Digraph;
   5.130 -    /// The type of the length map
   5.131 -    typedef typename TR::LengthMap LengthMap;
   5.132 -    /// The type of the arc lengths
   5.133 -    typedef typename TR::Value Value;
   5.134 -
   5.135 -    /// \brief The large value type
   5.136 -    ///
   5.137 -    /// The large value type used for internal computations.
   5.138 -    /// By default, it is \c long \c long if the \c Value type is integer,
   5.139 -    /// otherwise it is \c double.
   5.140 -    typedef typename TR::LargeValue LargeValue;
   5.141 -
   5.142 -    /// The tolerance type
   5.143 -    typedef typename TR::Tolerance Tolerance;
   5.144 -
   5.145 -    /// \brief The path type of the found cycles
   5.146 -    ///
   5.147 -    /// The path type of the found cycles.
   5.148 -    /// Using the \ref HowardDefaultTraits "default traits class",
   5.149 -    /// it is \ref lemon::Path "Path<Digraph>".
   5.150 -    typedef typename TR::Path Path;
   5.151 -
   5.152 -    /// The \ref HowardDefaultTraits "traits class" of the algorithm
   5.153 -    typedef TR Traits;
   5.154 -
   5.155 -  private:
   5.156 -
   5.157 -    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   5.158 -  
   5.159 -    // The digraph the algorithm runs on
   5.160 -    const Digraph &_gr;
   5.161 -    // The length of the arcs
   5.162 -    const LengthMap &_length;
   5.163 -
   5.164 -    // Data for the found cycles
   5.165 -    bool _curr_found, _best_found;
   5.166 -    LargeValue _curr_length, _best_length;
   5.167 -    int _curr_size, _best_size;
   5.168 -    Node _curr_node, _best_node;
   5.169 -
   5.170 -    Path *_cycle_path;
   5.171 -    bool _local_path;
   5.172 -
   5.173 -    // Internal data used by the algorithm
   5.174 -    typename Digraph::template NodeMap<Arc> _policy;
   5.175 -    typename Digraph::template NodeMap<bool> _reached;
   5.176 -    typename Digraph::template NodeMap<int> _level;
   5.177 -    typename Digraph::template NodeMap<LargeValue> _dist;
   5.178 -
   5.179 -    // Data for storing the strongly connected components
   5.180 -    int _comp_num;
   5.181 -    typename Digraph::template NodeMap<int> _comp;
   5.182 -    std::vector<std::vector<Node> > _comp_nodes;
   5.183 -    std::vector<Node>* _nodes;
   5.184 -    typename Digraph::template NodeMap<std::vector<Arc> > _in_arcs;
   5.185 -    
   5.186 -    // Queue used for BFS search
   5.187 -    std::vector<Node> _queue;
   5.188 -    int _qfront, _qback;
   5.189 -
   5.190 -    Tolerance _tolerance;
   5.191 -  
   5.192 -    // Infinite constant
   5.193 -    const LargeValue INF;
   5.194 -
   5.195 -  public:
   5.196 -  
   5.197 -    /// \name Named Template Parameters
   5.198 -    /// @{
   5.199 -
   5.200 -    template <typename T>
   5.201 -    struct SetLargeValueTraits : public Traits {
   5.202 -      typedef T LargeValue;
   5.203 -      typedef lemon::Tolerance<T> Tolerance;
   5.204 -    };
   5.205 -
   5.206 -    /// \brief \ref named-templ-param "Named parameter" for setting
   5.207 -    /// \c LargeValue type.
   5.208 -    ///
   5.209 -    /// \ref named-templ-param "Named parameter" for setting \c LargeValue
   5.210 -    /// type. It is used for internal computations in the algorithm.
   5.211 -    template <typename T>
   5.212 -    struct SetLargeValue
   5.213 -      : public Howard<GR, LEN, SetLargeValueTraits<T> > {
   5.214 -      typedef Howard<GR, LEN, SetLargeValueTraits<T> > Create;
   5.215 -    };
   5.216 -
   5.217 -    template <typename T>
   5.218 -    struct SetPathTraits : public Traits {
   5.219 -      typedef T Path;
   5.220 -    };
   5.221 -
   5.222 -    /// \brief \ref named-templ-param "Named parameter" for setting
   5.223 -    /// \c %Path type.
   5.224 -    ///
   5.225 -    /// \ref named-templ-param "Named parameter" for setting the \c %Path
   5.226 -    /// type of the found cycles.
   5.227 -    /// It must conform to the \ref lemon::concepts::Path "Path" concept
   5.228 -    /// and it must have an \c addBack() function.
   5.229 -    template <typename T>
   5.230 -    struct SetPath
   5.231 -      : public Howard<GR, LEN, SetPathTraits<T> > {
   5.232 -      typedef Howard<GR, LEN, SetPathTraits<T> > Create;
   5.233 -    };
   5.234 -    
   5.235 -    /// @}
   5.236 -
   5.237 -  protected:
   5.238 -
   5.239 -    Howard() {}
   5.240 -
   5.241 -  public:
   5.242 -
   5.243 -    /// \brief Constructor.
   5.244 -    ///
   5.245 -    /// The constructor of the class.
   5.246 -    ///
   5.247 -    /// \param digraph The digraph the algorithm runs on.
   5.248 -    /// \param length The lengths (costs) of the arcs.
   5.249 -    Howard( const Digraph &digraph,
   5.250 -            const LengthMap &length ) :
   5.251 -      _gr(digraph), _length(length), _best_found(false),
   5.252 -      _best_length(0), _best_size(1), _cycle_path(NULL), _local_path(false),
   5.253 -      _policy(digraph), _reached(digraph), _level(digraph), _dist(digraph),
   5.254 -      _comp(digraph), _in_arcs(digraph),
   5.255 -      INF(std::numeric_limits<LargeValue>::has_infinity ?
   5.256 -          std::numeric_limits<LargeValue>::infinity() :
   5.257 -          std::numeric_limits<LargeValue>::max())
   5.258 -    {}
   5.259 -
   5.260 -    /// Destructor.
   5.261 -    ~Howard() {
   5.262 -      if (_local_path) delete _cycle_path;
   5.263 -    }
   5.264 -
   5.265 -    /// \brief Set the path structure for storing the found cycle.
   5.266 -    ///
   5.267 -    /// This function sets an external path structure for storing the
   5.268 -    /// found cycle.
   5.269 -    ///
   5.270 -    /// If you don't call this function before calling \ref run() or
   5.271 -    /// \ref findMinMean(), it will allocate a local \ref Path "path"
   5.272 -    /// structure. The destuctor deallocates this automatically
   5.273 -    /// allocated object, of course.
   5.274 -    ///
   5.275 -    /// \note The algorithm calls only the \ref lemon::Path::addBack()
   5.276 -    /// "addBack()" function of the given path structure.
   5.277 -    ///
   5.278 -    /// \return <tt>(*this)</tt>
   5.279 -    Howard& cycle(Path &path) {
   5.280 -      if (_local_path) {
   5.281 -        delete _cycle_path;
   5.282 -        _local_path = false;
   5.283 -      }
   5.284 -      _cycle_path = &path;
   5.285 -      return *this;
   5.286 -    }
   5.287 -
   5.288 -    /// \brief Set the tolerance used by the algorithm.
   5.289 -    ///
   5.290 -    /// This function sets the tolerance object used by the algorithm.
   5.291 -    ///
   5.292 -    /// \return <tt>(*this)</tt>
   5.293 -    Howard& tolerance(const Tolerance& tolerance) {
   5.294 -      _tolerance = tolerance;
   5.295 -      return *this;
   5.296 -    }
   5.297 -
   5.298 -    /// \brief Return a const reference to the tolerance.
   5.299 -    ///
   5.300 -    /// This function returns a const reference to the tolerance object
   5.301 -    /// used by the algorithm.
   5.302 -    const Tolerance& tolerance() const {
   5.303 -      return _tolerance;
   5.304 -    }
   5.305 -
   5.306 -    /// \name Execution control
   5.307 -    /// The simplest way to execute the algorithm is to call the \ref run()
   5.308 -    /// function.\n
   5.309 -    /// If you only need the minimum mean length, you may call
   5.310 -    /// \ref findMinMean().
   5.311 -
   5.312 -    /// @{
   5.313 -
   5.314 -    /// \brief Run the algorithm.
   5.315 -    ///
   5.316 -    /// This function runs the algorithm.
   5.317 -    /// It can be called more than once (e.g. if the underlying digraph
   5.318 -    /// and/or the arc lengths have been modified).
   5.319 -    ///
   5.320 -    /// \return \c true if a directed cycle exists in the digraph.
   5.321 -    ///
   5.322 -    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.
   5.323 -    /// \code
   5.324 -    ///   return mmc.findMinMean() && mmc.findCycle();
   5.325 -    /// \endcode
   5.326 -    bool run() {
   5.327 -      return findMinMean() && findCycle();
   5.328 -    }
   5.329 -
   5.330 -    /// \brief Find the minimum cycle mean.
   5.331 -    ///
   5.332 -    /// This function finds the minimum mean length of the directed
   5.333 -    /// cycles in the digraph.
   5.334 -    ///
   5.335 -    /// \return \c true if a directed cycle exists in the digraph.
   5.336 -    bool findMinMean() {
   5.337 -      // Initialize and find strongly connected components
   5.338 -      init();
   5.339 -      findComponents();
   5.340 -      
   5.341 -      // Find the minimum cycle mean in the components
   5.342 -      for (int comp = 0; comp < _comp_num; ++comp) {
   5.343 -        // Find the minimum mean cycle in the current component
   5.344 -        if (!buildPolicyGraph(comp)) continue;
   5.345 -        while (true) {
   5.346 -          findPolicyCycle();
   5.347 -          if (!computeNodeDistances()) break;
   5.348 -        }
   5.349 -        // Update the best cycle (global minimum mean cycle)
   5.350 -        if ( _curr_found && (!_best_found ||
   5.351 -             _curr_length * _best_size < _best_length * _curr_size) ) {
   5.352 -          _best_found = true;
   5.353 -          _best_length = _curr_length;
   5.354 -          _best_size = _curr_size;
   5.355 -          _best_node = _curr_node;
   5.356 -        }
   5.357 -      }
   5.358 -      return _best_found;
   5.359 -    }
   5.360 -
   5.361 -    /// \brief Find a minimum mean directed cycle.
   5.362 -    ///
   5.363 -    /// This function finds a directed cycle of minimum mean length
   5.364 -    /// in the digraph using the data computed by findMinMean().
   5.365 -    ///
   5.366 -    /// \return \c true if a directed cycle exists in the digraph.
   5.367 -    ///
   5.368 -    /// \pre \ref findMinMean() must be called before using this function.
   5.369 -    bool findCycle() {
   5.370 -      if (!_best_found) return false;
   5.371 -      _cycle_path->addBack(_policy[_best_node]);
   5.372 -      for ( Node v = _best_node;
   5.373 -            (v = _gr.target(_policy[v])) != _best_node; ) {
   5.374 -        _cycle_path->addBack(_policy[v]);
   5.375 -      }
   5.376 -      return true;
   5.377 -    }
   5.378 -
   5.379 -    /// @}
   5.380 -
   5.381 -    /// \name Query Functions
   5.382 -    /// The results of the algorithm can be obtained using these
   5.383 -    /// functions.\n
   5.384 -    /// The algorithm should be executed before using them.
   5.385 -
   5.386 -    /// @{
   5.387 -
   5.388 -    /// \brief Return the total length of the found cycle.
   5.389 -    ///
   5.390 -    /// This function returns the total length of the found cycle.
   5.391 -    ///
   5.392 -    /// \pre \ref run() or \ref findMinMean() must be called before
   5.393 -    /// using this function.
   5.394 -    Value cycleLength() const {
   5.395 -      return static_cast<Value>(_best_length);
   5.396 -    }
   5.397 -
   5.398 -    /// \brief Return the number of arcs on the found cycle.
   5.399 -    ///
   5.400 -    /// This function returns the number of arcs on the found cycle.
   5.401 -    ///
   5.402 -    /// \pre \ref run() or \ref findMinMean() must be called before
   5.403 -    /// using this function.
   5.404 -    int cycleArcNum() const {
   5.405 -      return _best_size;
   5.406 -    }
   5.407 -
   5.408 -    /// \brief Return the mean length of the found cycle.
   5.409 -    ///
   5.410 -    /// This function returns the mean length of the found cycle.
   5.411 -    ///
   5.412 -    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the
   5.413 -    /// following code.
   5.414 -    /// \code
   5.415 -    ///   return static_cast<double>(alg.cycleLength()) / alg.cycleArcNum();
   5.416 -    /// \endcode
   5.417 -    ///
   5.418 -    /// \pre \ref run() or \ref findMinMean() must be called before
   5.419 -    /// using this function.
   5.420 -    double cycleMean() const {
   5.421 -      return static_cast<double>(_best_length) / _best_size;
   5.422 -    }
   5.423 -
   5.424 -    /// \brief Return the found cycle.
   5.425 -    ///
   5.426 -    /// This function returns a const reference to the path structure
   5.427 -    /// storing the found cycle.
   5.428 -    ///
   5.429 -    /// \pre \ref run() or \ref findCycle() must be called before using
   5.430 -    /// this function.
   5.431 -    const Path& cycle() const {
   5.432 -      return *_cycle_path;
   5.433 -    }
   5.434 -
   5.435 -    ///@}
   5.436 -
   5.437 -  private:
   5.438 -
   5.439 -    // Initialize
   5.440 -    void init() {
   5.441 -      if (!_cycle_path) {
   5.442 -        _local_path = true;
   5.443 -        _cycle_path = new Path;
   5.444 -      }
   5.445 -      _queue.resize(countNodes(_gr));
   5.446 -      _best_found = false;
   5.447 -      _best_length = 0;
   5.448 -      _best_size = 1;
   5.449 -      _cycle_path->clear();
   5.450 -    }
   5.451 -    
   5.452 -    // Find strongly connected components and initialize _comp_nodes
   5.453 -    // and _in_arcs
   5.454 -    void findComponents() {
   5.455 -      _comp_num = stronglyConnectedComponents(_gr, _comp);
   5.456 -      _comp_nodes.resize(_comp_num);
   5.457 -      if (_comp_num == 1) {
   5.458 -        _comp_nodes[0].clear();
   5.459 -        for (NodeIt n(_gr); n != INVALID; ++n) {
   5.460 -          _comp_nodes[0].push_back(n);
   5.461 -          _in_arcs[n].clear();
   5.462 -          for (InArcIt a(_gr, n); a != INVALID; ++a) {
   5.463 -            _in_arcs[n].push_back(a);
   5.464 -          }
   5.465 -        }
   5.466 -      } else {
   5.467 -        for (int i = 0; i < _comp_num; ++i)
   5.468 -          _comp_nodes[i].clear();
   5.469 -        for (NodeIt n(_gr); n != INVALID; ++n) {
   5.470 -          int k = _comp[n];
   5.471 -          _comp_nodes[k].push_back(n);
   5.472 -          _in_arcs[n].clear();
   5.473 -          for (InArcIt a(_gr, n); a != INVALID; ++a) {
   5.474 -            if (_comp[_gr.source(a)] == k) _in_arcs[n].push_back(a);
   5.475 -          }
   5.476 -        }
   5.477 -      }
   5.478 -    }
   5.479 -
   5.480 -    // Build the policy graph in the given strongly connected component
   5.481 -    // (the out-degree of every node is 1)
   5.482 -    bool buildPolicyGraph(int comp) {
   5.483 -      _nodes = &(_comp_nodes[comp]);
   5.484 -      if (_nodes->size() < 1 ||
   5.485 -          (_nodes->size() == 1 && _in_arcs[(*_nodes)[0]].size() == 0)) {
   5.486 -        return false;
   5.487 -      }
   5.488 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   5.489 -        _dist[(*_nodes)[i]] = INF;
   5.490 -      }
   5.491 -      Node u, v;
   5.492 -      Arc e;
   5.493 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   5.494 -        v = (*_nodes)[i];
   5.495 -        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   5.496 -          e = _in_arcs[v][j];
   5.497 -          u = _gr.source(e);
   5.498 -          if (_length[e] < _dist[u]) {
   5.499 -            _dist[u] = _length[e];
   5.500 -            _policy[u] = e;
   5.501 -          }
   5.502 -        }
   5.503 -      }
   5.504 -      return true;
   5.505 -    }
   5.506 -
   5.507 -    // Find the minimum mean cycle in the policy graph
   5.508 -    void findPolicyCycle() {
   5.509 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   5.510 -        _level[(*_nodes)[i]] = -1;
   5.511 -      }
   5.512 -      LargeValue clength;
   5.513 -      int csize;
   5.514 -      Node u, v;
   5.515 -      _curr_found = false;
   5.516 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   5.517 -        u = (*_nodes)[i];
   5.518 -        if (_level[u] >= 0) continue;
   5.519 -        for (; _level[u] < 0; u = _gr.target(_policy[u])) {
   5.520 -          _level[u] = i;
   5.521 -        }
   5.522 -        if (_level[u] == i) {
   5.523 -          // A cycle is found
   5.524 -          clength = _length[_policy[u]];
   5.525 -          csize = 1;
   5.526 -          for (v = u; (v = _gr.target(_policy[v])) != u; ) {
   5.527 -            clength += _length[_policy[v]];
   5.528 -            ++csize;
   5.529 -          }
   5.530 -          if ( !_curr_found ||
   5.531 -               (clength * _curr_size < _curr_length * csize) ) {
   5.532 -            _curr_found = true;
   5.533 -            _curr_length = clength;
   5.534 -            _curr_size = csize;
   5.535 -            _curr_node = u;
   5.536 -          }
   5.537 -        }
   5.538 -      }
   5.539 -    }
   5.540 -
   5.541 -    // Contract the policy graph and compute node distances
   5.542 -    bool computeNodeDistances() {
   5.543 -      // Find the component of the main cycle and compute node distances
   5.544 -      // using reverse BFS
   5.545 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   5.546 -        _reached[(*_nodes)[i]] = false;
   5.547 -      }
   5.548 -      _qfront = _qback = 0;
   5.549 -      _queue[0] = _curr_node;
   5.550 -      _reached[_curr_node] = true;
   5.551 -      _dist[_curr_node] = 0;
   5.552 -      Node u, v;
   5.553 -      Arc e;
   5.554 -      while (_qfront <= _qback) {
   5.555 -        v = _queue[_qfront++];
   5.556 -        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   5.557 -          e = _in_arcs[v][j];
   5.558 -          u = _gr.source(e);
   5.559 -          if (_policy[u] == e && !_reached[u]) {
   5.560 -            _reached[u] = true;
   5.561 -            _dist[u] = _dist[v] + _length[e] * _curr_size - _curr_length;
   5.562 -            _queue[++_qback] = u;
   5.563 -          }
   5.564 -        }
   5.565 -      }
   5.566 -
   5.567 -      // Connect all other nodes to this component and compute node
   5.568 -      // distances using reverse BFS
   5.569 -      _qfront = 0;
   5.570 -      while (_qback < int(_nodes->size())-1) {
   5.571 -        v = _queue[_qfront++];
   5.572 -        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   5.573 -          e = _in_arcs[v][j];
   5.574 -          u = _gr.source(e);
   5.575 -          if (!_reached[u]) {
   5.576 -            _reached[u] = true;
   5.577 -            _policy[u] = e;
   5.578 -            _dist[u] = _dist[v] + _length[e] * _curr_size - _curr_length;
   5.579 -            _queue[++_qback] = u;
   5.580 -          }
   5.581 -        }
   5.582 -      }
   5.583 -
   5.584 -      // Improve node distances
   5.585 -      bool improved = false;
   5.586 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   5.587 -        v = (*_nodes)[i];
   5.588 -        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   5.589 -          e = _in_arcs[v][j];
   5.590 -          u = _gr.source(e);
   5.591 -          LargeValue delta = _dist[v] + _length[e] * _curr_size - _curr_length;
   5.592 -          if (_tolerance.less(delta, _dist[u])) {
   5.593 -            _dist[u] = delta;
   5.594 -            _policy[u] = e;
   5.595 -            improved = true;
   5.596 -          }
   5.597 -        }
   5.598 -      }
   5.599 -      return improved;
   5.600 -    }
   5.601 -
   5.602 -  }; //class Howard
   5.603 -
   5.604 -  ///@}
   5.605 -
   5.606 -} //namespace lemon
   5.607 -
   5.608 -#endif //LEMON_HOWARD_H
     6.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     6.2 +++ b/lemon/howard_mmc.h	Sat Mar 13 22:01:38 2010 +0100
     6.3 @@ -0,0 +1,605 @@
     6.4 +/* -*- C++ -*-
     6.5 + *
     6.6 + * This file is a part of LEMON, a generic C++ optimization library
     6.7 + *
     6.8 + * Copyright (C) 2003-2008
     6.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    6.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    6.11 + *
    6.12 + * Permission to use, modify and distribute this software is granted
    6.13 + * provided that this copyright notice appears in all copies. For
    6.14 + * precise terms see the accompanying LICENSE file.
    6.15 + *
    6.16 + * This software is provided "AS IS" with no warranty of any kind,
    6.17 + * express or implied, and with no claim as to its suitability for any
    6.18 + * purpose.
    6.19 + *
    6.20 + */
    6.21 +
    6.22 +#ifndef LEMON_HOWARD_MMC_H
    6.23 +#define LEMON_HOWARD_MMC_H
    6.24 +
    6.25 +/// \ingroup min_mean_cycle
    6.26 +///
    6.27 +/// \file
    6.28 +/// \brief Howard's algorithm for finding a minimum mean cycle.
    6.29 +
    6.30 +#include <vector>
    6.31 +#include <limits>
    6.32 +#include <lemon/core.h>
    6.33 +#include <lemon/path.h>
    6.34 +#include <lemon/tolerance.h>
    6.35 +#include <lemon/connectivity.h>
    6.36 +
    6.37 +namespace lemon {
    6.38 +
    6.39 +  /// \brief Default traits class of HowardMmc class.
    6.40 +  ///
    6.41 +  /// Default traits class of HowardMmc class.
    6.42 +  /// \tparam GR The type of the digraph.
    6.43 +  /// \tparam CM The type of the cost map.
    6.44 +  /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
    6.45 +#ifdef DOXYGEN
    6.46 +  template <typename GR, typename CM>
    6.47 +#else
    6.48 +  template <typename GR, typename CM,
    6.49 +    bool integer = std::numeric_limits<typename CM::Value>::is_integer>
    6.50 +#endif
    6.51 +  struct HowardMmcDefaultTraits
    6.52 +  {
    6.53 +    /// The type of the digraph
    6.54 +    typedef GR Digraph;
    6.55 +    /// The type of the cost map
    6.56 +    typedef CM CostMap;
    6.57 +    /// The type of the arc costs
    6.58 +    typedef typename CostMap::Value Cost;
    6.59 +
    6.60 +    /// \brief The large cost type used for internal computations
    6.61 +    ///
    6.62 +    /// The large cost type used for internal computations.
    6.63 +    /// It is \c long \c long if the \c Cost type is integer,
    6.64 +    /// otherwise it is \c double.
    6.65 +    /// \c Cost must be convertible to \c LargeCost.
    6.66 +    typedef double LargeCost;
    6.67 +
    6.68 +    /// The tolerance type used for internal computations
    6.69 +    typedef lemon::Tolerance<LargeCost> Tolerance;
    6.70 +
    6.71 +    /// \brief The path type of the found cycles
    6.72 +    ///
    6.73 +    /// The path type of the found cycles.
    6.74 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
    6.75 +    /// and it must have an \c addBack() function.
    6.76 +    typedef lemon::Path<Digraph> Path;
    6.77 +  };
    6.78 +
    6.79 +  // Default traits class for integer cost types
    6.80 +  template <typename GR, typename CM>
    6.81 +  struct HowardMmcDefaultTraits<GR, CM, true>
    6.82 +  {
    6.83 +    typedef GR Digraph;
    6.84 +    typedef CM CostMap;
    6.85 +    typedef typename CostMap::Value Cost;
    6.86 +#ifdef LEMON_HAVE_LONG_LONG
    6.87 +    typedef long long LargeCost;
    6.88 +#else
    6.89 +    typedef long LargeCost;
    6.90 +#endif
    6.91 +    typedef lemon::Tolerance<LargeCost> Tolerance;
    6.92 +    typedef lemon::Path<Digraph> Path;
    6.93 +  };
    6.94 +
    6.95 +
    6.96 +  /// \addtogroup min_mean_cycle
    6.97 +  /// @{
    6.98 +
    6.99 +  /// \brief Implementation of Howard's algorithm for finding a minimum
   6.100 +  /// mean cycle.
   6.101 +  ///
   6.102 +  /// This class implements Howard's policy iteration algorithm for finding
   6.103 +  /// a directed cycle of minimum mean cost in a digraph
   6.104 +  /// \ref amo93networkflows, \ref dasdan98minmeancycle.
   6.105 +  /// This class provides the most efficient algorithm for the
   6.106 +  /// minimum mean cycle problem, though the best known theoretical
   6.107 +  /// bound on its running time is exponential.
   6.108 +  ///
   6.109 +  /// \tparam GR The type of the digraph the algorithm runs on.
   6.110 +  /// \tparam CM The type of the cost map. The default
   6.111 +  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
   6.112 +  /// \tparam TR The traits class that defines various types used by the
   6.113 +  /// algorithm. By default, it is \ref HowardMmcDefaultTraits
   6.114 +  /// "HowardMmcDefaultTraits<GR, CM>".
   6.115 +  /// In most cases, this parameter should not be set directly,
   6.116 +  /// consider to use the named template parameters instead.
   6.117 +#ifdef DOXYGEN
   6.118 +  template <typename GR, typename CM, typename TR>
   6.119 +#else
   6.120 +  template < typename GR,
   6.121 +             typename CM = typename GR::template ArcMap<int>,
   6.122 +             typename TR = HowardMmcDefaultTraits<GR, CM> >
   6.123 +#endif
   6.124 +  class HowardMmc
   6.125 +  {
   6.126 +  public:
   6.127 +  
   6.128 +    /// The type of the digraph
   6.129 +    typedef typename TR::Digraph Digraph;
   6.130 +    /// The type of the cost map
   6.131 +    typedef typename TR::CostMap CostMap;
   6.132 +    /// The type of the arc costs
   6.133 +    typedef typename TR::Cost Cost;
   6.134 +
   6.135 +    /// \brief The large cost type
   6.136 +    ///
   6.137 +    /// The large cost type used for internal computations.
   6.138 +    /// By default, it is \c long \c long if the \c Cost type is integer,
   6.139 +    /// otherwise it is \c double.
   6.140 +    typedef typename TR::LargeCost LargeCost;
   6.141 +
   6.142 +    /// The tolerance type
   6.143 +    typedef typename TR::Tolerance Tolerance;
   6.144 +
   6.145 +    /// \brief The path type of the found cycles
   6.146 +    ///
   6.147 +    /// The path type of the found cycles.
   6.148 +    /// Using the \ref HowardMmcDefaultTraits "default traits class",
   6.149 +    /// it is \ref lemon::Path "Path<Digraph>".
   6.150 +    typedef typename TR::Path Path;
   6.151 +
   6.152 +    /// The \ref HowardMmcDefaultTraits "traits class" of the algorithm
   6.153 +    typedef TR Traits;
   6.154 +
   6.155 +  private:
   6.156 +
   6.157 +    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   6.158 +  
   6.159 +    // The digraph the algorithm runs on
   6.160 +    const Digraph &_gr;
   6.161 +    // The cost of the arcs
   6.162 +    const CostMap &_cost;
   6.163 +
   6.164 +    // Data for the found cycles
   6.165 +    bool _curr_found, _best_found;
   6.166 +    LargeCost _curr_cost, _best_cost;
   6.167 +    int _curr_size, _best_size;
   6.168 +    Node _curr_node, _best_node;
   6.169 +
   6.170 +    Path *_cycle_path;
   6.171 +    bool _local_path;
   6.172 +
   6.173 +    // Internal data used by the algorithm
   6.174 +    typename Digraph::template NodeMap<Arc> _policy;
   6.175 +    typename Digraph::template NodeMap<bool> _reached;
   6.176 +    typename Digraph::template NodeMap<int> _level;
   6.177 +    typename Digraph::template NodeMap<LargeCost> _dist;
   6.178 +
   6.179 +    // Data for storing the strongly connected components
   6.180 +    int _comp_num;
   6.181 +    typename Digraph::template NodeMap<int> _comp;
   6.182 +    std::vector<std::vector<Node> > _comp_nodes;
   6.183 +    std::vector<Node>* _nodes;
   6.184 +    typename Digraph::template NodeMap<std::vector<Arc> > _in_arcs;
   6.185 +    
   6.186 +    // Queue used for BFS search
   6.187 +    std::vector<Node> _queue;
   6.188 +    int _qfront, _qback;
   6.189 +
   6.190 +    Tolerance _tolerance;
   6.191 +  
   6.192 +    // Infinite constant
   6.193 +    const LargeCost INF;
   6.194 +
   6.195 +  public:
   6.196 +  
   6.197 +    /// \name Named Template Parameters
   6.198 +    /// @{
   6.199 +
   6.200 +    template <typename T>
   6.201 +    struct SetLargeCostTraits : public Traits {
   6.202 +      typedef T LargeCost;
   6.203 +      typedef lemon::Tolerance<T> Tolerance;
   6.204 +    };
   6.205 +
   6.206 +    /// \brief \ref named-templ-param "Named parameter" for setting
   6.207 +    /// \c LargeCost type.
   6.208 +    ///
   6.209 +    /// \ref named-templ-param "Named parameter" for setting \c LargeCost
   6.210 +    /// type. It is used for internal computations in the algorithm.
   6.211 +    template <typename T>
   6.212 +    struct SetLargeCost
   6.213 +      : public HowardMmc<GR, CM, SetLargeCostTraits<T> > {
   6.214 +      typedef HowardMmc<GR, CM, SetLargeCostTraits<T> > Create;
   6.215 +    };
   6.216 +
   6.217 +    template <typename T>
   6.218 +    struct SetPathTraits : public Traits {
   6.219 +      typedef T Path;
   6.220 +    };
   6.221 +
   6.222 +    /// \brief \ref named-templ-param "Named parameter" for setting
   6.223 +    /// \c %Path type.
   6.224 +    ///
   6.225 +    /// \ref named-templ-param "Named parameter" for setting the \c %Path
   6.226 +    /// type of the found cycles.
   6.227 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
   6.228 +    /// and it must have an \c addBack() function.
   6.229 +    template <typename T>
   6.230 +    struct SetPath
   6.231 +      : public HowardMmc<GR, CM, SetPathTraits<T> > {
   6.232 +      typedef HowardMmc<GR, CM, SetPathTraits<T> > Create;
   6.233 +    };
   6.234 +    
   6.235 +    /// @}
   6.236 +
   6.237 +  protected:
   6.238 +
   6.239 +    HowardMmc() {}
   6.240 +
   6.241 +  public:
   6.242 +
   6.243 +    /// \brief Constructor.
   6.244 +    ///
   6.245 +    /// The constructor of the class.
   6.246 +    ///
   6.247 +    /// \param digraph The digraph the algorithm runs on.
   6.248 +    /// \param cost The costs of the arcs.
   6.249 +    HowardMmc( const Digraph &digraph,
   6.250 +               const CostMap &cost ) :
   6.251 +      _gr(digraph), _cost(cost), _best_found(false),
   6.252 +      _best_cost(0), _best_size(1), _cycle_path(NULL), _local_path(false),
   6.253 +      _policy(digraph), _reached(digraph), _level(digraph), _dist(digraph),
   6.254 +      _comp(digraph), _in_arcs(digraph),
   6.255 +      INF(std::numeric_limits<LargeCost>::has_infinity ?
   6.256 +          std::numeric_limits<LargeCost>::infinity() :
   6.257 +          std::numeric_limits<LargeCost>::max())
   6.258 +    {}
   6.259 +
   6.260 +    /// Destructor.
   6.261 +    ~HowardMmc() {
   6.262 +      if (_local_path) delete _cycle_path;
   6.263 +    }
   6.264 +
   6.265 +    /// \brief Set the path structure for storing the found cycle.
   6.266 +    ///
   6.267 +    /// This function sets an external path structure for storing the
   6.268 +    /// found cycle.
   6.269 +    ///
   6.270 +    /// If you don't call this function before calling \ref run() or
   6.271 +    /// \ref findCycleMean(), it will allocate a local \ref Path "path"
   6.272 +    /// structure. The destuctor deallocates this automatically
   6.273 +    /// allocated object, of course.
   6.274 +    ///
   6.275 +    /// \note The algorithm calls only the \ref lemon::Path::addBack()
   6.276 +    /// "addBack()" function of the given path structure.
   6.277 +    ///
   6.278 +    /// \return <tt>(*this)</tt>
   6.279 +    HowardMmc& cycle(Path &path) {
   6.280 +      if (_local_path) {
   6.281 +        delete _cycle_path;
   6.282 +        _local_path = false;
   6.283 +      }
   6.284 +      _cycle_path = &path;
   6.285 +      return *this;
   6.286 +    }
   6.287 +
   6.288 +    /// \brief Set the tolerance used by the algorithm.
   6.289 +    ///
   6.290 +    /// This function sets the tolerance object used by the algorithm.
   6.291 +    ///
   6.292 +    /// \return <tt>(*this)</tt>
   6.293 +    HowardMmc& tolerance(const Tolerance& tolerance) {
   6.294 +      _tolerance = tolerance;
   6.295 +      return *this;
   6.296 +    }
   6.297 +
   6.298 +    /// \brief Return a const reference to the tolerance.
   6.299 +    ///
   6.300 +    /// This function returns a const reference to the tolerance object
   6.301 +    /// used by the algorithm.
   6.302 +    const Tolerance& tolerance() const {
   6.303 +      return _tolerance;
   6.304 +    }
   6.305 +
   6.306 +    /// \name Execution control
   6.307 +    /// The simplest way to execute the algorithm is to call the \ref run()
   6.308 +    /// function.\n
   6.309 +    /// If you only need the minimum mean cost, you may call
   6.310 +    /// \ref findCycleMean().
   6.311 +
   6.312 +    /// @{
   6.313 +
   6.314 +    /// \brief Run the algorithm.
   6.315 +    ///
   6.316 +    /// This function runs the algorithm.
   6.317 +    /// It can be called more than once (e.g. if the underlying digraph
   6.318 +    /// and/or the arc costs have been modified).
   6.319 +    ///
   6.320 +    /// \return \c true if a directed cycle exists in the digraph.
   6.321 +    ///
   6.322 +    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.
   6.323 +    /// \code
   6.324 +    ///   return mmc.findCycleMean() && mmc.findCycle();
   6.325 +    /// \endcode
   6.326 +    bool run() {
   6.327 +      return findCycleMean() && findCycle();
   6.328 +    }
   6.329 +
   6.330 +    /// \brief Find the minimum cycle mean.
   6.331 +    ///
   6.332 +    /// This function finds the minimum mean cost of the directed
   6.333 +    /// cycles in the digraph.
   6.334 +    ///
   6.335 +    /// \return \c true if a directed cycle exists in the digraph.
   6.336 +    bool findCycleMean() {
   6.337 +      // Initialize and find strongly connected components
   6.338 +      init();
   6.339 +      findComponents();
   6.340 +      
   6.341 +      // Find the minimum cycle mean in the components
   6.342 +      for (int comp = 0; comp < _comp_num; ++comp) {
   6.343 +        // Find the minimum mean cycle in the current component
   6.344 +        if (!buildPolicyGraph(comp)) continue;
   6.345 +        while (true) {
   6.346 +          findPolicyCycle();
   6.347 +          if (!computeNodeDistances()) break;
   6.348 +        }
   6.349 +        // Update the best cycle (global minimum mean cycle)
   6.350 +        if ( _curr_found && (!_best_found ||
   6.351 +             _curr_cost * _best_size < _best_cost * _curr_size) ) {
   6.352 +          _best_found = true;
   6.353 +          _best_cost = _curr_cost;
   6.354 +          _best_size = _curr_size;
   6.355 +          _best_node = _curr_node;
   6.356 +        }
   6.357 +      }
   6.358 +      return _best_found;
   6.359 +    }
   6.360 +
   6.361 +    /// \brief Find a minimum mean directed cycle.
   6.362 +    ///
   6.363 +    /// This function finds a directed cycle of minimum mean cost
   6.364 +    /// in the digraph using the data computed by findCycleMean().
   6.365 +    ///
   6.366 +    /// \return \c true if a directed cycle exists in the digraph.
   6.367 +    ///
   6.368 +    /// \pre \ref findCycleMean() must be called before using this function.
   6.369 +    bool findCycle() {
   6.370 +      if (!_best_found) return false;
   6.371 +      _cycle_path->addBack(_policy[_best_node]);
   6.372 +      for ( Node v = _best_node;
   6.373 +            (v = _gr.target(_policy[v])) != _best_node; ) {
   6.374 +        _cycle_path->addBack(_policy[v]);
   6.375 +      }
   6.376 +      return true;
   6.377 +    }
   6.378 +
   6.379 +    /// @}
   6.380 +
   6.381 +    /// \name Query Functions
   6.382 +    /// The results of the algorithm can be obtained using these
   6.383 +    /// functions.\n
   6.384 +    /// The algorithm should be executed before using them.
   6.385 +
   6.386 +    /// @{
   6.387 +
   6.388 +    /// \brief Return the total cost of the found cycle.
   6.389 +    ///
   6.390 +    /// This function returns the total cost of the found cycle.
   6.391 +    ///
   6.392 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   6.393 +    /// using this function.
   6.394 +    Cost cycleCost() const {
   6.395 +      return static_cast<Cost>(_best_cost);
   6.396 +    }
   6.397 +
   6.398 +    /// \brief Return the number of arcs on the found cycle.
   6.399 +    ///
   6.400 +    /// This function returns the number of arcs on the found cycle.
   6.401 +    ///
   6.402 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   6.403 +    /// using this function.
   6.404 +    int cycleSize() const {
   6.405 +      return _best_size;
   6.406 +    }
   6.407 +
   6.408 +    /// \brief Return the mean cost of the found cycle.
   6.409 +    ///
   6.410 +    /// This function returns the mean cost of the found cycle.
   6.411 +    ///
   6.412 +    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the
   6.413 +    /// following code.
   6.414 +    /// \code
   6.415 +    ///   return static_cast<double>(alg.cycleCost()) / alg.cycleSize();
   6.416 +    /// \endcode
   6.417 +    ///
   6.418 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   6.419 +    /// using this function.
   6.420 +    double cycleMean() const {
   6.421 +      return static_cast<double>(_best_cost) / _best_size;
   6.422 +    }
   6.423 +
   6.424 +    /// \brief Return the found cycle.
   6.425 +    ///
   6.426 +    /// This function returns a const reference to the path structure
   6.427 +    /// storing the found cycle.
   6.428 +    ///
   6.429 +    /// \pre \ref run() or \ref findCycle() must be called before using
   6.430 +    /// this function.
   6.431 +    const Path& cycle() const {
   6.432 +      return *_cycle_path;
   6.433 +    }
   6.434 +
   6.435 +    ///@}
   6.436 +
   6.437 +  private:
   6.438 +
   6.439 +    // Initialize
   6.440 +    void init() {
   6.441 +      if (!_cycle_path) {
   6.442 +        _local_path = true;
   6.443 +        _cycle_path = new Path;
   6.444 +      }
   6.445 +      _queue.resize(countNodes(_gr));
   6.446 +      _best_found = false;
   6.447 +      _best_cost = 0;
   6.448 +      _best_size = 1;
   6.449 +      _cycle_path->clear();
   6.450 +    }
   6.451 +    
   6.452 +    // Find strongly connected components and initialize _comp_nodes
   6.453 +    // and _in_arcs
   6.454 +    void findComponents() {
   6.455 +      _comp_num = stronglyConnectedComponents(_gr, _comp);
   6.456 +      _comp_nodes.resize(_comp_num);
   6.457 +      if (_comp_num == 1) {
   6.458 +        _comp_nodes[0].clear();
   6.459 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   6.460 +          _comp_nodes[0].push_back(n);
   6.461 +          _in_arcs[n].clear();
   6.462 +          for (InArcIt a(_gr, n); a != INVALID; ++a) {
   6.463 +            _in_arcs[n].push_back(a);
   6.464 +          }
   6.465 +        }
   6.466 +      } else {
   6.467 +        for (int i = 0; i < _comp_num; ++i)
   6.468 +          _comp_nodes[i].clear();
   6.469 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   6.470 +          int k = _comp[n];
   6.471 +          _comp_nodes[k].push_back(n);
   6.472 +          _in_arcs[n].clear();
   6.473 +          for (InArcIt a(_gr, n); a != INVALID; ++a) {
   6.474 +            if (_comp[_gr.source(a)] == k) _in_arcs[n].push_back(a);
   6.475 +          }
   6.476 +        }
   6.477 +      }
   6.478 +    }
   6.479 +
   6.480 +    // Build the policy graph in the given strongly connected component
   6.481 +    // (the out-degree of every node is 1)
   6.482 +    bool buildPolicyGraph(int comp) {
   6.483 +      _nodes = &(_comp_nodes[comp]);
   6.484 +      if (_nodes->size() < 1 ||
   6.485 +          (_nodes->size() == 1 && _in_arcs[(*_nodes)[0]].size() == 0)) {
   6.486 +        return false;
   6.487 +      }
   6.488 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   6.489 +        _dist[(*_nodes)[i]] = INF;
   6.490 +      }
   6.491 +      Node u, v;
   6.492 +      Arc e;
   6.493 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   6.494 +        v = (*_nodes)[i];
   6.495 +        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   6.496 +          e = _in_arcs[v][j];
   6.497 +          u = _gr.source(e);
   6.498 +          if (_cost[e] < _dist[u]) {
   6.499 +            _dist[u] = _cost[e];
   6.500 +            _policy[u] = e;
   6.501 +          }
   6.502 +        }
   6.503 +      }
   6.504 +      return true;
   6.505 +    }
   6.506 +
   6.507 +    // Find the minimum mean cycle in the policy graph
   6.508 +    void findPolicyCycle() {
   6.509 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   6.510 +        _level[(*_nodes)[i]] = -1;
   6.511 +      }
   6.512 +      LargeCost ccost;
   6.513 +      int csize;
   6.514 +      Node u, v;
   6.515 +      _curr_found = false;
   6.516 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   6.517 +        u = (*_nodes)[i];
   6.518 +        if (_level[u] >= 0) continue;
   6.519 +        for (; _level[u] < 0; u = _gr.target(_policy[u])) {
   6.520 +          _level[u] = i;
   6.521 +        }
   6.522 +        if (_level[u] == i) {
   6.523 +          // A cycle is found
   6.524 +          ccost = _cost[_policy[u]];
   6.525 +          csize = 1;
   6.526 +          for (v = u; (v = _gr.target(_policy[v])) != u; ) {
   6.527 +            ccost += _cost[_policy[v]];
   6.528 +            ++csize;
   6.529 +          }
   6.530 +          if ( !_curr_found ||
   6.531 +               (ccost * _curr_size < _curr_cost * csize) ) {
   6.532 +            _curr_found = true;
   6.533 +            _curr_cost = ccost;
   6.534 +            _curr_size = csize;
   6.535 +            _curr_node = u;
   6.536 +          }
   6.537 +        }
   6.538 +      }
   6.539 +    }
   6.540 +
   6.541 +    // Contract the policy graph and compute node distances
   6.542 +    bool computeNodeDistances() {
   6.543 +      // Find the component of the main cycle and compute node distances
   6.544 +      // using reverse BFS
   6.545 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   6.546 +        _reached[(*_nodes)[i]] = false;
   6.547 +      }
   6.548 +      _qfront = _qback = 0;
   6.549 +      _queue[0] = _curr_node;
   6.550 +      _reached[_curr_node] = true;
   6.551 +      _dist[_curr_node] = 0;
   6.552 +      Node u, v;
   6.553 +      Arc e;
   6.554 +      while (_qfront <= _qback) {
   6.555 +        v = _queue[_qfront++];
   6.556 +        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   6.557 +          e = _in_arcs[v][j];
   6.558 +          u = _gr.source(e);
   6.559 +          if (_policy[u] == e && !_reached[u]) {
   6.560 +            _reached[u] = true;
   6.561 +            _dist[u] = _dist[v] + _cost[e] * _curr_size - _curr_cost;
   6.562 +            _queue[++_qback] = u;
   6.563 +          }
   6.564 +        }
   6.565 +      }
   6.566 +
   6.567 +      // Connect all other nodes to this component and compute node
   6.568 +      // distances using reverse BFS
   6.569 +      _qfront = 0;
   6.570 +      while (_qback < int(_nodes->size())-1) {
   6.571 +        v = _queue[_qfront++];
   6.572 +        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   6.573 +          e = _in_arcs[v][j];
   6.574 +          u = _gr.source(e);
   6.575 +          if (!_reached[u]) {
   6.576 +            _reached[u] = true;
   6.577 +            _policy[u] = e;
   6.578 +            _dist[u] = _dist[v] + _cost[e] * _curr_size - _curr_cost;
   6.579 +            _queue[++_qback] = u;
   6.580 +          }
   6.581 +        }
   6.582 +      }
   6.583 +
   6.584 +      // Improve node distances
   6.585 +      bool improved = false;
   6.586 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   6.587 +        v = (*_nodes)[i];
   6.588 +        for (int j = 0; j < int(_in_arcs[v].size()); ++j) {
   6.589 +          e = _in_arcs[v][j];
   6.590 +          u = _gr.source(e);
   6.591 +          LargeCost delta = _dist[v] + _cost[e] * _curr_size - _curr_cost;
   6.592 +          if (_tolerance.less(delta, _dist[u])) {
   6.593 +            _dist[u] = delta;
   6.594 +            _policy[u] = e;
   6.595 +            improved = true;
   6.596 +          }
   6.597 +        }
   6.598 +      }
   6.599 +      return improved;
   6.600 +    }
   6.601 +
   6.602 +  }; //class HowardMmc
   6.603 +
   6.604 +  ///@}
   6.605 +
   6.606 +} //namespace lemon
   6.607 +
   6.608 +#endif //LEMON_HOWARD_MMC_H
     7.1 --- a/lemon/karp.h	Mon Mar 08 08:33:41 2010 +0100
     7.2 +++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
     7.3 @@ -1,590 +0,0 @@
     7.4 -/* -*- C++ -*-
     7.5 - *
     7.6 - * This file is a part of LEMON, a generic C++ optimization library
     7.7 - *
     7.8 - * Copyright (C) 2003-2008
     7.9 - * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    7.10 - * (Egervary Research Group on Combinatorial Optimization, EGRES).
    7.11 - *
    7.12 - * Permission to use, modify and distribute this software is granted
    7.13 - * provided that this copyright notice appears in all copies. For
    7.14 - * precise terms see the accompanying LICENSE file.
    7.15 - *
    7.16 - * This software is provided "AS IS" with no warranty of any kind,
    7.17 - * express or implied, and with no claim as to its suitability for any
    7.18 - * purpose.
    7.19 - *
    7.20 - */
    7.21 -
    7.22 -#ifndef LEMON_KARP_H
    7.23 -#define LEMON_KARP_H
    7.24 -
    7.25 -/// \ingroup min_mean_cycle
    7.26 -///
    7.27 -/// \file
    7.28 -/// \brief Karp's algorithm for finding a minimum mean cycle.
    7.29 -
    7.30 -#include <vector>
    7.31 -#include <limits>
    7.32 -#include <lemon/core.h>
    7.33 -#include <lemon/path.h>
    7.34 -#include <lemon/tolerance.h>
    7.35 -#include <lemon/connectivity.h>
    7.36 -
    7.37 -namespace lemon {
    7.38 -
    7.39 -  /// \brief Default traits class of Karp algorithm.
    7.40 -  ///
    7.41 -  /// Default traits class of Karp algorithm.
    7.42 -  /// \tparam GR The type of the digraph.
    7.43 -  /// \tparam LEN The type of the length map.
    7.44 -  /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
    7.45 -#ifdef DOXYGEN
    7.46 -  template <typename GR, typename LEN>
    7.47 -#else
    7.48 -  template <typename GR, typename LEN,
    7.49 -    bool integer = std::numeric_limits<typename LEN::Value>::is_integer>
    7.50 -#endif
    7.51 -  struct KarpDefaultTraits
    7.52 -  {
    7.53 -    /// The type of the digraph
    7.54 -    typedef GR Digraph;
    7.55 -    /// The type of the length map
    7.56 -    typedef LEN LengthMap;
    7.57 -    /// The type of the arc lengths
    7.58 -    typedef typename LengthMap::Value Value;
    7.59 -
    7.60 -    /// \brief The large value type used for internal computations
    7.61 -    ///
    7.62 -    /// The large value type used for internal computations.
    7.63 -    /// It is \c long \c long if the \c Value type is integer,
    7.64 -    /// otherwise it is \c double.
    7.65 -    /// \c Value must be convertible to \c LargeValue.
    7.66 -    typedef double LargeValue;
    7.67 -
    7.68 -    /// The tolerance type used for internal computations
    7.69 -    typedef lemon::Tolerance<LargeValue> Tolerance;
    7.70 -
    7.71 -    /// \brief The path type of the found cycles
    7.72 -    ///
    7.73 -    /// The path type of the found cycles.
    7.74 -    /// It must conform to the \ref lemon::concepts::Path "Path" concept
    7.75 -    /// and it must have an \c addFront() function.
    7.76 -    typedef lemon::Path<Digraph> Path;
    7.77 -  };
    7.78 -
    7.79 -  // Default traits class for integer value types
    7.80 -  template <typename GR, typename LEN>
    7.81 -  struct KarpDefaultTraits<GR, LEN, true>
    7.82 -  {
    7.83 -    typedef GR Digraph;
    7.84 -    typedef LEN LengthMap;
    7.85 -    typedef typename LengthMap::Value Value;
    7.86 -#ifdef LEMON_HAVE_LONG_LONG
    7.87 -    typedef long long LargeValue;
    7.88 -#else
    7.89 -    typedef long LargeValue;
    7.90 -#endif
    7.91 -    typedef lemon::Tolerance<LargeValue> Tolerance;
    7.92 -    typedef lemon::Path<Digraph> Path;
    7.93 -  };
    7.94 -
    7.95 -
    7.96 -  /// \addtogroup min_mean_cycle
    7.97 -  /// @{
    7.98 -
    7.99 -  /// \brief Implementation of Karp's algorithm for finding a minimum
   7.100 -  /// mean cycle.
   7.101 -  ///
   7.102 -  /// This class implements Karp's algorithm for finding a directed
   7.103 -  /// cycle of minimum mean length (cost) in a digraph
   7.104 -  /// \ref amo93networkflows, \ref dasdan98minmeancycle.
   7.105 -  /// It runs in time O(ne) and uses space O(n<sup>2</sup>+e).
   7.106 -  ///
   7.107 -  /// \tparam GR The type of the digraph the algorithm runs on.
   7.108 -  /// \tparam LEN The type of the length map. The default
   7.109 -  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
   7.110 -  /// \tparam TR The traits class that defines various types used by the
   7.111 -  /// algorithm. By default, it is \ref KarpDefaultTraits
   7.112 -  /// "KarpDefaultTraits<GR, LEN>".
   7.113 -  /// In most cases, this parameter should not be set directly,
   7.114 -  /// consider to use the named template parameters instead.
   7.115 -#ifdef DOXYGEN
   7.116 -  template <typename GR, typename LEN, typename TR>
   7.117 -#else
   7.118 -  template < typename GR,
   7.119 -             typename LEN = typename GR::template ArcMap<int>,
   7.120 -             typename TR = KarpDefaultTraits<GR, LEN> >
   7.121 -#endif
   7.122 -  class Karp
   7.123 -  {
   7.124 -  public:
   7.125 -
   7.126 -    /// The type of the digraph
   7.127 -    typedef typename TR::Digraph Digraph;
   7.128 -    /// The type of the length map
   7.129 -    typedef typename TR::LengthMap LengthMap;
   7.130 -    /// The type of the arc lengths
   7.131 -    typedef typename TR::Value Value;
   7.132 -
   7.133 -    /// \brief The large value type
   7.134 -    ///
   7.135 -    /// The large value type used for internal computations.
   7.136 -    /// By default, it is \c long \c long if the \c Value type is integer,
   7.137 -    /// otherwise it is \c double.
   7.138 -    typedef typename TR::LargeValue LargeValue;
   7.139 -
   7.140 -    /// The tolerance type
   7.141 -    typedef typename TR::Tolerance Tolerance;
   7.142 -
   7.143 -    /// \brief The path type of the found cycles
   7.144 -    ///
   7.145 -    /// The path type of the found cycles.
   7.146 -    /// Using the \ref KarpDefaultTraits "default traits class",
   7.147 -    /// it is \ref lemon::Path "Path<Digraph>".
   7.148 -    typedef typename TR::Path Path;
   7.149 -
   7.150 -    /// The \ref KarpDefaultTraits "traits class" of the algorithm
   7.151 -    typedef TR Traits;
   7.152 -
   7.153 -  private:
   7.154 -
   7.155 -    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   7.156 -
   7.157 -    // Data sturcture for path data
   7.158 -    struct PathData
   7.159 -    {
   7.160 -      LargeValue dist;
   7.161 -      Arc pred;
   7.162 -      PathData(LargeValue d, Arc p = INVALID) :
   7.163 -        dist(d), pred(p) {}
   7.164 -    };
   7.165 -
   7.166 -    typedef typename Digraph::template NodeMap<std::vector<PathData> >
   7.167 -      PathDataNodeMap;
   7.168 -
   7.169 -  private:
   7.170 -
   7.171 -    // The digraph the algorithm runs on
   7.172 -    const Digraph &_gr;
   7.173 -    // The length of the arcs
   7.174 -    const LengthMap &_length;
   7.175 -
   7.176 -    // Data for storing the strongly connected components
   7.177 -    int _comp_num;
   7.178 -    typename Digraph::template NodeMap<int> _comp;
   7.179 -    std::vector<std::vector<Node> > _comp_nodes;
   7.180 -    std::vector<Node>* _nodes;
   7.181 -    typename Digraph::template NodeMap<std::vector<Arc> > _out_arcs;
   7.182 -
   7.183 -    // Data for the found cycle
   7.184 -    LargeValue _cycle_length;
   7.185 -    int _cycle_size;
   7.186 -    Node _cycle_node;
   7.187 -
   7.188 -    Path *_cycle_path;
   7.189 -    bool _local_path;
   7.190 -
   7.191 -    // Node map for storing path data
   7.192 -    PathDataNodeMap _data;
   7.193 -    // The processed nodes in the last round
   7.194 -    std::vector<Node> _process;
   7.195 -
   7.196 -    Tolerance _tolerance;
   7.197 -    
   7.198 -    // Infinite constant
   7.199 -    const LargeValue INF;
   7.200 -
   7.201 -  public:
   7.202 -
   7.203 -    /// \name Named Template Parameters
   7.204 -    /// @{
   7.205 -
   7.206 -    template <typename T>
   7.207 -    struct SetLargeValueTraits : public Traits {
   7.208 -      typedef T LargeValue;
   7.209 -      typedef lemon::Tolerance<T> Tolerance;
   7.210 -    };
   7.211 -
   7.212 -    /// \brief \ref named-templ-param "Named parameter" for setting
   7.213 -    /// \c LargeValue type.
   7.214 -    ///
   7.215 -    /// \ref named-templ-param "Named parameter" for setting \c LargeValue
   7.216 -    /// type. It is used for internal computations in the algorithm.
   7.217 -    template <typename T>
   7.218 -    struct SetLargeValue
   7.219 -      : public Karp<GR, LEN, SetLargeValueTraits<T> > {
   7.220 -      typedef Karp<GR, LEN, SetLargeValueTraits<T> > Create;
   7.221 -    };
   7.222 -
   7.223 -    template <typename T>
   7.224 -    struct SetPathTraits : public Traits {
   7.225 -      typedef T Path;
   7.226 -    };
   7.227 -
   7.228 -    /// \brief \ref named-templ-param "Named parameter" for setting
   7.229 -    /// \c %Path type.
   7.230 -    ///
   7.231 -    /// \ref named-templ-param "Named parameter" for setting the \c %Path
   7.232 -    /// type of the found cycles.
   7.233 -    /// It must conform to the \ref lemon::concepts::Path "Path" concept
   7.234 -    /// and it must have an \c addFront() function.
   7.235 -    template <typename T>
   7.236 -    struct SetPath
   7.237 -      : public Karp<GR, LEN, SetPathTraits<T> > {
   7.238 -      typedef Karp<GR, LEN, SetPathTraits<T> > Create;
   7.239 -    };
   7.240 -
   7.241 -    /// @}
   7.242 -
   7.243 -  protected:
   7.244 -
   7.245 -    Karp() {}
   7.246 -
   7.247 -  public:
   7.248 -
   7.249 -    /// \brief Constructor.
   7.250 -    ///
   7.251 -    /// The constructor of the class.
   7.252 -    ///
   7.253 -    /// \param digraph The digraph the algorithm runs on.
   7.254 -    /// \param length The lengths (costs) of the arcs.
   7.255 -    Karp( const Digraph &digraph,
   7.256 -          const LengthMap &length ) :
   7.257 -      _gr(digraph), _length(length), _comp(digraph), _out_arcs(digraph),
   7.258 -      _cycle_length(0), _cycle_size(1), _cycle_node(INVALID),
   7.259 -      _cycle_path(NULL), _local_path(false), _data(digraph),
   7.260 -      INF(std::numeric_limits<LargeValue>::has_infinity ?
   7.261 -          std::numeric_limits<LargeValue>::infinity() :
   7.262 -          std::numeric_limits<LargeValue>::max())
   7.263 -    {}
   7.264 -
   7.265 -    /// Destructor.
   7.266 -    ~Karp() {
   7.267 -      if (_local_path) delete _cycle_path;
   7.268 -    }
   7.269 -
   7.270 -    /// \brief Set the path structure for storing the found cycle.
   7.271 -    ///
   7.272 -    /// This function sets an external path structure for storing the
   7.273 -    /// found cycle.
   7.274 -    ///
   7.275 -    /// If you don't call this function before calling \ref run() or
   7.276 -    /// \ref findMinMean(), it will allocate a local \ref Path "path"
   7.277 -    /// structure. The destuctor deallocates this automatically
   7.278 -    /// allocated object, of course.
   7.279 -    ///
   7.280 -    /// \note The algorithm calls only the \ref lemon::Path::addFront()
   7.281 -    /// "addFront()" function of the given path structure.
   7.282 -    ///
   7.283 -    /// \return <tt>(*this)</tt>
   7.284 -    Karp& cycle(Path &path) {
   7.285 -      if (_local_path) {
   7.286 -        delete _cycle_path;
   7.287 -        _local_path = false;
   7.288 -      }
   7.289 -      _cycle_path = &path;
   7.290 -      return *this;
   7.291 -    }
   7.292 -
   7.293 -    /// \brief Set the tolerance used by the algorithm.
   7.294 -    ///
   7.295 -    /// This function sets the tolerance object used by the algorithm.
   7.296 -    ///
   7.297 -    /// \return <tt>(*this)</tt>
   7.298 -    Karp& tolerance(const Tolerance& tolerance) {
   7.299 -      _tolerance = tolerance;
   7.300 -      return *this;
   7.301 -    }
   7.302 -
   7.303 -    /// \brief Return a const reference to the tolerance.
   7.304 -    ///
   7.305 -    /// This function returns a const reference to the tolerance object
   7.306 -    /// used by the algorithm.
   7.307 -    const Tolerance& tolerance() const {
   7.308 -      return _tolerance;
   7.309 -    }
   7.310 -
   7.311 -    /// \name Execution control
   7.312 -    /// The simplest way to execute the algorithm is to call the \ref run()
   7.313 -    /// function.\n
   7.314 -    /// If you only need the minimum mean length, you may call
   7.315 -    /// \ref findMinMean().
   7.316 -
   7.317 -    /// @{
   7.318 -
   7.319 -    /// \brief Run the algorithm.
   7.320 -    ///
   7.321 -    /// This function runs the algorithm.
   7.322 -    /// It can be called more than once (e.g. if the underlying digraph
   7.323 -    /// and/or the arc lengths have been modified).
   7.324 -    ///
   7.325 -    /// \return \c true if a directed cycle exists in the digraph.
   7.326 -    ///
   7.327 -    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.
   7.328 -    /// \code
   7.329 -    ///   return mmc.findMinMean() && mmc.findCycle();
   7.330 -    /// \endcode
   7.331 -    bool run() {
   7.332 -      return findMinMean() && findCycle();
   7.333 -    }
   7.334 -
   7.335 -    /// \brief Find the minimum cycle mean.
   7.336 -    ///
   7.337 -    /// This function finds the minimum mean length of the directed
   7.338 -    /// cycles in the digraph.
   7.339 -    ///
   7.340 -    /// \return \c true if a directed cycle exists in the digraph.
   7.341 -    bool findMinMean() {
   7.342 -      // Initialization and find strongly connected components
   7.343 -      init();
   7.344 -      findComponents();
   7.345 -      
   7.346 -      // Find the minimum cycle mean in the components
   7.347 -      for (int comp = 0; comp < _comp_num; ++comp) {
   7.348 -        if (!initComponent(comp)) continue;
   7.349 -        processRounds();
   7.350 -        updateMinMean();
   7.351 -      }
   7.352 -      return (_cycle_node != INVALID);
   7.353 -    }
   7.354 -
   7.355 -    /// \brief Find a minimum mean directed cycle.
   7.356 -    ///
   7.357 -    /// This function finds a directed cycle of minimum mean length
   7.358 -    /// in the digraph using the data computed by findMinMean().
   7.359 -    ///
   7.360 -    /// \return \c true if a directed cycle exists in the digraph.
   7.361 -    ///
   7.362 -    /// \pre \ref findMinMean() must be called before using this function.
   7.363 -    bool findCycle() {
   7.364 -      if (_cycle_node == INVALID) return false;
   7.365 -      IntNodeMap reached(_gr, -1);
   7.366 -      int r = _data[_cycle_node].size();
   7.367 -      Node u = _cycle_node;
   7.368 -      while (reached[u] < 0) {
   7.369 -        reached[u] = --r;
   7.370 -        u = _gr.source(_data[u][r].pred);
   7.371 -      }
   7.372 -      r = reached[u];
   7.373 -      Arc e = _data[u][r].pred;
   7.374 -      _cycle_path->addFront(e);
   7.375 -      _cycle_length = _length[e];
   7.376 -      _cycle_size = 1;
   7.377 -      Node v;
   7.378 -      while ((v = _gr.source(e)) != u) {
   7.379 -        e = _data[v][--r].pred;
   7.380 -        _cycle_path->addFront(e);
   7.381 -        _cycle_length += _length[e];
   7.382 -        ++_cycle_size;
   7.383 -      }
   7.384 -      return true;
   7.385 -    }
   7.386 -
   7.387 -    /// @}
   7.388 -
   7.389 -    /// \name Query Functions
   7.390 -    /// The results of the algorithm can be obtained using these
   7.391 -    /// functions.\n
   7.392 -    /// The algorithm should be executed before using them.
   7.393 -
   7.394 -    /// @{
   7.395 -
   7.396 -    /// \brief Return the total length of the found cycle.
   7.397 -    ///
   7.398 -    /// This function returns the total length of the found cycle.
   7.399 -    ///
   7.400 -    /// \pre \ref run() or \ref findMinMean() must be called before
   7.401 -    /// using this function.
   7.402 -    Value cycleLength() const {
   7.403 -      return static_cast<Value>(_cycle_length);
   7.404 -    }
   7.405 -
   7.406 -    /// \brief Return the number of arcs on the found cycle.
   7.407 -    ///
   7.408 -    /// This function returns the number of arcs on the found cycle.
   7.409 -    ///
   7.410 -    /// \pre \ref run() or \ref findMinMean() must be called before
   7.411 -    /// using this function.
   7.412 -    int cycleArcNum() const {
   7.413 -      return _cycle_size;
   7.414 -    }
   7.415 -
   7.416 -    /// \brief Return the mean length of the found cycle.
   7.417 -    ///
   7.418 -    /// This function returns the mean length of the found cycle.
   7.419 -    ///
   7.420 -    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the
   7.421 -    /// following code.
   7.422 -    /// \code
   7.423 -    ///   return static_cast<double>(alg.cycleLength()) / alg.cycleArcNum();
   7.424 -    /// \endcode
   7.425 -    ///
   7.426 -    /// \pre \ref run() or \ref findMinMean() must be called before
   7.427 -    /// using this function.
   7.428 -    double cycleMean() const {
   7.429 -      return static_cast<double>(_cycle_length) / _cycle_size;
   7.430 -    }
   7.431 -
   7.432 -    /// \brief Return the found cycle.
   7.433 -    ///
   7.434 -    /// This function returns a const reference to the path structure
   7.435 -    /// storing the found cycle.
   7.436 -    ///
   7.437 -    /// \pre \ref run() or \ref findCycle() must be called before using
   7.438 -    /// this function.
   7.439 -    const Path& cycle() const {
   7.440 -      return *_cycle_path;
   7.441 -    }
   7.442 -
   7.443 -    ///@}
   7.444 -
   7.445 -  private:
   7.446 -
   7.447 -    // Initialization
   7.448 -    void init() {
   7.449 -      if (!_cycle_path) {
   7.450 -        _local_path = true;
   7.451 -        _cycle_path = new Path;
   7.452 -      }
   7.453 -      _cycle_path->clear();
   7.454 -      _cycle_length = 0;
   7.455 -      _cycle_size = 1;
   7.456 -      _cycle_node = INVALID;
   7.457 -      for (NodeIt u(_gr); u != INVALID; ++u)
   7.458 -        _data[u].clear();
   7.459 -    }
   7.460 -
   7.461 -    // Find strongly connected components and initialize _comp_nodes
   7.462 -    // and _out_arcs
   7.463 -    void findComponents() {
   7.464 -      _comp_num = stronglyConnectedComponents(_gr, _comp);
   7.465 -      _comp_nodes.resize(_comp_num);
   7.466 -      if (_comp_num == 1) {
   7.467 -        _comp_nodes[0].clear();
   7.468 -        for (NodeIt n(_gr); n != INVALID; ++n) {
   7.469 -          _comp_nodes[0].push_back(n);
   7.470 -          _out_arcs[n].clear();
   7.471 -          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   7.472 -            _out_arcs[n].push_back(a);
   7.473 -          }
   7.474 -        }
   7.475 -      } else {
   7.476 -        for (int i = 0; i < _comp_num; ++i)
   7.477 -          _comp_nodes[i].clear();
   7.478 -        for (NodeIt n(_gr); n != INVALID; ++n) {
   7.479 -          int k = _comp[n];
   7.480 -          _comp_nodes[k].push_back(n);
   7.481 -          _out_arcs[n].clear();
   7.482 -          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   7.483 -            if (_comp[_gr.target(a)] == k) _out_arcs[n].push_back(a);
   7.484 -          }
   7.485 -        }
   7.486 -      }
   7.487 -    }
   7.488 -
   7.489 -    // Initialize path data for the current component
   7.490 -    bool initComponent(int comp) {
   7.491 -      _nodes = &(_comp_nodes[comp]);
   7.492 -      int n = _nodes->size();
   7.493 -      if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) {
   7.494 -        return false;
   7.495 -      }      
   7.496 -      for (int i = 0; i < n; ++i) {
   7.497 -        _data[(*_nodes)[i]].resize(n + 1, PathData(INF));
   7.498 -      }
   7.499 -      return true;
   7.500 -    }
   7.501 -
   7.502 -    // Process all rounds of computing path data for the current component.
   7.503 -    // _data[v][k] is the length of a shortest directed walk from the root
   7.504 -    // node to node v containing exactly k arcs.
   7.505 -    void processRounds() {
   7.506 -      Node start = (*_nodes)[0];
   7.507 -      _data[start][0] = PathData(0);
   7.508 -      _process.clear();
   7.509 -      _process.push_back(start);
   7.510 -
   7.511 -      int k, n = _nodes->size();
   7.512 -      for (k = 1; k <= n && int(_process.size()) < n; ++k) {
   7.513 -        processNextBuildRound(k);
   7.514 -      }
   7.515 -      for ( ; k <= n; ++k) {
   7.516 -        processNextFullRound(k);
   7.517 -      }
   7.518 -    }
   7.519 -
   7.520 -    // Process one round and rebuild _process
   7.521 -    void processNextBuildRound(int k) {
   7.522 -      std::vector<Node> next;
   7.523 -      Node u, v;
   7.524 -      Arc e;
   7.525 -      LargeValue d;
   7.526 -      for (int i = 0; i < int(_process.size()); ++i) {
   7.527 -        u = _process[i];
   7.528 -        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   7.529 -          e = _out_arcs[u][j];
   7.530 -          v = _gr.target(e);
   7.531 -          d = _data[u][k-1].dist + _length[e];
   7.532 -          if (_tolerance.less(d, _data[v][k].dist)) {
   7.533 -            if (_data[v][k].dist == INF) next.push_back(v);
   7.534 -            _data[v][k] = PathData(d, e);
   7.535 -          }
   7.536 -        }
   7.537 -      }
   7.538 -      _process.swap(next);
   7.539 -    }
   7.540 -
   7.541 -    // Process one round using _nodes instead of _process
   7.542 -    void processNextFullRound(int k) {
   7.543 -      Node u, v;
   7.544 -      Arc e;
   7.545 -      LargeValue d;
   7.546 -      for (int i = 0; i < int(_nodes->size()); ++i) {
   7.547 -        u = (*_nodes)[i];
   7.548 -        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   7.549 -          e = _out_arcs[u][j];
   7.550 -          v = _gr.target(e);
   7.551 -          d = _data[u][k-1].dist + _length[e];
   7.552 -          if (_tolerance.less(d, _data[v][k].dist)) {
   7.553 -            _data[v][k] = PathData(d, e);
   7.554 -          }
   7.555 -        }
   7.556 -      }
   7.557 -    }
   7.558 -
   7.559 -    // Update the minimum cycle mean
   7.560 -    void updateMinMean() {
   7.561 -      int n = _nodes->size();
   7.562 -      for (int i = 0; i < n; ++i) {
   7.563 -        Node u = (*_nodes)[i];
   7.564 -        if (_data[u][n].dist == INF) continue;
   7.565 -        LargeValue length, max_length = 0;
   7.566 -        int size, max_size = 1;
   7.567 -        bool found_curr = false;
   7.568 -        for (int k = 0; k < n; ++k) {
   7.569 -          if (_data[u][k].dist == INF) continue;
   7.570 -          length = _data[u][n].dist - _data[u][k].dist;
   7.571 -          size = n - k;
   7.572 -          if (!found_curr || length * max_size > max_length * size) {
   7.573 -            found_curr = true;
   7.574 -            max_length = length;
   7.575 -            max_size = size;
   7.576 -          }
   7.577 -        }
   7.578 -        if ( found_curr && (_cycle_node == INVALID ||
   7.579 -             max_length * _cycle_size < _cycle_length * max_size) ) {
   7.580 -          _cycle_length = max_length;
   7.581 -          _cycle_size = max_size;
   7.582 -          _cycle_node = u;
   7.583 -        }
   7.584 -      }
   7.585 -    }
   7.586 -
   7.587 -  }; //class Karp
   7.588 -
   7.589 -  ///@}
   7.590 -
   7.591 -} //namespace lemon
   7.592 -
   7.593 -#endif //LEMON_KARP_H
     8.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     8.2 +++ b/lemon/karp_mmc.h	Sat Mar 13 22:01:38 2010 +0100
     8.3 @@ -0,0 +1,590 @@
     8.4 +/* -*- C++ -*-
     8.5 + *
     8.6 + * This file is a part of LEMON, a generic C++ optimization library
     8.7 + *
     8.8 + * Copyright (C) 2003-2008
     8.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    8.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    8.11 + *
    8.12 + * Permission to use, modify and distribute this software is granted
    8.13 + * provided that this copyright notice appears in all copies. For
    8.14 + * precise terms see the accompanying LICENSE file.
    8.15 + *
    8.16 + * This software is provided "AS IS" with no warranty of any kind,
    8.17 + * express or implied, and with no claim as to its suitability for any
    8.18 + * purpose.
    8.19 + *
    8.20 + */
    8.21 +
    8.22 +#ifndef LEMON_KARP_MMC_H
    8.23 +#define LEMON_KARP_MMC_H
    8.24 +
    8.25 +/// \ingroup min_mean_cycle
    8.26 +///
    8.27 +/// \file
    8.28 +/// \brief Karp's algorithm for finding a minimum mean cycle.
    8.29 +
    8.30 +#include <vector>
    8.31 +#include <limits>
    8.32 +#include <lemon/core.h>
    8.33 +#include <lemon/path.h>
    8.34 +#include <lemon/tolerance.h>
    8.35 +#include <lemon/connectivity.h>
    8.36 +
    8.37 +namespace lemon {
    8.38 +
    8.39 +  /// \brief Default traits class of KarpMmc class.
    8.40 +  ///
    8.41 +  /// Default traits class of KarpMmc class.
    8.42 +  /// \tparam GR The type of the digraph.
    8.43 +  /// \tparam CM The type of the cost map.
    8.44 +  /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
    8.45 +#ifdef DOXYGEN
    8.46 +  template <typename GR, typename CM>
    8.47 +#else
    8.48 +  template <typename GR, typename CM,
    8.49 +    bool integer = std::numeric_limits<typename CM::Value>::is_integer>
    8.50 +#endif
    8.51 +  struct KarpMmcDefaultTraits
    8.52 +  {
    8.53 +    /// The type of the digraph
    8.54 +    typedef GR Digraph;
    8.55 +    /// The type of the cost map
    8.56 +    typedef CM CostMap;
    8.57 +    /// The type of the arc costs
    8.58 +    typedef typename CostMap::Value Cost;
    8.59 +
    8.60 +    /// \brief The large cost type used for internal computations
    8.61 +    ///
    8.62 +    /// The large cost type used for internal computations.
    8.63 +    /// It is \c long \c long if the \c Cost type is integer,
    8.64 +    /// otherwise it is \c double.
    8.65 +    /// \c Cost must be convertible to \c LargeCost.
    8.66 +    typedef double LargeCost;
    8.67 +
    8.68 +    /// The tolerance type used for internal computations
    8.69 +    typedef lemon::Tolerance<LargeCost> Tolerance;
    8.70 +
    8.71 +    /// \brief The path type of the found cycles
    8.72 +    ///
    8.73 +    /// The path type of the found cycles.
    8.74 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
    8.75 +    /// and it must have an \c addFront() function.
    8.76 +    typedef lemon::Path<Digraph> Path;
    8.77 +  };
    8.78 +
    8.79 +  // Default traits class for integer cost types
    8.80 +  template <typename GR, typename CM>
    8.81 +  struct KarpMmcDefaultTraits<GR, CM, true>
    8.82 +  {
    8.83 +    typedef GR Digraph;
    8.84 +    typedef CM CostMap;
    8.85 +    typedef typename CostMap::Value Cost;
    8.86 +#ifdef LEMON_HAVE_LONG_LONG
    8.87 +    typedef long long LargeCost;
    8.88 +#else
    8.89 +    typedef long LargeCost;
    8.90 +#endif
    8.91 +    typedef lemon::Tolerance<LargeCost> Tolerance;
    8.92 +    typedef lemon::Path<Digraph> Path;
    8.93 +  };
    8.94 +
    8.95 +
    8.96 +  /// \addtogroup min_mean_cycle
    8.97 +  /// @{
    8.98 +
    8.99 +  /// \brief Implementation of Karp's algorithm for finding a minimum
   8.100 +  /// mean cycle.
   8.101 +  ///
   8.102 +  /// This class implements Karp's algorithm for finding a directed
   8.103 +  /// cycle of minimum mean cost in a digraph
   8.104 +  /// \ref amo93networkflows, \ref dasdan98minmeancycle.
   8.105 +  /// It runs in time O(ne) and uses space O(n<sup>2</sup>+e).
   8.106 +  ///
   8.107 +  /// \tparam GR The type of the digraph the algorithm runs on.
   8.108 +  /// \tparam CM The type of the cost map. The default
   8.109 +  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
   8.110 +  /// \tparam TR The traits class that defines various types used by the
   8.111 +  /// algorithm. By default, it is \ref KarpMmcDefaultTraits
   8.112 +  /// "KarpMmcDefaultTraits<GR, CM>".
   8.113 +  /// In most cases, this parameter should not be set directly,
   8.114 +  /// consider to use the named template parameters instead.
   8.115 +#ifdef DOXYGEN
   8.116 +  template <typename GR, typename CM, typename TR>
   8.117 +#else
   8.118 +  template < typename GR,
   8.119 +             typename CM = typename GR::template ArcMap<int>,
   8.120 +             typename TR = KarpMmcDefaultTraits<GR, CM> >
   8.121 +#endif
   8.122 +  class KarpMmc
   8.123 +  {
   8.124 +  public:
   8.125 +
   8.126 +    /// The type of the digraph
   8.127 +    typedef typename TR::Digraph Digraph;
   8.128 +    /// The type of the cost map
   8.129 +    typedef typename TR::CostMap CostMap;
   8.130 +    /// The type of the arc costs
   8.131 +    typedef typename TR::Cost Cost;
   8.132 +
   8.133 +    /// \brief The large cost type
   8.134 +    ///
   8.135 +    /// The large cost type used for internal computations.
   8.136 +    /// By default, it is \c long \c long if the \c Cost type is integer,
   8.137 +    /// otherwise it is \c double.
   8.138 +    typedef typename TR::LargeCost LargeCost;
   8.139 +
   8.140 +    /// The tolerance type
   8.141 +    typedef typename TR::Tolerance Tolerance;
   8.142 +
   8.143 +    /// \brief The path type of the found cycles
   8.144 +    ///
   8.145 +    /// The path type of the found cycles.
   8.146 +    /// Using the \ref KarpMmcDefaultTraits "default traits class",
   8.147 +    /// it is \ref lemon::Path "Path<Digraph>".
   8.148 +    typedef typename TR::Path Path;
   8.149 +
   8.150 +    /// The \ref KarpMmcDefaultTraits "traits class" of the algorithm
   8.151 +    typedef TR Traits;
   8.152 +
   8.153 +  private:
   8.154 +
   8.155 +    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   8.156 +
   8.157 +    // Data sturcture for path data
   8.158 +    struct PathData
   8.159 +    {
   8.160 +      LargeCost dist;
   8.161 +      Arc pred;
   8.162 +      PathData(LargeCost d, Arc p = INVALID) :
   8.163 +        dist(d), pred(p) {}
   8.164 +    };
   8.165 +
   8.166 +    typedef typename Digraph::template NodeMap<std::vector<PathData> >
   8.167 +      PathDataNodeMap;
   8.168 +
   8.169 +  private:
   8.170 +
   8.171 +    // The digraph the algorithm runs on
   8.172 +    const Digraph &_gr;
   8.173 +    // The cost of the arcs
   8.174 +    const CostMap &_cost;
   8.175 +
   8.176 +    // Data for storing the strongly connected components
   8.177 +    int _comp_num;
   8.178 +    typename Digraph::template NodeMap<int> _comp;
   8.179 +    std::vector<std::vector<Node> > _comp_nodes;
   8.180 +    std::vector<Node>* _nodes;
   8.181 +    typename Digraph::template NodeMap<std::vector<Arc> > _out_arcs;
   8.182 +
   8.183 +    // Data for the found cycle
   8.184 +    LargeCost _cycle_cost;
   8.185 +    int _cycle_size;
   8.186 +    Node _cycle_node;
   8.187 +
   8.188 +    Path *_cycle_path;
   8.189 +    bool _local_path;
   8.190 +
   8.191 +    // Node map for storing path data
   8.192 +    PathDataNodeMap _data;
   8.193 +    // The processed nodes in the last round
   8.194 +    std::vector<Node> _process;
   8.195 +
   8.196 +    Tolerance _tolerance;
   8.197 +    
   8.198 +    // Infinite constant
   8.199 +    const LargeCost INF;
   8.200 +
   8.201 +  public:
   8.202 +
   8.203 +    /// \name Named Template Parameters
   8.204 +    /// @{
   8.205 +
   8.206 +    template <typename T>
   8.207 +    struct SetLargeCostTraits : public Traits {
   8.208 +      typedef T LargeCost;
   8.209 +      typedef lemon::Tolerance<T> Tolerance;
   8.210 +    };
   8.211 +
   8.212 +    /// \brief \ref named-templ-param "Named parameter" for setting
   8.213 +    /// \c LargeCost type.
   8.214 +    ///
   8.215 +    /// \ref named-templ-param "Named parameter" for setting \c LargeCost
   8.216 +    /// type. It is used for internal computations in the algorithm.
   8.217 +    template <typename T>
   8.218 +    struct SetLargeCost
   8.219 +      : public KarpMmc<GR, CM, SetLargeCostTraits<T> > {
   8.220 +      typedef KarpMmc<GR, CM, SetLargeCostTraits<T> > Create;
   8.221 +    };
   8.222 +
   8.223 +    template <typename T>
   8.224 +    struct SetPathTraits : public Traits {
   8.225 +      typedef T Path;
   8.226 +    };
   8.227 +
   8.228 +    /// \brief \ref named-templ-param "Named parameter" for setting
   8.229 +    /// \c %Path type.
   8.230 +    ///
   8.231 +    /// \ref named-templ-param "Named parameter" for setting the \c %Path
   8.232 +    /// type of the found cycles.
   8.233 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
   8.234 +    /// and it must have an \c addFront() function.
   8.235 +    template <typename T>
   8.236 +    struct SetPath
   8.237 +      : public KarpMmc<GR, CM, SetPathTraits<T> > {
   8.238 +      typedef KarpMmc<GR, CM, SetPathTraits<T> > Create;
   8.239 +    };
   8.240 +
   8.241 +    /// @}
   8.242 +
   8.243 +  protected:
   8.244 +
   8.245 +    KarpMmc() {}
   8.246 +
   8.247 +  public:
   8.248 +
   8.249 +    /// \brief Constructor.
   8.250 +    ///
   8.251 +    /// The constructor of the class.
   8.252 +    ///
   8.253 +    /// \param digraph The digraph the algorithm runs on.
   8.254 +    /// \param cost The costs of the arcs.
   8.255 +    KarpMmc( const Digraph &digraph,
   8.256 +             const CostMap &cost ) :
   8.257 +      _gr(digraph), _cost(cost), _comp(digraph), _out_arcs(digraph),
   8.258 +      _cycle_cost(0), _cycle_size(1), _cycle_node(INVALID),
   8.259 +      _cycle_path(NULL), _local_path(false), _data(digraph),
   8.260 +      INF(std::numeric_limits<LargeCost>::has_infinity ?
   8.261 +          std::numeric_limits<LargeCost>::infinity() :
   8.262 +          std::numeric_limits<LargeCost>::max())
   8.263 +    {}
   8.264 +
   8.265 +    /// Destructor.
   8.266 +    ~KarpMmc() {
   8.267 +      if (_local_path) delete _cycle_path;
   8.268 +    }
   8.269 +
   8.270 +    /// \brief Set the path structure for storing the found cycle.
   8.271 +    ///
   8.272 +    /// This function sets an external path structure for storing the
   8.273 +    /// found cycle.
   8.274 +    ///
   8.275 +    /// If you don't call this function before calling \ref run() or
   8.276 +    /// \ref findCycleMean(), it will allocate a local \ref Path "path"
   8.277 +    /// structure. The destuctor deallocates this automatically
   8.278 +    /// allocated object, of course.
   8.279 +    ///
   8.280 +    /// \note The algorithm calls only the \ref lemon::Path::addFront()
   8.281 +    /// "addFront()" function of the given path structure.
   8.282 +    ///
   8.283 +    /// \return <tt>(*this)</tt>
   8.284 +    KarpMmc& cycle(Path &path) {
   8.285 +      if (_local_path) {
   8.286 +        delete _cycle_path;
   8.287 +        _local_path = false;
   8.288 +      }
   8.289 +      _cycle_path = &path;
   8.290 +      return *this;
   8.291 +    }
   8.292 +
   8.293 +    /// \brief Set the tolerance used by the algorithm.
   8.294 +    ///
   8.295 +    /// This function sets the tolerance object used by the algorithm.
   8.296 +    ///
   8.297 +    /// \return <tt>(*this)</tt>
   8.298 +    KarpMmc& tolerance(const Tolerance& tolerance) {
   8.299 +      _tolerance = tolerance;
   8.300 +      return *this;
   8.301 +    }
   8.302 +
   8.303 +    /// \brief Return a const reference to the tolerance.
   8.304 +    ///
   8.305 +    /// This function returns a const reference to the tolerance object
   8.306 +    /// used by the algorithm.
   8.307 +    const Tolerance& tolerance() const {
   8.308 +      return _tolerance;
   8.309 +    }
   8.310 +
   8.311 +    /// \name Execution control
   8.312 +    /// The simplest way to execute the algorithm is to call the \ref run()
   8.313 +    /// function.\n
   8.314 +    /// If you only need the minimum mean cost, you may call
   8.315 +    /// \ref findCycleMean().
   8.316 +
   8.317 +    /// @{
   8.318 +
   8.319 +    /// \brief Run the algorithm.
   8.320 +    ///
   8.321 +    /// This function runs the algorithm.
   8.322 +    /// It can be called more than once (e.g. if the underlying digraph
   8.323 +    /// and/or the arc costs have been modified).
   8.324 +    ///
   8.325 +    /// \return \c true if a directed cycle exists in the digraph.
   8.326 +    ///
   8.327 +    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.
   8.328 +    /// \code
   8.329 +    ///   return mmc.findCycleMean() && mmc.findCycle();
   8.330 +    /// \endcode
   8.331 +    bool run() {
   8.332 +      return findCycleMean() && findCycle();
   8.333 +    }
   8.334 +
   8.335 +    /// \brief Find the minimum cycle mean.
   8.336 +    ///
   8.337 +    /// This function finds the minimum mean cost of the directed
   8.338 +    /// cycles in the digraph.
   8.339 +    ///
   8.340 +    /// \return \c true if a directed cycle exists in the digraph.
   8.341 +    bool findCycleMean() {
   8.342 +      // Initialization and find strongly connected components
   8.343 +      init();
   8.344 +      findComponents();
   8.345 +      
   8.346 +      // Find the minimum cycle mean in the components
   8.347 +      for (int comp = 0; comp < _comp_num; ++comp) {
   8.348 +        if (!initComponent(comp)) continue;
   8.349 +        processRounds();
   8.350 +        updateMinMean();
   8.351 +      }
   8.352 +      return (_cycle_node != INVALID);
   8.353 +    }
   8.354 +
   8.355 +    /// \brief Find a minimum mean directed cycle.
   8.356 +    ///
   8.357 +    /// This function finds a directed cycle of minimum mean cost
   8.358 +    /// in the digraph using the data computed by findCycleMean().
   8.359 +    ///
   8.360 +    /// \return \c true if a directed cycle exists in the digraph.
   8.361 +    ///
   8.362 +    /// \pre \ref findCycleMean() must be called before using this function.
   8.363 +    bool findCycle() {
   8.364 +      if (_cycle_node == INVALID) return false;
   8.365 +      IntNodeMap reached(_gr, -1);
   8.366 +      int r = _data[_cycle_node].size();
   8.367 +      Node u = _cycle_node;
   8.368 +      while (reached[u] < 0) {
   8.369 +        reached[u] = --r;
   8.370 +        u = _gr.source(_data[u][r].pred);
   8.371 +      }
   8.372 +      r = reached[u];
   8.373 +      Arc e = _data[u][r].pred;
   8.374 +      _cycle_path->addFront(e);
   8.375 +      _cycle_cost = _cost[e];
   8.376 +      _cycle_size = 1;
   8.377 +      Node v;
   8.378 +      while ((v = _gr.source(e)) != u) {
   8.379 +        e = _data[v][--r].pred;
   8.380 +        _cycle_path->addFront(e);
   8.381 +        _cycle_cost += _cost[e];
   8.382 +        ++_cycle_size;
   8.383 +      }
   8.384 +      return true;
   8.385 +    }
   8.386 +
   8.387 +    /// @}
   8.388 +
   8.389 +    /// \name Query Functions
   8.390 +    /// The results of the algorithm can be obtained using these
   8.391 +    /// functions.\n
   8.392 +    /// The algorithm should be executed before using them.
   8.393 +
   8.394 +    /// @{
   8.395 +
   8.396 +    /// \brief Return the total cost of the found cycle.
   8.397 +    ///
   8.398 +    /// This function returns the total cost of the found cycle.
   8.399 +    ///
   8.400 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   8.401 +    /// using this function.
   8.402 +    Cost cycleCost() const {
   8.403 +      return static_cast<Cost>(_cycle_cost);
   8.404 +    }
   8.405 +
   8.406 +    /// \brief Return the number of arcs on the found cycle.
   8.407 +    ///
   8.408 +    /// This function returns the number of arcs on the found cycle.
   8.409 +    ///
   8.410 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   8.411 +    /// using this function.
   8.412 +    int cycleSize() const {
   8.413 +      return _cycle_size;
   8.414 +    }
   8.415 +
   8.416 +    /// \brief Return the mean cost of the found cycle.
   8.417 +    ///
   8.418 +    /// This function returns the mean cost of the found cycle.
   8.419 +    ///
   8.420 +    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the
   8.421 +    /// following code.
   8.422 +    /// \code
   8.423 +    ///   return static_cast<double>(alg.cycleCost()) / alg.cycleSize();
   8.424 +    /// \endcode
   8.425 +    ///
   8.426 +    /// \pre \ref run() or \ref findCycleMean() must be called before
   8.427 +    /// using this function.
   8.428 +    double cycleMean() const {
   8.429 +      return static_cast<double>(_cycle_cost) / _cycle_size;
   8.430 +    }
   8.431 +
   8.432 +    /// \brief Return the found cycle.
   8.433 +    ///
   8.434 +    /// This function returns a const reference to the path structure
   8.435 +    /// storing the found cycle.
   8.436 +    ///
   8.437 +    /// \pre \ref run() or \ref findCycle() must be called before using
   8.438 +    /// this function.
   8.439 +    const Path& cycle() const {
   8.440 +      return *_cycle_path;
   8.441 +    }
   8.442 +
   8.443 +    ///@}
   8.444 +
   8.445 +  private:
   8.446 +
   8.447 +    // Initialization
   8.448 +    void init() {
   8.449 +      if (!_cycle_path) {
   8.450 +        _local_path = true;
   8.451 +        _cycle_path = new Path;
   8.452 +      }
   8.453 +      _cycle_path->clear();
   8.454 +      _cycle_cost = 0;
   8.455 +      _cycle_size = 1;
   8.456 +      _cycle_node = INVALID;
   8.457 +      for (NodeIt u(_gr); u != INVALID; ++u)
   8.458 +        _data[u].clear();
   8.459 +    }
   8.460 +
   8.461 +    // Find strongly connected components and initialize _comp_nodes
   8.462 +    // and _out_arcs
   8.463 +    void findComponents() {
   8.464 +      _comp_num = stronglyConnectedComponents(_gr, _comp);
   8.465 +      _comp_nodes.resize(_comp_num);
   8.466 +      if (_comp_num == 1) {
   8.467 +        _comp_nodes[0].clear();
   8.468 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   8.469 +          _comp_nodes[0].push_back(n);
   8.470 +          _out_arcs[n].clear();
   8.471 +          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   8.472 +            _out_arcs[n].push_back(a);
   8.473 +          }
   8.474 +        }
   8.475 +      } else {
   8.476 +        for (int i = 0; i < _comp_num; ++i)
   8.477 +          _comp_nodes[i].clear();
   8.478 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   8.479 +          int k = _comp[n];
   8.480 +          _comp_nodes[k].push_back(n);
   8.481 +          _out_arcs[n].clear();
   8.482 +          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   8.483 +            if (_comp[_gr.target(a)] == k) _out_arcs[n].push_back(a);
   8.484 +          }
   8.485 +        }
   8.486 +      }
   8.487 +    }
   8.488 +
   8.489 +    // Initialize path data for the current component
   8.490 +    bool initComponent(int comp) {
   8.491 +      _nodes = &(_comp_nodes[comp]);
   8.492 +      int n = _nodes->size();
   8.493 +      if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) {
   8.494 +        return false;
   8.495 +      }      
   8.496 +      for (int i = 0; i < n; ++i) {
   8.497 +        _data[(*_nodes)[i]].resize(n + 1, PathData(INF));
   8.498 +      }
   8.499 +      return true;
   8.500 +    }
   8.501 +
   8.502 +    // Process all rounds of computing path data for the current component.
   8.503 +    // _data[v][k] is the cost of a shortest directed walk from the root
   8.504 +    // node to node v containing exactly k arcs.
   8.505 +    void processRounds() {
   8.506 +      Node start = (*_nodes)[0];
   8.507 +      _data[start][0] = PathData(0);
   8.508 +      _process.clear();
   8.509 +      _process.push_back(start);
   8.510 +
   8.511 +      int k, n = _nodes->size();
   8.512 +      for (k = 1; k <= n && int(_process.size()) < n; ++k) {
   8.513 +        processNextBuildRound(k);
   8.514 +      }
   8.515 +      for ( ; k <= n; ++k) {
   8.516 +        processNextFullRound(k);
   8.517 +      }
   8.518 +    }
   8.519 +
   8.520 +    // Process one round and rebuild _process
   8.521 +    void processNextBuildRound(int k) {
   8.522 +      std::vector<Node> next;
   8.523 +      Node u, v;
   8.524 +      Arc e;
   8.525 +      LargeCost d;
   8.526 +      for (int i = 0; i < int(_process.size()); ++i) {
   8.527 +        u = _process[i];
   8.528 +        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   8.529 +          e = _out_arcs[u][j];
   8.530 +          v = _gr.target(e);
   8.531 +          d = _data[u][k-1].dist + _cost[e];
   8.532 +          if (_tolerance.less(d, _data[v][k].dist)) {
   8.533 +            if (_data[v][k].dist == INF) next.push_back(v);
   8.534 +            _data[v][k] = PathData(d, e);
   8.535 +          }
   8.536 +        }
   8.537 +      }
   8.538 +      _process.swap(next);
   8.539 +    }
   8.540 +
   8.541 +    // Process one round using _nodes instead of _process
   8.542 +    void processNextFullRound(int k) {
   8.543 +      Node u, v;
   8.544 +      Arc e;
   8.545 +      LargeCost d;
   8.546 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   8.547 +        u = (*_nodes)[i];
   8.548 +        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   8.549 +          e = _out_arcs[u][j];
   8.550 +          v = _gr.target(e);
   8.551 +          d = _data[u][k-1].dist + _cost[e];
   8.552 +          if (_tolerance.less(d, _data[v][k].dist)) {
   8.553 +            _data[v][k] = PathData(d, e);
   8.554 +          }
   8.555 +        }
   8.556 +      }
   8.557 +    }
   8.558 +
   8.559 +    // Update the minimum cycle mean
   8.560 +    void updateMinMean() {
   8.561 +      int n = _nodes->size();
   8.562 +      for (int i = 0; i < n; ++i) {
   8.563 +        Node u = (*_nodes)[i];
   8.564 +        if (_data[u][n].dist == INF) continue;
   8.565 +        LargeCost cost, max_cost = 0;
   8.566 +        int size, max_size = 1;
   8.567 +        bool found_curr = false;
   8.568 +        for (int k = 0; k < n; ++k) {
   8.569 +          if (_data[u][k].dist == INF) continue;
   8.570 +          cost = _data[u][n].dist - _data[u][k].dist;
   8.571 +          size = n - k;
   8.572 +          if (!found_curr || cost * max_size > max_cost * size) {
   8.573 +            found_curr = true;
   8.574 +            max_cost = cost;
   8.575 +            max_size = size;
   8.576 +          }
   8.577 +        }
   8.578 +        if ( found_curr && (_cycle_node == INVALID ||
   8.579 +             max_cost * _cycle_size < _cycle_cost * max_size) ) {
   8.580 +          _cycle_cost = max_cost;
   8.581 +          _cycle_size = max_size;
   8.582 +          _cycle_node = u;
   8.583 +        }
   8.584 +      }
   8.585 +    }
   8.586 +
   8.587 +  }; //class KarpMmc
   8.588 +
   8.589 +  ///@}
   8.590 +
   8.591 +} //namespace lemon
   8.592 +
   8.593 +#endif //LEMON_KARP_MMC_H
     9.1 --- a/test/min_mean_cycle_test.cc	Mon Mar 08 08:33:41 2010 +0100
     9.2 +++ b/test/min_mean_cycle_test.cc	Sat Mar 13 22:01:38 2010 +0100
     9.3 @@ -25,9 +25,9 @@
     9.4  #include <lemon/concepts/digraph.h>
     9.5  #include <lemon/concept_check.h>
     9.6  
     9.7 -#include <lemon/karp.h>
     9.8 -#include <lemon/hartmann_orlin.h>
     9.9 -#include <lemon/howard.h>
    9.10 +#include <lemon/karp_mmc.h>
    9.11 +#include <lemon/hartmann_orlin_mmc.h>
    9.12 +#include <lemon/howard_mmc.h>
    9.13  
    9.14  #include "test_tools.h"
    9.15  
    9.16 @@ -63,7 +63,7 @@
    9.17  
    9.18                          
    9.19  // Check the interface of an MMC algorithm
    9.20 -template <typename GR, typename Value>
    9.21 +template <typename GR, typename Cost>
    9.22  struct MmcClassConcept
    9.23  {
    9.24    template <typename MMC>
    9.25 @@ -73,30 +73,30 @@
    9.26  
    9.27        typedef typename MMC
    9.28          ::template SetPath<ListPath<GR> >
    9.29 -        ::template SetLargeValue<Value>
    9.30 +        ::template SetLargeCost<Cost>
    9.31          ::Create MmcAlg;
    9.32 -      MmcAlg mmc(me.g, me.length);
    9.33 +      MmcAlg mmc(me.g, me.cost);
    9.34        const MmcAlg& const_mmc = mmc;
    9.35        
    9.36        typename MmcAlg::Tolerance tol = const_mmc.tolerance();
    9.37        mmc.tolerance(tol);
    9.38        
    9.39        b = mmc.cycle(p).run();
    9.40 -      b = mmc.findMinMean();
    9.41 +      b = mmc.findCycleMean();
    9.42        b = mmc.findCycle();
    9.43  
    9.44 -      v = const_mmc.cycleLength();
    9.45 -      i = const_mmc.cycleArcNum();
    9.46 +      v = const_mmc.cycleCost();
    9.47 +      i = const_mmc.cycleSize();
    9.48        d = const_mmc.cycleMean();
    9.49        p = const_mmc.cycle();
    9.50      }
    9.51  
    9.52 -    typedef concepts::ReadMap<typename GR::Arc, Value> LM;
    9.53 +    typedef concepts::ReadMap<typename GR::Arc, Cost> CM;
    9.54    
    9.55      GR g;
    9.56 -    LM length;
    9.57 +    CM cost;
    9.58      ListPath<GR> p;
    9.59 -    Value v;
    9.60 +    Cost v;
    9.61      int i;
    9.62      double d;
    9.63      bool b;
    9.64 @@ -108,13 +108,13 @@
    9.65  void checkMmcAlg(const SmartDigraph& gr,
    9.66                   const SmartDigraph::ArcMap<int>& lm,
    9.67                   const SmartDigraph::ArcMap<int>& cm,
    9.68 -                 int length, int size) {
    9.69 +                 int cost, int size) {
    9.70    MMC alg(gr, lm);
    9.71 -  alg.findMinMean();
    9.72 -  check(alg.cycleMean() == static_cast<double>(length) / size,
    9.73 +  alg.findCycleMean();
    9.74 +  check(alg.cycleMean() == static_cast<double>(cost) / size,
    9.75          "Wrong cycle mean");
    9.76    alg.findCycle();
    9.77 -  check(alg.cycleLength() == length && alg.cycleArcNum() == size,
    9.78 +  check(alg.cycleCost() == cost && alg.cycleSize() == size,
    9.79          "Wrong path");
    9.80    SmartDigraph::ArcMap<int> cycle(gr, 0);
    9.81    for (typename MMC::Path::ArcIt a(alg.cycle()); a != INVALID; ++a) {
    9.82 @@ -148,28 +148,28 @@
    9.83    {
    9.84      typedef concepts::Digraph GR;
    9.85  
    9.86 -    // Karp
    9.87 +    // KarpMmc
    9.88      checkConcept< MmcClassConcept<GR, int>,
    9.89 -                  Karp<GR, concepts::ReadMap<GR::Arc, int> > >();
    9.90 +                  KarpMmc<GR, concepts::ReadMap<GR::Arc, int> > >();
    9.91      checkConcept< MmcClassConcept<GR, float>,
    9.92 -                  Karp<GR, concepts::ReadMap<GR::Arc, float> > >();
    9.93 +                  KarpMmc<GR, concepts::ReadMap<GR::Arc, float> > >();
    9.94      
    9.95 -    // HartmannOrlin
    9.96 +    // HartmannOrlinMmc
    9.97      checkConcept< MmcClassConcept<GR, int>,
    9.98 -                  HartmannOrlin<GR, concepts::ReadMap<GR::Arc, int> > >();
    9.99 +                  HartmannOrlinMmc<GR, concepts::ReadMap<GR::Arc, int> > >();
   9.100      checkConcept< MmcClassConcept<GR, float>,
   9.101 -                  HartmannOrlin<GR, concepts::ReadMap<GR::Arc, float> > >();
   9.102 +                  HartmannOrlinMmc<GR, concepts::ReadMap<GR::Arc, float> > >();
   9.103      
   9.104 -    // Howard
   9.105 +    // HowardMmc
   9.106      checkConcept< MmcClassConcept<GR, int>,
   9.107 -                  Howard<GR, concepts::ReadMap<GR::Arc, int> > >();
   9.108 +                  HowardMmc<GR, concepts::ReadMap<GR::Arc, int> > >();
   9.109      checkConcept< MmcClassConcept<GR, float>,
   9.110 -                  Howard<GR, concepts::ReadMap<GR::Arc, float> > >();
   9.111 +                  HowardMmc<GR, concepts::ReadMap<GR::Arc, float> > >();
   9.112  
   9.113 -    if (IsSameType<Howard<GR, concepts::ReadMap<GR::Arc, int> >::LargeValue,
   9.114 -          long_int>::result == 0) check(false, "Wrong LargeValue type");
   9.115 -    if (IsSameType<Howard<GR, concepts::ReadMap<GR::Arc, float> >::LargeValue,
   9.116 -          double>::result == 0) check(false, "Wrong LargeValue type");
   9.117 +    check((IsSameType<HowardMmc<GR, concepts::ReadMap<GR::Arc, int> >
   9.118 +           ::LargeCost, long_int>::result == 1), "Wrong LargeCost type");
   9.119 +    check((IsSameType<HowardMmc<GR, concepts::ReadMap<GR::Arc, float> >
   9.120 +           ::LargeCost, double>::result == 1), "Wrong LargeCost type");
   9.121    }
   9.122  
   9.123    // Run various tests
   9.124 @@ -194,22 +194,22 @@
   9.125        run();
   9.126  
   9.127      // Karp
   9.128 -    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l1, c1,  6, 3);
   9.129 -    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l2, c2,  5, 2);
   9.130 -    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l3, c3,  0, 1);
   9.131 -    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l4, c4, -1, 1);
   9.132 +    checkMmcAlg<KarpMmc<GR, IntArcMap> >(gr, l1, c1,  6, 3);
   9.133 +    checkMmcAlg<KarpMmc<GR, IntArcMap> >(gr, l2, c2,  5, 2);
   9.134 +    checkMmcAlg<KarpMmc<GR, IntArcMap> >(gr, l3, c3,  0, 1);
   9.135 +    checkMmcAlg<KarpMmc<GR, IntArcMap> >(gr, l4, c4, -1, 1);
   9.136  
   9.137      // HartmannOrlin
   9.138 -    checkMmcAlg<HartmannOrlin<GR, IntArcMap> >(gr, l1, c1,  6, 3);
   9.139 -    checkMmcAlg<HartmannOrlin<GR, IntArcMap> >(gr, l2, c2,  5, 2);
   9.140 -    checkMmcAlg<HartmannOrlin<GR, IntArcMap> >(gr, l3, c3,  0, 1);
   9.141 -    checkMmcAlg<HartmannOrlin<GR, IntArcMap> >(gr, l4, c4, -1, 1);
   9.142 +    checkMmcAlg<HartmannOrlinMmc<GR, IntArcMap> >(gr, l1, c1,  6, 3);
   9.143 +    checkMmcAlg<HartmannOrlinMmc<GR, IntArcMap> >(gr, l2, c2,  5, 2);
   9.144 +    checkMmcAlg<HartmannOrlinMmc<GR, IntArcMap> >(gr, l3, c3,  0, 1);
   9.145 +    checkMmcAlg<HartmannOrlinMmc<GR, IntArcMap> >(gr, l4, c4, -1, 1);
   9.146  
   9.147      // Howard
   9.148 -    checkMmcAlg<Howard<GR, IntArcMap> >(gr, l1, c1,  6, 3);
   9.149 -    checkMmcAlg<Howard<GR, IntArcMap> >(gr, l2, c2,  5, 2);
   9.150 -    checkMmcAlg<Howard<GR, IntArcMap> >(gr, l3, c3,  0, 1);
   9.151 -    checkMmcAlg<Howard<GR, IntArcMap> >(gr, l4, c4, -1, 1);
   9.152 +    checkMmcAlg<HowardMmc<GR, IntArcMap> >(gr, l1, c1,  6, 3);
   9.153 +    checkMmcAlg<HowardMmc<GR, IntArcMap> >(gr, l2, c2,  5, 2);
   9.154 +    checkMmcAlg<HowardMmc<GR, IntArcMap> >(gr, l3, c3,  0, 1);
   9.155 +    checkMmcAlg<HowardMmc<GR, IntArcMap> >(gr, l4, c4, -1, 1);
   9.156    }
   9.157  
   9.158    return 0;