Add Karp algorithm class (#179)
authorPeter Kovacs <kpeter@inf.elte.hu>
Tue, 11 Aug 2009 20:55:40 +0200
changeset 7653b544a9c92db
parent 764 1fac515a59c1
child 766 97744b6dabf8
Add Karp algorithm class (#179)
based on the MinMeanCycle implementation in SVN -r3436.
The interface is reworked to be the same as Howard's interface.
lemon/Makefile.am
lemon/karp.h
test/min_mean_cycle_test.cc
     1.1 --- a/lemon/Makefile.am	Mon Aug 10 14:50:57 2009 +0200
     1.2 +++ b/lemon/Makefile.am	Tue Aug 11 20:55:40 2009 +0200
     1.3 @@ -85,6 +85,7 @@
     1.4  	lemon/grid_graph.h \
     1.5  	lemon/howard.h \
     1.6  	lemon/hypercube_graph.h \
     1.7 +	lemon/karp.h \
     1.8  	lemon/kruskal.h \
     1.9  	lemon/hao_orlin.h \
    1.10  	lemon/lgf_reader.h \
     2.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     2.2 +++ b/lemon/karp.h	Tue Aug 11 20:55:40 2009 +0200
     2.3 @@ -0,0 +1,560 @@
     2.4 +/* -*- C++ -*-
     2.5 + *
     2.6 + * This file is a part of LEMON, a generic C++ optimization library
     2.7 + *
     2.8 + * Copyright (C) 2003-2008
     2.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    2.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    2.11 + *
    2.12 + * Permission to use, modify and distribute this software is granted
    2.13 + * provided that this copyright notice appears in all copies. For
    2.14 + * precise terms see the accompanying LICENSE file.
    2.15 + *
    2.16 + * This software is provided "AS IS" with no warranty of any kind,
    2.17 + * express or implied, and with no claim as to its suitability for any
    2.18 + * purpose.
    2.19 + *
    2.20 + */
    2.21 +
    2.22 +#ifndef LEMON_KARP_H
    2.23 +#define LEMON_KARP_H
    2.24 +
    2.25 +/// \ingroup shortest_path
    2.26 +///
    2.27 +/// \file
    2.28 +/// \brief Karp's algorithm for finding a minimum mean cycle.
    2.29 +
    2.30 +#include <vector>
    2.31 +#include <limits>
    2.32 +#include <lemon/core.h>
    2.33 +#include <lemon/path.h>
    2.34 +#include <lemon/tolerance.h>
    2.35 +#include <lemon/connectivity.h>
    2.36 +
    2.37 +namespace lemon {
    2.38 +
    2.39 +  /// \brief Default traits class of Karp algorithm.
    2.40 +  ///
    2.41 +  /// Default traits class of Karp algorithm.
    2.42 +  /// \tparam GR The type of the digraph.
    2.43 +  /// \tparam LEN The type of the length map.
    2.44 +  /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
    2.45 +#ifdef DOXYGEN
    2.46 +  template <typename GR, typename LEN>
    2.47 +#else
    2.48 +  template <typename GR, typename LEN,
    2.49 +    bool integer = std::numeric_limits<typename LEN::Value>::is_integer>
    2.50 +#endif
    2.51 +  struct KarpDefaultTraits
    2.52 +  {
    2.53 +    /// The type of the digraph
    2.54 +    typedef GR Digraph;
    2.55 +    /// The type of the length map
    2.56 +    typedef LEN LengthMap;
    2.57 +    /// The type of the arc lengths
    2.58 +    typedef typename LengthMap::Value Value;
    2.59 +
    2.60 +    /// \brief The large value type used for internal computations
    2.61 +    ///
    2.62 +    /// The large value type used for internal computations.
    2.63 +    /// It is \c long \c long if the \c Value type is integer,
    2.64 +    /// otherwise it is \c double.
    2.65 +    /// \c Value must be convertible to \c LargeValue.
    2.66 +    typedef double LargeValue;
    2.67 +
    2.68 +    /// The tolerance type used for internal computations
    2.69 +    typedef lemon::Tolerance<LargeValue> Tolerance;
    2.70 +
    2.71 +    /// \brief The path type of the found cycles
    2.72 +    ///
    2.73 +    /// The path type of the found cycles.
    2.74 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
    2.75 +    /// and it must have an \c addBack() function.
    2.76 +    typedef lemon::Path<Digraph> Path;
    2.77 +  };
    2.78 +
    2.79 +  // Default traits class for integer value types
    2.80 +  template <typename GR, typename LEN>
    2.81 +  struct KarpDefaultTraits<GR, LEN, true>
    2.82 +  {
    2.83 +    typedef GR Digraph;
    2.84 +    typedef LEN LengthMap;
    2.85 +    typedef typename LengthMap::Value Value;
    2.86 +#ifdef LEMON_HAVE_LONG_LONG
    2.87 +    typedef long long LargeValue;
    2.88 +#else
    2.89 +    typedef long LargeValue;
    2.90 +#endif
    2.91 +    typedef lemon::Tolerance<LargeValue> Tolerance;
    2.92 +    typedef lemon::Path<Digraph> Path;
    2.93 +  };
    2.94 +
    2.95 +
    2.96 +  /// \addtogroup shortest_path
    2.97 +  /// @{
    2.98 +
    2.99 +  /// \brief Implementation of Karp's algorithm for finding a minimum
   2.100 +  /// mean cycle.
   2.101 +  ///
   2.102 +  /// This class implements Karp's algorithm for finding a directed
   2.103 +  /// cycle of minimum mean length (cost) in a digraph.
   2.104 +  ///
   2.105 +  /// \tparam GR The type of the digraph the algorithm runs on.
   2.106 +  /// \tparam LEN The type of the length map. The default
   2.107 +  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
   2.108 +#ifdef DOXYGEN
   2.109 +  template <typename GR, typename LEN, typename TR>
   2.110 +#else
   2.111 +  template < typename GR,
   2.112 +             typename LEN = typename GR::template ArcMap<int>,
   2.113 +             typename TR = KarpDefaultTraits<GR, LEN> >
   2.114 +#endif
   2.115 +  class Karp
   2.116 +  {
   2.117 +  public:
   2.118 +
   2.119 +    /// The type of the digraph
   2.120 +    typedef typename TR::Digraph Digraph;
   2.121 +    /// The type of the length map
   2.122 +    typedef typename TR::LengthMap LengthMap;
   2.123 +    /// The type of the arc lengths
   2.124 +    typedef typename TR::Value Value;
   2.125 +
   2.126 +    /// \brief The large value type
   2.127 +    ///
   2.128 +    /// The large value type used for internal computations.
   2.129 +    /// Using the \ref KarpDefaultTraits "default traits class",
   2.130 +    /// it is \c long \c long if the \c Value type is integer,
   2.131 +    /// otherwise it is \c double.
   2.132 +    typedef typename TR::LargeValue LargeValue;
   2.133 +
   2.134 +    /// The tolerance type
   2.135 +    typedef typename TR::Tolerance Tolerance;
   2.136 +
   2.137 +    /// \brief The path type of the found cycles
   2.138 +    ///
   2.139 +    /// The path type of the found cycles.
   2.140 +    /// Using the \ref KarpDefaultTraits "default traits class",
   2.141 +    /// it is \ref lemon::Path "Path<Digraph>".
   2.142 +    typedef typename TR::Path Path;
   2.143 +
   2.144 +    /// The \ref KarpDefaultTraits "traits class" of the algorithm
   2.145 +    typedef TR Traits;
   2.146 +
   2.147 +  private:
   2.148 +
   2.149 +    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   2.150 +
   2.151 +    // Data sturcture for path data
   2.152 +    struct PathData
   2.153 +    {
   2.154 +      bool found;
   2.155 +      LargeValue dist;
   2.156 +      Arc pred;
   2.157 +      PathData(bool f = false, LargeValue d = 0, Arc p = INVALID) :
   2.158 +        found(f), dist(d), pred(p) {}
   2.159 +    };
   2.160 +
   2.161 +    typedef typename Digraph::template NodeMap<std::vector<PathData> >
   2.162 +      PathDataNodeMap;
   2.163 +
   2.164 +  private:
   2.165 +
   2.166 +    // The digraph the algorithm runs on
   2.167 +    const Digraph &_gr;
   2.168 +    // The length of the arcs
   2.169 +    const LengthMap &_length;
   2.170 +
   2.171 +    // Data for storing the strongly connected components
   2.172 +    int _comp_num;
   2.173 +    typename Digraph::template NodeMap<int> _comp;
   2.174 +    std::vector<std::vector<Node> > _comp_nodes;
   2.175 +    std::vector<Node>* _nodes;
   2.176 +    typename Digraph::template NodeMap<std::vector<Arc> > _out_arcs;
   2.177 +
   2.178 +    // Data for the found cycle
   2.179 +    LargeValue _cycle_length;
   2.180 +    int _cycle_size;
   2.181 +    Node _cycle_node;
   2.182 +
   2.183 +    Path *_cycle_path;
   2.184 +    bool _local_path;
   2.185 +
   2.186 +    // Node map for storing path data
   2.187 +    PathDataNodeMap _data;
   2.188 +    // The processed nodes in the last round
   2.189 +    std::vector<Node> _process;
   2.190 +
   2.191 +    Tolerance _tolerance;
   2.192 +
   2.193 +  public:
   2.194 +
   2.195 +    /// \name Named Template Parameters
   2.196 +    /// @{
   2.197 +
   2.198 +    template <typename T>
   2.199 +    struct SetLargeValueTraits : public Traits {
   2.200 +      typedef T LargeValue;
   2.201 +      typedef lemon::Tolerance<T> Tolerance;
   2.202 +    };
   2.203 +
   2.204 +    /// \brief \ref named-templ-param "Named parameter" for setting
   2.205 +    /// \c LargeValue type.
   2.206 +    ///
   2.207 +    /// \ref named-templ-param "Named parameter" for setting \c LargeValue
   2.208 +    /// type. It is used for internal computations in the algorithm.
   2.209 +    template <typename T>
   2.210 +    struct SetLargeValue
   2.211 +      : public Karp<GR, LEN, SetLargeValueTraits<T> > {
   2.212 +      typedef Karp<GR, LEN, SetLargeValueTraits<T> > Create;
   2.213 +    };
   2.214 +
   2.215 +    template <typename T>
   2.216 +    struct SetPathTraits : public Traits {
   2.217 +      typedef T Path;
   2.218 +    };
   2.219 +
   2.220 +    /// \brief \ref named-templ-param "Named parameter" for setting
   2.221 +    /// \c %Path type.
   2.222 +    ///
   2.223 +    /// \ref named-templ-param "Named parameter" for setting the \c %Path
   2.224 +    /// type of the found cycles.
   2.225 +    /// It must conform to the \ref lemon::concepts::Path "Path" concept
   2.226 +    /// and it must have an \c addFront() function.
   2.227 +    template <typename T>
   2.228 +    struct SetPath
   2.229 +      : public Karp<GR, LEN, SetPathTraits<T> > {
   2.230 +      typedef Karp<GR, LEN, SetPathTraits<T> > Create;
   2.231 +    };
   2.232 +
   2.233 +    /// @}
   2.234 +
   2.235 +  public:
   2.236 +
   2.237 +    /// \brief Constructor.
   2.238 +    ///
   2.239 +    /// The constructor of the class.
   2.240 +    ///
   2.241 +    /// \param digraph The digraph the algorithm runs on.
   2.242 +    /// \param length The lengths (costs) of the arcs.
   2.243 +    Karp( const Digraph &digraph,
   2.244 +          const LengthMap &length ) :
   2.245 +      _gr(digraph), _length(length), _comp(digraph), _out_arcs(digraph),
   2.246 +      _cycle_length(0), _cycle_size(1), _cycle_node(INVALID),
   2.247 +      _cycle_path(NULL), _local_path(false), _data(digraph)
   2.248 +    {}
   2.249 +
   2.250 +    /// Destructor.
   2.251 +    ~Karp() {
   2.252 +      if (_local_path) delete _cycle_path;
   2.253 +    }
   2.254 +
   2.255 +    /// \brief Set the path structure for storing the found cycle.
   2.256 +    ///
   2.257 +    /// This function sets an external path structure for storing the
   2.258 +    /// found cycle.
   2.259 +    ///
   2.260 +    /// If you don't call this function before calling \ref run() or
   2.261 +    /// \ref findMinMean(), it will allocate a local \ref Path "path"
   2.262 +    /// structure. The destuctor deallocates this automatically
   2.263 +    /// allocated object, of course.
   2.264 +    ///
   2.265 +    /// \note The algorithm calls only the \ref lemon::Path::addFront()
   2.266 +    /// "addFront()" function of the given path structure.
   2.267 +    ///
   2.268 +    /// \return <tt>(*this)</tt>
   2.269 +    Karp& cycle(Path &path) {
   2.270 +      if (_local_path) {
   2.271 +        delete _cycle_path;
   2.272 +        _local_path = false;
   2.273 +      }
   2.274 +      _cycle_path = &path;
   2.275 +      return *this;
   2.276 +    }
   2.277 +
   2.278 +    /// \name Execution control
   2.279 +    /// The simplest way to execute the algorithm is to call the \ref run()
   2.280 +    /// function.\n
   2.281 +    /// If you only need the minimum mean length, you may call
   2.282 +    /// \ref findMinMean().
   2.283 +
   2.284 +    /// @{
   2.285 +
   2.286 +    /// \brief Run the algorithm.
   2.287 +    ///
   2.288 +    /// This function runs the algorithm.
   2.289 +    /// It can be called more than once (e.g. if the underlying digraph
   2.290 +    /// and/or the arc lengths have been modified).
   2.291 +    ///
   2.292 +    /// \return \c true if a directed cycle exists in the digraph.
   2.293 +    ///
   2.294 +    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.
   2.295 +    /// \code
   2.296 +    ///   return mmc.findMinMean() && mmc.findCycle();
   2.297 +    /// \endcode
   2.298 +    bool run() {
   2.299 +      return findMinMean() && findCycle();
   2.300 +    }
   2.301 +
   2.302 +    /// \brief Find the minimum cycle mean.
   2.303 +    ///
   2.304 +    /// This function finds the minimum mean length of the directed
   2.305 +    /// cycles in the digraph.
   2.306 +    ///
   2.307 +    /// \return \c true if a directed cycle exists in the digraph.
   2.308 +    bool findMinMean() {
   2.309 +      // Initialization and find strongly connected components
   2.310 +      init();
   2.311 +      findComponents();
   2.312 +      
   2.313 +      // Find the minimum cycle mean in the components
   2.314 +      for (int comp = 0; comp < _comp_num; ++comp) {
   2.315 +        if (!initComponent(comp)) continue;
   2.316 +        processRounds();
   2.317 +        updateMinMean();
   2.318 +      }
   2.319 +      return (_cycle_node != INVALID);
   2.320 +    }
   2.321 +
   2.322 +    /// \brief Find a minimum mean directed cycle.
   2.323 +    ///
   2.324 +    /// This function finds a directed cycle of minimum mean length
   2.325 +    /// in the digraph using the data computed by findMinMean().
   2.326 +    ///
   2.327 +    /// \return \c true if a directed cycle exists in the digraph.
   2.328 +    ///
   2.329 +    /// \pre \ref findMinMean() must be called before using this function.
   2.330 +    bool findCycle() {
   2.331 +      if (_cycle_node == INVALID) return false;
   2.332 +      IntNodeMap reached(_gr, -1);
   2.333 +      int r = _data[_cycle_node].size();
   2.334 +      Node u = _cycle_node;
   2.335 +      while (reached[u] < 0) {
   2.336 +        reached[u] = --r;
   2.337 +        u = _gr.source(_data[u][r].pred);
   2.338 +      }
   2.339 +      r = reached[u];
   2.340 +      Arc e = _data[u][r].pred;
   2.341 +      _cycle_path->addFront(e);
   2.342 +      _cycle_length = _length[e];
   2.343 +      _cycle_size = 1;
   2.344 +      Node v;
   2.345 +      while ((v = _gr.source(e)) != u) {
   2.346 +        e = _data[v][--r].pred;
   2.347 +        _cycle_path->addFront(e);
   2.348 +        _cycle_length += _length[e];
   2.349 +        ++_cycle_size;
   2.350 +      }
   2.351 +      return true;
   2.352 +    }
   2.353 +
   2.354 +    /// @}
   2.355 +
   2.356 +    /// \name Query Functions
   2.357 +    /// The results of the algorithm can be obtained using these
   2.358 +    /// functions.\n
   2.359 +    /// The algorithm should be executed before using them.
   2.360 +
   2.361 +    /// @{
   2.362 +
   2.363 +    /// \brief Return the total length of the found cycle.
   2.364 +    ///
   2.365 +    /// This function returns the total length of the found cycle.
   2.366 +    ///
   2.367 +    /// \pre \ref run() or \ref findMinMean() must be called before
   2.368 +    /// using this function.
   2.369 +    LargeValue cycleLength() const {
   2.370 +      return _cycle_length;
   2.371 +    }
   2.372 +
   2.373 +    /// \brief Return the number of arcs on the found cycle.
   2.374 +    ///
   2.375 +    /// This function returns the number of arcs on the found cycle.
   2.376 +    ///
   2.377 +    /// \pre \ref run() or \ref findMinMean() must be called before
   2.378 +    /// using this function.
   2.379 +    int cycleArcNum() const {
   2.380 +      return _cycle_size;
   2.381 +    }
   2.382 +
   2.383 +    /// \brief Return the mean length of the found cycle.
   2.384 +    ///
   2.385 +    /// This function returns the mean length of the found cycle.
   2.386 +    ///
   2.387 +    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the
   2.388 +    /// following code.
   2.389 +    /// \code
   2.390 +    ///   return static_cast<double>(alg.cycleLength()) / alg.cycleArcNum();
   2.391 +    /// \endcode
   2.392 +    ///
   2.393 +    /// \pre \ref run() or \ref findMinMean() must be called before
   2.394 +    /// using this function.
   2.395 +    double cycleMean() const {
   2.396 +      return static_cast<double>(_cycle_length) / _cycle_size;
   2.397 +    }
   2.398 +
   2.399 +    /// \brief Return the found cycle.
   2.400 +    ///
   2.401 +    /// This function returns a const reference to the path structure
   2.402 +    /// storing the found cycle.
   2.403 +    ///
   2.404 +    /// \pre \ref run() or \ref findCycle() must be called before using
   2.405 +    /// this function.
   2.406 +    const Path& cycle() const {
   2.407 +      return *_cycle_path;
   2.408 +    }
   2.409 +
   2.410 +    ///@}
   2.411 +
   2.412 +  private:
   2.413 +
   2.414 +    // Initialization
   2.415 +    void init() {
   2.416 +      if (!_cycle_path) {
   2.417 +        _local_path = true;
   2.418 +        _cycle_path = new Path;
   2.419 +      }
   2.420 +      _cycle_path->clear();
   2.421 +      _cycle_length = 0;
   2.422 +      _cycle_size = 1;
   2.423 +      _cycle_node = INVALID;
   2.424 +      for (NodeIt u(_gr); u != INVALID; ++u)
   2.425 +        _data[u].clear();
   2.426 +    }
   2.427 +
   2.428 +    // Find strongly connected components and initialize _comp_nodes
   2.429 +    // and _out_arcs
   2.430 +    void findComponents() {
   2.431 +      _comp_num = stronglyConnectedComponents(_gr, _comp);
   2.432 +      _comp_nodes.resize(_comp_num);
   2.433 +      if (_comp_num == 1) {
   2.434 +        _comp_nodes[0].clear();
   2.435 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   2.436 +          _comp_nodes[0].push_back(n);
   2.437 +          _out_arcs[n].clear();
   2.438 +          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   2.439 +            _out_arcs[n].push_back(a);
   2.440 +          }
   2.441 +        }
   2.442 +      } else {
   2.443 +        for (int i = 0; i < _comp_num; ++i)
   2.444 +          _comp_nodes[i].clear();
   2.445 +        for (NodeIt n(_gr); n != INVALID; ++n) {
   2.446 +          int k = _comp[n];
   2.447 +          _comp_nodes[k].push_back(n);
   2.448 +          _out_arcs[n].clear();
   2.449 +          for (OutArcIt a(_gr, n); a != INVALID; ++a) {
   2.450 +            if (_comp[_gr.target(a)] == k) _out_arcs[n].push_back(a);
   2.451 +          }
   2.452 +        }
   2.453 +      }
   2.454 +    }
   2.455 +
   2.456 +    // Initialize path data for the current component
   2.457 +    bool initComponent(int comp) {
   2.458 +      _nodes = &(_comp_nodes[comp]);
   2.459 +      int n = _nodes->size();
   2.460 +      if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) {
   2.461 +        return false;
   2.462 +      }      
   2.463 +      for (int i = 0; i < n; ++i) {
   2.464 +        _data[(*_nodes)[i]].resize(n + 1);
   2.465 +      }
   2.466 +      return true;
   2.467 +    }
   2.468 +
   2.469 +    // Process all rounds of computing path data for the current component.
   2.470 +    // _data[v][k] is the length of a shortest directed walk from the root
   2.471 +    // node to node v containing exactly k arcs.
   2.472 +    void processRounds() {
   2.473 +      Node start = (*_nodes)[0];
   2.474 +      _data[start][0] = PathData(true, 0);
   2.475 +      _process.clear();
   2.476 +      _process.push_back(start);
   2.477 +
   2.478 +      int k, n = _nodes->size();
   2.479 +      for (k = 1; k <= n && int(_process.size()) < n; ++k) {
   2.480 +        processNextBuildRound(k);
   2.481 +      }
   2.482 +      for ( ; k <= n; ++k) {
   2.483 +        processNextFullRound(k);
   2.484 +      }
   2.485 +    }
   2.486 +
   2.487 +    // Process one round and rebuild _process
   2.488 +    void processNextBuildRound(int k) {
   2.489 +      std::vector<Node> next;
   2.490 +      Node u, v;
   2.491 +      Arc e;
   2.492 +      LargeValue d;
   2.493 +      for (int i = 0; i < int(_process.size()); ++i) {
   2.494 +        u = _process[i];
   2.495 +        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   2.496 +          e = _out_arcs[u][j];
   2.497 +          v = _gr.target(e);
   2.498 +          d = _data[u][k-1].dist + _length[e];
   2.499 +          if (!_data[v][k].found) {
   2.500 +            next.push_back(v);
   2.501 +            _data[v][k] = PathData(true, _data[u][k-1].dist + _length[e], e);
   2.502 +          }
   2.503 +          else if (_tolerance.less(d, _data[v][k].dist)) {
   2.504 +            _data[v][k] = PathData(true, d, e);
   2.505 +          }
   2.506 +        }
   2.507 +      }
   2.508 +      _process.swap(next);
   2.509 +    }
   2.510 +
   2.511 +    // Process one round using _nodes instead of _process
   2.512 +    void processNextFullRound(int k) {
   2.513 +      Node u, v;
   2.514 +      Arc e;
   2.515 +      LargeValue d;
   2.516 +      for (int i = 0; i < int(_nodes->size()); ++i) {
   2.517 +        u = (*_nodes)[i];
   2.518 +        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {
   2.519 +          e = _out_arcs[u][j];
   2.520 +          v = _gr.target(e);
   2.521 +          d = _data[u][k-1].dist + _length[e];
   2.522 +          if (!_data[v][k].found || _tolerance.less(d, _data[v][k].dist)) {
   2.523 +            _data[v][k] = PathData(true, d, e);
   2.524 +          }
   2.525 +        }
   2.526 +      }
   2.527 +    }
   2.528 +
   2.529 +    // Update the minimum cycle mean
   2.530 +    void updateMinMean() {
   2.531 +      int n = _nodes->size();
   2.532 +      for (int i = 0; i < n; ++i) {
   2.533 +        Node u = (*_nodes)[i];
   2.534 +        if (!_data[u][n].found) continue;
   2.535 +        LargeValue length, max_length = 0;
   2.536 +        int size, max_size = 1;
   2.537 +        bool found_curr = false;
   2.538 +        for (int k = 0; k < n; ++k) {
   2.539 +          if (!_data[u][k].found) continue;
   2.540 +          length = _data[u][n].dist - _data[u][k].dist;
   2.541 +          size = n - k;
   2.542 +          if (!found_curr || length * max_size > max_length * size) {
   2.543 +            found_curr = true;
   2.544 +            max_length = length;
   2.545 +            max_size = size;
   2.546 +          }
   2.547 +        }
   2.548 +        if ( found_curr && (_cycle_node == INVALID ||
   2.549 +             max_length * _cycle_size < _cycle_length * max_size) ) {
   2.550 +          _cycle_length = max_length;
   2.551 +          _cycle_size = max_size;
   2.552 +          _cycle_node = u;
   2.553 +        }
   2.554 +      }
   2.555 +    }
   2.556 +
   2.557 +  }; //class Karp
   2.558 +
   2.559 +  ///@}
   2.560 +
   2.561 +} //namespace lemon
   2.562 +
   2.563 +#endif //LEMON_KARP_H
     3.1 --- a/test/min_mean_cycle_test.cc	Mon Aug 10 14:50:57 2009 +0200
     3.2 +++ b/test/min_mean_cycle_test.cc	Tue Aug 11 20:55:40 2009 +0200
     3.3 @@ -21,11 +21,13 @@
     3.4  
     3.5  #include <lemon/smart_graph.h>
     3.6  #include <lemon/lgf_reader.h>
     3.7 -#include <lemon/howard.h>
     3.8  #include <lemon/path.h>
     3.9  #include <lemon/concepts/digraph.h>
    3.10  #include <lemon/concept_check.h>
    3.11  
    3.12 +#include <lemon/karp.h>
    3.13 +#include <lemon/howard.h>
    3.14 +
    3.15  #include "test_tools.h"
    3.16  
    3.17  using namespace lemon;
    3.18 @@ -141,16 +143,23 @@
    3.19    // Check the interface
    3.20    {
    3.21      typedef concepts::Digraph GR;
    3.22 -    typedef Howard<GR, concepts::ReadMap<GR::Arc, int> > IntMmcAlg;
    3.23 -    typedef Howard<GR, concepts::ReadMap<GR::Arc, float> > FloatMmcAlg;
    3.24 +
    3.25 +    // Karp
    3.26 +    checkConcept< MmcClassConcept<GR, int>,
    3.27 +                  Karp<GR, concepts::ReadMap<GR::Arc, int> > >();
    3.28 +    checkConcept< MmcClassConcept<GR, float>,
    3.29 +                  Karp<GR, concepts::ReadMap<GR::Arc, float> > >();
    3.30      
    3.31 -    checkConcept<MmcClassConcept<GR, int>, IntMmcAlg>();
    3.32 -    checkConcept<MmcClassConcept<GR, float>, FloatMmcAlg>();
    3.33 -  
    3.34 -    if (IsSameType<IntMmcAlg::LargeValue, long_int>::result == 0)
    3.35 -      check(false, "Wrong LargeValue type");
    3.36 -    if (IsSameType<FloatMmcAlg::LargeValue, double>::result == 0)
    3.37 -      check(false, "Wrong LargeValue type");
    3.38 +    // Howard
    3.39 +    checkConcept< MmcClassConcept<GR, int>,
    3.40 +                  Howard<GR, concepts::ReadMap<GR::Arc, int> > >();
    3.41 +    checkConcept< MmcClassConcept<GR, float>,
    3.42 +                  Howard<GR, concepts::ReadMap<GR::Arc, float> > >();
    3.43 +
    3.44 +    if (IsSameType<Howard<GR, concepts::ReadMap<GR::Arc, int> >::LargeValue,
    3.45 +          long_int>::result == 0) check(false, "Wrong LargeValue type");
    3.46 +    if (IsSameType<Howard<GR, concepts::ReadMap<GR::Arc, float> >::LargeValue,
    3.47 +          double>::result == 0) check(false, "Wrong LargeValue type");
    3.48    }
    3.49  
    3.50    // Run various tests
    3.51 @@ -174,6 +183,13 @@
    3.52        arcMap("c4", c4).
    3.53        run();
    3.54  
    3.55 +    // Karp
    3.56 +    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l1, c1,  6, 3);
    3.57 +    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l2, c2,  5, 2);
    3.58 +    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l3, c3,  0, 1);
    3.59 +    checkMmcAlg<Karp<GR, IntArcMap> >(gr, l4, c4, -1, 1);
    3.60 +
    3.61 +    // Howard
    3.62      checkMmcAlg<Howard<GR, IntArcMap> >(gr, l1, c1,  6, 3);
    3.63      checkMmcAlg<Howard<GR, IntArcMap> >(gr, l2, c2,  5, 2);
    3.64      checkMmcAlg<Howard<GR, IntArcMap> >(gr, l3, c3,  0, 1);