COIN-OR::LEMON - Graph Library

source: lemon-0.x/lemon/min_mean_cycle.h @ 2583:7216b6a52ab9

Last change on this file since 2583:7216b6a52ab9 was 2583:7216b6a52ab9, checked in by Peter Kovacs, 11 years ago

Small fixes and doc improvements in MinMeanCycle?.

File size: 13.9 KB
Line 
1/* -*- C++ -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library
4 *
5 * Copyright (C) 2003-2008
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 *
9 * Permission to use, modify and distribute this software is granted
10 * provided that this copyright notice appears in all copies. For
11 * precise terms see the accompanying LICENSE file.
12 *
13 * This software is provided "AS IS" with no warranty of any kind,
14 * express or implied, and with no claim as to its suitability for any
15 * purpose.
16 *
17 */
18
19#ifndef LEMON_MIN_MEAN_CYCLE_H
20#define LEMON_MIN_MEAN_CYCLE_H
21
22/// \ingroup shortest_path
23///
24/// \file
25/// \brief Howard's algorithm for finding a minimum mean directed cycle.
26
27#include <vector>
28#include <lemon/graph_utils.h>
29#include <lemon/path.h>
30#include <lemon/tolerance.h>
31#include <lemon/topology.h>
32
33namespace lemon {
34
35  /// \addtogroup shortest_path
36  /// @{
37
38  /// \brief Implementation of Howard's algorithm for finding a minimum
39  /// mean directed cycle.
40  ///
41  /// \ref MinMeanCycle implements Howard's algorithm for finding a
42  /// minimum mean directed cycle.
43  ///
44  /// \tparam Graph The directed graph type the algorithm runs on.
45  /// \tparam LengthMap The type of the length (cost) map.
46  ///
47  /// \warning \c LengthMap::Value must be convertible to \c double.
48  ///
49  /// \author Peter Kovacs
50
51  template < typename Graph,
52             typename LengthMap = typename Graph::template EdgeMap<int> >
53  class MinMeanCycle
54  {
55    GRAPH_TYPEDEFS(typename Graph);
56
57    typedef typename LengthMap::Value Length;
58    typedef Path<Graph> Path;
59
60  private:
61
62    // The directed graph the algorithm runs on
63    const Graph &_graph;
64    // The length of the edges
65    const LengthMap &_length;
66
67    // The total length of the found cycle
68    Length _cycle_length;
69    // The number of edges on the found cycle
70    int _cycle_size;
71    // The found cycle
72    Path *_cycle_path;
73
74    bool _local_path;
75    bool _cycle_found;
76    Node _cycle_node;
77
78    typename Graph::template NodeMap<bool> _reached;
79    typename Graph::template NodeMap<double> _dist;
80    typename Graph::template NodeMap<Edge> _policy;
81
82    typename Graph::template NodeMap<int> _component;
83    int _component_num;
84
85    std::vector<Node> _nodes;
86    std::vector<Edge> _edges;
87    Tolerance<double> _tolerance;
88
89  public:
90
91    /// \brief The constructor of the class.
92    ///
93    /// The constructor of the class.
94    ///
95    /// \param graph The directed graph the algorithm runs on.
96    /// \param length The length (cost) of the edges.
97    MinMeanCycle( const Graph &graph,
98                  const LengthMap &length ) :
99      _graph(graph), _length(length), _cycle_length(0), _cycle_size(-1),
100      _cycle_path(NULL), _local_path(false), _reached(graph),
101      _dist(graph), _policy(graph), _component(graph)
102    {}
103
104    /// The destructor of the class.
105    ~MinMeanCycle() {
106      if (_local_path) delete _cycle_path;
107    }
108
109    /// \brief Sets the \ref Path "path" structure for storing the found
110    /// cycle.
111    ///
112    /// Sets an external \ref Path "path" structure for storing the
113    /// found cycle.
114    ///
115    /// If you don't call this function before calling \ref run() or
116    /// \ref init(), it will allocate a local \ref Path "path"
117    /// structure.
118    /// The destuctor deallocates this automatically allocated map,
119    /// of course.
120    ///
121    /// \note The algorithm calls only the \ref lemon::Path::addBack()
122    /// "addBack()" function of the given \ref Path "path" structure.
123    ///
124    /// \return <tt>(*this)</tt>
125    ///
126    /// \sa cycle()
127    MinMeanCycle& cyclePath(Path &path) {
128      if (_local_path) {
129        delete _cycle_path;
130        _local_path = false;
131      }
132      _cycle_path = &path;
133      return *this;
134    }
135
136    /// \name Execution control
137    /// The simplest way to execute the algorithm is to call run().
138    /// \n
139    /// If you only need the minimum mean value, you may call init()
140    /// and findMinMean().
141    /// \n
142    /// If you would like to run the algorithm again (e.g. the
143    /// underlaying graph and/or the edge costs were modified), you may
144    /// not create a new instance of the class, rather call reset(),
145    /// findMinMean(), and findCycle() instead.
146
147    /// @{
148
149    /// \brief Runs the algorithm.
150    ///
151    /// Runs the algorithm.
152    ///
153    /// \return Returns \c true if a directed cycle exists in the graph.
154    ///
155    /// \note Apart from the return value, <tt>mmc.run()</tt> is just a
156    /// shortcut of the following code.
157    /// \code
158    ///   mmc.init();
159    ///   mmc.findMinMean();
160    ///   mmc.findCycle();
161    /// \endcode
162    bool run() {
163      init();
164      return findMinMean() && findCycle();
165    }
166
167    /// \brief Initializes the internal data structures.
168    ///
169    /// Initializes the internal data structures.
170    ///
171    /// \sa reset()
172    void init() {
173      _tolerance.epsilon(1e-6);
174      if (!_cycle_path) {
175        _local_path = true;
176        _cycle_path = new Path;
177      }
178      _cycle_found = false;
179      _component_num = stronglyConnectedComponents(_graph, _component);
180    }
181
182    /// \brief Resets the internal data structures.
183    ///
184    /// Resets the internal data structures so that \ref findMinMean()
185    /// and \ref findCycle() can be called again (e.g. when the
186    /// underlaying graph has been modified).
187    ///
188    /// \sa init()
189    void reset() {
190      if (_cycle_path) _cycle_path->clear();
191      _cycle_found = false;
192      _component_num = stronglyConnectedComponents(_graph, _component);
193    }
194
195    /// \brief Finds the minimum cycle mean length in the graph.
196    ///
197    /// Computes all the required data and finds the minimum cycle mean
198    /// length in the graph.
199    ///
200    /// \return Returns \c true if a directed cycle exists in the graph.
201    ///
202    /// \pre \ref init() must be called before using this function.
203    bool findMinMean() {
204      // Finding the minimum mean cycle in the components
205      for (int comp = 0; comp < _component_num; ++comp) {
206        if (!initCurrentComponent(comp)) continue;
207        while (true) {
208          if (!findPolicyCycles()) break;
209          contractPolicyGraph(comp);
210          if (!computeNodeDistances(comp)) break;
211        }
212      }
213      return _cycle_found;
214    }
215
216    /// \brief Finds a critical (minimum mean) directed cycle.
217    ///
218    /// Finds a critical (minimum mean) directed cycle using the data
219    /// computed in the \ref findMinMean() function.
220    ///
221    /// \return Returns \c true if a directed cycle exists in the graph.
222    ///
223    /// \pre \ref init() and \ref findMinMean() must be called before
224    /// using this function.
225    bool findCycle() {
226      if (!_cycle_found) return false;
227      _cycle_path->addBack(_policy[_cycle_node]);
228      for ( Node v = _cycle_node;
229            (v = _graph.target(_policy[v])) != _cycle_node; ) {
230        _cycle_path->addBack(_policy[v]);
231      }
232      return true;
233    }
234   
235    /// @}
236
237    /// \name Query Functions
238    /// The result of the algorithm can be obtained using these
239    /// functions.
240    /// \n run() must be called before using them.
241
242    /// @{
243   
244    /// \brief Returns the total length of the found cycle.
245    ///
246    /// Returns the total length of the found cycle.
247    ///
248    /// \pre \ref run() or \ref findMinMean() must be called before
249    /// using this function.
250    Length cycleLength() const {
251      return _cycle_length;
252    }
253
254    /// \brief Returns the number of edges on the found cycle.
255    ///
256    /// Returns the number of edges on the found cycle.
257    ///
258    /// \pre \ref run() or \ref findMinMean() must be called before
259    /// using this function.
260    int cycleEdgeNum() const {
261      return _cycle_size;
262    }
263
264    /// \brief Returns the mean length of the found cycle.
265    ///
266    /// Returns the mean length of the found cycle.
267    ///
268    /// \pre \ref run() or \ref findMinMean() must be called before
269    /// using this function.
270    ///
271    /// \note <tt>mmc.cycleMean()</tt> is just a shortcut of the
272    /// following code.
273    /// \code
274    ///   return double(mmc.cycleLength()) / mmc.cycleEdgeNum();
275    /// \endcode
276    double cycleMean() const {
277      return double(_cycle_length) / _cycle_size;
278    }
279
280    /// \brief Returns a const reference to the \ref Path "path"
281    /// structure storing the found cycle.
282    ///
283    /// Returns a const reference to the \ref Path "path"
284    /// structure storing the found cycle.
285    ///
286    /// \pre \ref run() or \ref findCycle() must be called before using
287    /// this function.
288    ///
289    /// \sa cyclePath()
290    const Path& cycle() const {
291      return *_cycle_path;
292    }
293   
294    ///@}
295   
296  private:
297
298    // Initializes the internal data structures for the current strongly
299    // connected component and creating the policy graph.
300    // The policy graph can be represented by the _policy map because
301    // the out degree of every node is 1.
302    bool initCurrentComponent(int comp) {
303      // Finding the nodes of the current component
304      _nodes.clear();
305      for (NodeIt n(_graph); n != INVALID; ++n) {
306        if (_component[n] == comp) _nodes.push_back(n);
307      }
308      if (_nodes.size() <= 1) return false;
309      // Finding the edges of the current component
310      _edges.clear();
311      for (EdgeIt e(_graph); e != INVALID; ++e) {
312        if ( _component[_graph.source(e)] == comp &&
313             _component[_graph.target(e)] == comp )
314          _edges.push_back(e);
315      }
316      // Initializing _reached, _dist, _policy maps
317      for (int i = 0; i < int(_nodes.size()); ++i) {
318        _reached[_nodes[i]] = false;
319        _policy[_nodes[i]] = INVALID;
320      }
321      Node u; Edge e;
322      for (int j = 0; j < int(_edges.size()); ++j) {
323        e = _edges[j];
324        u = _graph.source(e);
325        if (!_reached[u] || _length[e] < _dist[u]) {
326          _dist[u] = _length[e];
327          _policy[u] = e;
328          _reached[u] = true;
329        }
330      }
331      return true;
332    }
333
334    // Finds all cycles in the policy graph.
335    // Sets _cycle_found to true if a cycle is found and sets
336    // _cycle_length, _cycle_size, _cycle_node to represent the minimum
337    // mean cycle in the policy graph.
338    bool findPolicyCycles() {
339      typename Graph::template NodeMap<int> level(_graph, -1);
340      bool curr_cycle_found = false;
341      Length clength;
342      int csize;
343      int path_cnt = 0;
344      Node u, v;
345      // Searching for cycles
346      for (int i = 0; i < int(_nodes.size()); ++i) {
347        if (level[_nodes[i]] < 0) {
348          u = _nodes[i];
349          level[u] = path_cnt;
350          while (level[u = _graph.target(_policy[u])] < 0)
351            level[u] = path_cnt;
352          if (level[u] == path_cnt) {
353            // A cycle is found
354            curr_cycle_found = true;
355            clength = _length[_policy[u]];
356            csize = 1;
357            for (v = u; (v = _graph.target(_policy[v])) != u; ) {
358              clength += _length[_policy[v]];
359              ++csize;
360            }
361            if ( !_cycle_found ||
362                 clength * _cycle_size < _cycle_length * csize ) {
363              _cycle_found = true;
364              _cycle_length = clength;
365              _cycle_size = csize;
366              _cycle_node = u;
367            }
368          }
369          ++path_cnt;
370        }
371      }
372      return curr_cycle_found;
373    }
374
375    // Contracts the policy graph to be connected by cutting all cycles
376    // except for the main cycle (i.e. the minimum mean cycle).
377    void contractPolicyGraph(int comp) {
378      // Finding the component of the main cycle using
379      // reverse BFS search
380      typename Graph::template NodeMap<int> found(_graph, false);
381      std::deque<Node> queue;
382      queue.push_back(_cycle_node);
383      found[_cycle_node] = true;
384      Node u, v;
385      while (!queue.empty()) {
386        v = queue.front(); queue.pop_front();
387        for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
388          u = _graph.source(e);
389          if (_component[u] == comp && !found[u] && _policy[u] == e) {
390            found[u] = true;
391            queue.push_back(u);
392          }
393        }
394      }
395      // Connecting all other nodes to this component using
396      // reverse BFS search
397      queue.clear();
398      for (int i = 0; i < int(_nodes.size()); ++i)
399        if (found[_nodes[i]]) queue.push_back(_nodes[i]);
400      int found_cnt = queue.size();
401      while (found_cnt < int(_nodes.size()) && !queue.empty()) {
402        v = queue.front(); queue.pop_front();
403        for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
404          u = _graph.source(e);
405          if (_component[u] == comp && !found[u]) {
406            found[u] = true;
407            ++found_cnt;
408            _policy[u] = e;
409            queue.push_back(u);
410          }
411        }
412      }
413    }
414
415    // Computes node distances in the policy graph and updates the
416    // policy graph if the node distances can be improved.
417    bool computeNodeDistances(int comp) {
418      // Computing node distances using reverse BFS search
419      double cycle_mean = double(_cycle_length) / _cycle_size;
420      typename Graph::template NodeMap<int> found(_graph, false);
421      std::deque<Node> queue;
422      queue.push_back(_cycle_node);
423      found[_cycle_node] = true;
424      _dist[_cycle_node] = 0;
425      Node u, v;
426      while (!queue.empty()) {
427        v = queue.front(); queue.pop_front();
428        for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
429          u = _graph.source(e);
430          if (_component[u] == comp && !found[u] && _policy[u] == e) {
431            found[u] = true;
432            _dist[u] = _dist[v] + _length[e] - cycle_mean;
433            queue.push_back(u);
434          }
435        }
436      }
437      // Improving node distances
438      bool improved = false;
439      for (int j = 0; j < int(_edges.size()); ++j) {
440        Edge e = _edges[j];
441        u = _graph.source(e); v = _graph.target(e);
442        double delta = _dist[v] + _length[e] - cycle_mean;
443        if (_tolerance.less(delta, _dist[u])) {
444          improved = true;
445          _dist[u] = delta;
446          _policy[u] = e;
447        }
448      }
449      return improved;
450    }
451
452  }; //class MinMeanCycle
453
454  ///@}
455
456} //namespace lemon
457
458#endif //LEMON_MIN_MEAN_CYCLE_H
Note: See TracBrowser for help on using the repository browser.