[Lemon-commits] kpeter: r3437 - lemon/trunk/lemon
Lemon SVN
svn at lemon.cs.elte.hu
Sun Jan 13 11:26:56 CET 2008
Author: kpeter
Date: Sun Jan 13 11:26:55 2008
New Revision: 3437
Modified:
lemon/trunk/lemon/min_mean_cycle.h
Log:
Reimplemented MinMeanCycle to be much more efficient.
The new version implements Howard's algorithm instead of Karp's algorithm and
it is at least 10-20 times faster on all the 40-50 random graphs we have tested.
Modified: lemon/trunk/lemon/min_mean_cycle.h
==============================================================================
--- lemon/trunk/lemon/min_mean_cycle.h (original)
+++ lemon/trunk/lemon/min_mean_cycle.h Sun Jan 13 11:26:55 2008
@@ -22,224 +22,242 @@
/// \ingroup shortest_path
///
/// \file
-/// \brief Karp's algorithm for finding a minimum mean (directed) cycle.
+/// \brief Howard's algorithm for finding a minimum mean directed cycle.
#include <vector>
#include <lemon/graph_utils.h>
#include <lemon/topology.h>
#include <lemon/path.h>
+#include <lemon/tolerance.h>
namespace lemon {
/// \addtogroup shortest_path
/// @{
- /// \brief Implementation of Karp's algorithm for finding a
- /// minimum mean (directed) cycle.
+ /// \brief Implementation of Howard's algorithm for finding a minimum
+ /// mean directed cycle.
///
- /// The \ref lemon::MinMeanCycle "MinMeanCycle" implements Karp's
- /// algorithm for finding a minimum mean (directed) cycle.
+ /// \ref MinMeanCycle implements Howard's algorithm for finding a
+ /// minimum mean directed cycle.
///
/// \param Graph The directed graph type the algorithm runs on.
/// \param LengthMap The type of the length (cost) map.
///
+ /// \warning \c LengthMap::Value must be convertible to \c double.
+ ///
/// \author Peter Kovacs
- template <typename Graph,
- typename LengthMap = typename Graph::template EdgeMap<int> >
+ template < typename Graph,
+ typename LengthMap = typename Graph::template EdgeMap<int> >
class MinMeanCycle
{
- typedef typename Graph::Node Node;
- typedef typename Graph::NodeIt NodeIt;
- typedef typename Graph::Edge Edge;
- typedef typename Graph::EdgeIt EdgeIt;
- typedef typename Graph::OutEdgeIt OutEdgeIt;
+ GRAPH_TYPEDEFS(typename Graph);
typedef typename LengthMap::Value Length;
-
- typedef typename Graph::template NodeMap<int> IntNodeMap;
- typedef typename Graph::template NodeMap<Edge> PredNodeMap;
typedef Path<Graph> Path;
- typedef std::vector<Node> NodeVector;
protected:
- /// \brief Data structure for path data.
- struct PathData
- {
- bool found;
- Length dist;
- Edge pred;
- PathData(bool _found = false, Length _dist = 0) :
- found(_found), dist(_dist), pred(INVALID) {}
- PathData(bool _found, Length _dist, Edge _pred) :
- found(_found), dist(_dist), pred(_pred) {}
- };
-
- private:
-
- typedef typename Graph::template NodeMap<std::vector<PathData> >
- PathDataNodeMap;
-
- protected:
+ typename Graph::template NodeMap<bool> _reached;
+ typename Graph::template NodeMap<double> _dist;
+ typename Graph::template NodeMap<Edge> _policy;
+
+ /// The directed graph the algorithm runs on.
+ const Graph &_graph;
+ /// The length of the edges.
+ const LengthMap &_length;
+
+ /// The total length of the found cycle.
+ Length _cycle_length;
+ /// The number of edges on the found cycle.
+ int _cycle_size;
+ /// The found cycle.
+ Path *_cycle_path;
+
+ bool _local_path;
+ bool _cycle_found;
+ Node _cycle_node;
+
+ typename Graph::template NodeMap<int> _component;
+ int _component_num;
+
+ std::vector<Node> _nodes;
+ std::vector<Edge> _edges;
+ Tolerance<double> _tolerance;
- /// \brief Node map for storing path data.
- ///
- /// Node map for storing path data of all nodes in the current
- /// component. dmap[v][k] is the length of a shortest directed walk
- /// to node v from the starting node containing exactly k edges.
- PathDataNodeMap dmap;
-
- /// \brief The directed graph the algorithm runs on.
- const Graph &graph;
- /// \brief The length of the edges.
- const LengthMap &length;
-
- /// \brief The total length of the found cycle.
- Length cycle_length;
- /// \brief The number of edges in the found cycle.
- int cycle_size;
- /// \brief A node for obtaining a minimum mean cycle.
- Node cycle_node;
-
- /// \brief The found cycle.
- Path *cycle_path;
- /// \brief The algorithm uses local \ref lemon::Path "Path"
- /// structure to store the found cycle.
- bool local_path;
-
- /// \brief Node map for identifying strongly connected components.
- IntNodeMap comp;
- /// \brief The number of strongly connected components.
- int comp_num;
- /// \brief Counter for identifying the current component.
- int comp_cnt;
- /// \brief Nodes of the current component.
- NodeVector nodes;
- /// \brief The processed nodes in the last round.
- NodeVector process;
-
- public :
+ public:
/// \brief The constructor of the class.
///
/// The constructor of the class.
///
- /// \param _graph The directed graph the algorithm runs on.
- /// \param _length The length (cost) of the edges.
- MinMeanCycle( const Graph &_graph,
- const LengthMap &_length ) :
- graph(_graph), length(_length), dmap(_graph), comp(_graph),
- cycle_length(0), cycle_size(-1), cycle_node(INVALID),
- cycle_path(NULL), local_path(false)
+ /// \param graph The directed graph the algorithm runs on.
+ /// \param length The length (cost) of the edges.
+ MinMeanCycle( const Graph &graph,
+ const LengthMap &length ) :
+ _graph(graph), _length(length), _dist(graph), _reached(graph),
+ _policy(graph), _component(graph), _cycle_length(0),
+ _cycle_size(-1), _cycle_path(NULL), _local_path(false)
{ }
- /// \brief The destructor of the class.
+ /// The destructor of the class.
~MinMeanCycle() {
- if (local_path) delete cycle_path;
+ if (_local_path) delete _cycle_path;
}
protected:
- /// \brief Initializes the internal data structures for the current
- /// component.
- void initCurrent() {
- nodes.clear();
+ // Initializes the internal data structures for the current strongly
+ // connected component and creating the policy graph.
+ // The policy graph can be represented by the _policy map because
+ // the out degree of every node is 1.
+ bool initCurrentComponent(int comp) {
// Finding the nodes of the current component
- for (NodeIt v(graph); v != INVALID; ++v) {
- if (comp[v] == comp_cnt) nodes.push_back(v);
+ _nodes.clear();
+ for (NodeIt n(_graph); n != INVALID; ++n) {
+ if (_component[n] == comp) _nodes.push_back(n);
+ }
+ if (_nodes.size() <= 1) return false;
+ // Finding the edges of the current component
+ _edges.clear();
+ for (EdgeIt e(_graph); e != INVALID; ++e) {
+ if ( _component[_graph.source(e)] == comp &&
+ _component[_graph.target(e)] == comp )
+ _edges.push_back(e);
+ }
+ // Initializing _reached, _dist, _policy maps
+ for (int i = 0; i < _nodes.size(); ++i) {
+ _reached[_nodes[i]] = false;
+ _policy[_nodes[i]] = INVALID;
+ }
+ Node u; Edge e;
+ for (int j = 0; j < _edges.size(); ++j) {
+ e = _edges[j];
+ u = _graph.source(e);
+ if (!_reached[u] || _length[e] < _dist[u]) {
+ _dist[u] = _length[e];
+ _policy[u] = e;
+ _reached[u] = true;
+ }
}
- // Creating vectors for all nodes
- int n = nodes.size();
- for (int i = 0; i < n; ++i) {
- dmap[nodes[i]].resize(n + 1);
- }
- }
-
- /// \brief Processes all rounds of computing required path data for
- /// the current component.
- void processRounds() {
- dmap[nodes[0]][0] = PathData(true, 0);
- process.clear();
- // Processing the first round
- for (OutEdgeIt e(graph, nodes[0]); e != INVALID; ++e) {
- Node v = graph.target(e);
- if (comp[v] != comp_cnt || v == nodes[0]) continue;
- dmap[v][1] = PathData(true, length[e], e);
- process.push_back(v);
- }
- // Processing other rounds
- int n = nodes.size(), k;
- for (k = 2; k <= n && process.size() < n; ++k)
- processNextBuildRound(k);
- for ( ; k <= n; ++k)
- processNextFullRound(k);
- }
-
- /// \brief Processes one round of computing required path data and
- /// rebuilds \ref process vector.
- void processNextBuildRound(int k) {
- NodeVector next;
- for (int i = 0; i < process.size(); ++i) {
- for (OutEdgeIt e(graph, process[i]); e != INVALID; ++e) {
- Node v = graph.target(e);
- if (comp[v] != comp_cnt) continue;
- if (!dmap[v][k].found) {
- next.push_back(v);
- dmap[v][k] = PathData(true, dmap[process[i]][k-1].dist + length[e], e);
- }
- else if (dmap[process[i]][k-1].dist + length[e] < dmap[v][k].dist) {
- dmap[v][k] = PathData(true, dmap[process[i]][k-1].dist + length[e], e);
- }
- }
- }
- process.swap(next);
- }
-
- /// \brief Processes one round of computing required path data
- /// using \ref nodes vector instead of \ref process vector.
- void processNextFullRound(int k) {
- for (int i = 0; i < nodes.size(); ++i) {
- for (OutEdgeIt e(graph, nodes[i]); e != INVALID; ++e) {
- Node v = graph.target(e);
- if (comp[v] != comp_cnt) continue;
- if ( !dmap[v][k].found ||
- dmap[nodes[i]][k-1].dist + length[e] < dmap[v][k].dist ) {
- dmap[v][k] = PathData(true, dmap[nodes[i]][k-1].dist + length[e], e);
- }
- }
- }
- }
-
- /// \brief Finds the minimum cycle mean value in the current
- /// component.
- bool findCurrentMin(Length &min_length, int &min_size, Node &min_node) {
- bool found_min = false;
- for (int i = 0; i < nodes.size(); ++i) {
- int n = nodes.size();
- if (!dmap[nodes[i]][n].found) continue;
- Length len;
- int size;
- bool found_one = false;
- for (int k = 0; k < n; ++k) {
- if (!dmap[nodes[i]][k].found) continue;
- Length _len = dmap[nodes[i]][n].dist - dmap[nodes[i]][k].dist;
- int _size = n - k;
- if (!found_one || len * _size < _len * size) {
- found_one = true;
- len = _len;
- size = _size;
- }
- }
- if ( found_one &&
- (!found_min || len * min_size < min_length * size) ) {
- found_min = true;
- min_length = len;
- min_size = size;
- min_node = nodes[i];
- }
+ return true;
+ }
+
+ // Finds all cycles in the policy graph.
+ // Sets _cycle_found to true if a cycle is found and sets
+ // _cycle_length, _cycle_size, _cycle_node to represent the minimum
+ // mean cycle in the policy graph.
+ bool findPolicyCycles(int comp) {
+ typename Graph::template NodeMap<int> level(_graph, -1);
+ _cycle_found = false;
+ Length clength;
+ int csize;
+ int path_cnt = 0;
+ Node u, v;
+ // Searching for cycles
+ for (int i = 0; i < _nodes.size(); ++i) {
+ if (level[_nodes[i]] < 0) {
+ u = _nodes[i];
+ level[u] = path_cnt;
+ while (level[u = _graph.target(_policy[u])] < 0)
+ level[u] = path_cnt;
+ if (level[u] == path_cnt) {
+ // A cycle is found
+ clength = _length[_policy[u]];
+ csize = 1;
+ for (v = u; (v = _graph.target(_policy[v])) != u; ) {
+ clength += _length[_policy[v]];
+ ++csize;
+ }
+ if ( !_cycle_found ||
+ clength * _cycle_size < _cycle_length * csize ) {
+ _cycle_found = true;
+ _cycle_length = clength;
+ _cycle_size = csize;
+ _cycle_node = u;
+ }
+ }
+ ++path_cnt;
+ }
+ }
+ return _cycle_found;
+ }
+
+ // Contracts the policy graph to be connected by cutting all cycles
+ // except for the main cycle (i.e. the minimum mean cycle).
+ void contractPolicyGraph(int comp) {
+ // Finding the component of the main cycle using
+ // reverse BFS search
+ typename Graph::template NodeMap<int> found(_graph, false);
+ std::deque<Node> queue;
+ queue.push_back(_cycle_node);
+ found[_cycle_node] = true;
+ Node u, v;
+ while (!queue.empty()) {
+ v = queue.front(); queue.pop_front();
+ for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
+ u = _graph.source(e);
+ if (_component[u] == comp && !found[u] && _policy[u] == e) {
+ found[u] = true;
+ queue.push_back(u);
+ }
+ }
+ }
+ // Connecting all other nodes to this component using
+ // reverse BFS search
+ queue.clear();
+ for (int i = 0; i < _nodes.size(); ++i)
+ if (found[_nodes[i]]) queue.push_back(_nodes[i]);
+ int found_cnt = queue.size();
+ while (found_cnt < _nodes.size() && !queue.empty()) {
+ v = queue.front(); queue.pop_front();
+ for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
+ u = _graph.source(e);
+ if (_component[u] == comp && !found[u]) {
+ found[u] = true;
+ ++found_cnt;
+ _policy[u] = e;
+ queue.push_back(u);
+ }
+ }
+ }
+ }
+
+ // Computes node distances in the policy graph and updates the
+ // policy graph if the node distances can be improved.
+ bool computeNodeDistances(int comp) {
+ // Computing node distances using reverse BFS search
+ double cycle_mean = (double)_cycle_length / _cycle_size;
+ typename Graph::template NodeMap<int> found(_graph, false);
+ std::deque<Node> queue;
+ queue.push_back(_cycle_node);
+ found[_cycle_node] = true;
+ Node u, v;
+ while (!queue.empty()) {
+ v = queue.front(); queue.pop_front();
+ for (InEdgeIt e(_graph, v); e != INVALID; ++e) {
+ u = _graph.source(e);
+ if (_component[u] == comp && !found[u] && _policy[u] == e) {
+ found[u] = true;
+ _dist[u] = _dist[v] + _length[e] - cycle_mean;
+ queue.push_back(u);
+ }
+ }
+ }
+ // Improving node distances
+ bool improved = false;
+ for (int j = 0; j < _edges.size(); ++j) {
+ Edge e = _edges[j];
+ u = _graph.source(e); v = _graph.target(e);
+ double delta = _dist[v] + _length[e] - cycle_mean;
+ if (_tolerance.less(delta, _dist[u])) {
+ improved = true;
+ _dist[u] = delta;
+ _policy[u] = e;
+ }
}
- return found_min;
+ return improved;
}
public:
@@ -248,19 +266,18 @@
///
/// Runs the algorithm.
///
- /// \return \c true if a cycle exists in the graph.
+ /// \return Returns \c true if a directed cycle exists in the graph.
///
- /// \note Apart from the return value, <tt>m.run()</tt> is just a
+ /// \note Apart from the return value, <tt>mmc.run()</tt> is just a
/// shortcut of the following code.
/// \code
- /// m.init();
- /// m.findMinMean();
- /// m.findCycle();
+ /// mmc.init();
+ /// mmc.findMinMean();
+ /// mmc.findCycle();
/// \endcode
bool run() {
init();
- findMinMean();
- return findCycle();
+ return findMinMean() && findCycle();
}
/// \brief Initializes the internal data structures.
@@ -269,112 +286,85 @@
///
/// \sa reset()
void init() {
- comp_num = stronglyConnectedComponents(graph, comp);
- if (!cycle_path) {
- local_path = true;
- cycle_path = new Path;
+ _tolerance.epsilon(1e-8);
+ if (!_cycle_path) {
+ _local_path = true;
+ _cycle_path = new Path;
}
+ _cycle_found = false;
+ _component_num = stronglyConnectedComponents(_graph, _component);
+ }
+
+ /// \brief Resets the internal data structures.
+ ///
+ /// Resets the internal data structures so that \ref findMinMean()
+ /// and \ref findCycle() can be called again (e.g. when the
+ /// underlaying graph has been modified).
+ ///
+ /// \sa init()
+ void reset() {
+ if (_cycle_path) _cycle_path->clear();
+ _cycle_found = false;
+ _component_num = stronglyConnectedComponents(_graph, _component);
}
- /// \brief Finds the minimum cycle mean value in the graph.
+ /// \brief Finds the minimum cycle mean length in the graph.
///
- /// Computes all the required path data and finds the minimum cycle
- /// mean value in the graph.
+ /// Computes all the required data and finds the minimum cycle mean
+ /// length in the graph.
///
- /// \return \c true if a cycle exists in the graph.
+ /// \return Returns \c true if a directed cycle exists in the graph.
///
/// \pre \ref init() must be called before using this function.
bool findMinMean() {
- cycle_node = INVALID;
- for (comp_cnt = 0; comp_cnt < comp_num; ++comp_cnt) {
- initCurrent();
- processRounds();
-
- Length min_length;
- int min_size;
- Node min_node;
- bool found_min = findCurrentMin(min_length, min_size, min_node);
-
- if ( found_min && (cycle_node == INVALID ||
- min_length * cycle_size < cycle_length * min_size) ) {
- cycle_length = min_length;
- cycle_size = min_size;
- cycle_node = min_node;
- }
+ // Finding the minimum mean cycle in the components
+ for (int comp = 0; comp < _component_num; ++comp) {
+ if (!initCurrentComponent(comp)) continue;
+ while (true) {
+ if (!findPolicyCycles(comp)) return false;
+ contractPolicyGraph(comp);
+ if (!computeNodeDistances(comp)) return true;
+ }
}
- return (cycle_node != INVALID);
}
- /// \brief Finds a critical (minimum mean) cycle.
+ /// \brief Finds a critical (minimum mean) directed cycle.
///
- /// Finds a critical (minimum mean) cycle using the path data
- /// stored in \ref dmap.
+ /// Finds a critical (minimum mean) directed cycle using the data
+ /// computed in the \ref findMinMean() function.
///
- /// \return \c true if a cycle exists in the graph.
+ /// \return Returns \c true if a directed cycle exists in the graph.
///
/// \pre \ref init() and \ref findMinMean() must be called before
/// using this function.
bool findCycle() {
- if (cycle_node == INVALID) return false;
- cycle_length = 0;
- cycle_size = 0;
- IntNodeMap reached(graph, -1);
- int r = reached[cycle_node] = dmap[cycle_node].size() - 1;
- Node u = graph.source(dmap[cycle_node][r].pred);
- while (reached[u] < 0) {
- reached[u] = --r;
- u = graph.source(dmap[u][r].pred);
- }
- r = reached[u];
- Edge e = dmap[u][r].pred;
- cycle_path->addFront(e);
- cycle_length = cycle_length + length[e];
- ++cycle_size;
- Node v;
- while ((v = graph.source(e)) != u) {
- e = dmap[v][--r].pred;
- cycle_path->addFront(e);
- cycle_length = cycle_length + length[e];
- ++cycle_size;
+ if (!_cycle_found) return false;
+ _cycle_path->addBack(_policy[_cycle_node]);
+ for ( Node v = _cycle_node;
+ (v = _graph.target(_policy[v])) != _cycle_node; ) {
+ _cycle_path->addBack(_policy[v]);
}
return true;
}
- /// \brief Resets the internal data structures.
- ///
- /// Resets the internal data structures so that \ref findMinMean()
- /// and \ref findCycle() can be called again (e.g. when the
- /// underlaying graph has been modified).
- ///
- /// \sa init()
- void reset() {
- for (NodeIt u(graph); u != INVALID; ++u)
- dmap[u].clear();
- cycle_node = INVALID;
- if (cycle_path) cycle_path->clear();
- comp_num = stronglyConnectedComponents(graph, comp);
- }
-
/// \brief Returns the total length of the found cycle.
///
/// Returns the total length of the found cycle.
///
- /// \pre \ref run() or \ref findCycle() must be called before
- /// using this function. If only \ref findMinMean() is called,
- /// the return value is not valid.
+ /// \pre \ref run() or \ref findMinMean() must be called before
+ /// using this function.
Length cycleLength() const {
- return cycle_length;
+ return _cycle_length;
}
- /// \brief Returns the number of edges in the found cycle.
+ /// \brief Returns the number of edges on the found cycle.
///
- /// Returns the number of edges in the found cycle.
+ /// Returns the number of edges on the found cycle.
///
- /// \pre \ref run() or \ref findCycle() must be called before
- /// using this function. If only \ref findMinMean() is called,
- /// the return value is not valid.
+ /// \pre \ref run() or \ref findMinMean() must be called before
+ /// using this function.
int cycleEdgeNum() const {
- return cycle_size;
+ return _cycle_size;
}
/// \brief Returns the mean length of the found cycle.
@@ -384,52 +374,53 @@
/// \pre \ref run() or \ref findMinMean() must be called before
/// using this function.
///
- /// \warning LengthMap::Value must be convertible to double.
- ///
- /// \note <tt>m.minMean()</tt> is just a shortcut of the following
- /// code.
+ /// \note <tt>mmc.cycleMean()</tt> is just a shortcut of the
+ /// following code.
/// \code
- /// return m.cycleLength() / double(m.cycleEdgeNum());
+ /// return double(mmc.cycleLength()) / mmc.cycleEdgeNum();
/// \endcode
- /// However if only \ref findMinMean() is called before using this
- /// function, <tt>m.cycleLength()</tt> and <tt>m.cycleEdgeNum()</tt>
- /// are not valid alone, only their ratio is valid.
- double minMean() const {
- return cycle_length / (double)cycle_size;
+ double cycleMean() const {
+ return (double)_cycle_length / _cycle_size;
}
- /// \brief Returns a const reference to the \ref lemon::Path "Path"
- /// structure of the found cycle.
+ /// \brief Returns a const reference to the \ref Path "path"
+ /// structure storing the found cycle.
///
- /// Returns a const reference to the \ref lemon::Path "Path"
- /// structure of the found cycle.
+ /// Returns a const reference to the \ref Path "path"
+ /// structure storing the found cycle.
///
- /// \pre \ref run() must be called before using this function.
+ /// \pre \ref run() or \ref findCycle() must be called before using
+ /// this function.
///
- /// \sa \ref cyclePath()
+ /// \sa cyclePath()
const Path& cycle() const {
- return *cycle_path;
+ return *_cycle_path;
}
- /// \brief Sets the \ref lemon::Path "Path" structure storing the
- /// found cycle.
+ /// \brief Sets the \ref Path "path" structure for storing the found
+ /// cycle.
///
- /// Sets the \ref lemon::Path "Path" structure storing the found
- /// cycle. If you don't use this function before calling
- /// \ref run(), it will allocate one. The destuctor deallocates
- /// this automatically allocated map, of course.
+ /// Sets an external \ref Path "path" structure for storing the
+ /// found cycle.
///
- /// \note The algorithm calls only the \ref lemon::Path::addFront()
- /// "addFront()" method of the given \ref lemon::Path "Path"
+ /// If you don't call this function before calling \ref run() or
+ /// \ref init(), it will allocate a local \ref Path "path"
/// structure.
+ /// The destuctor deallocates this automatically allocated map,
+ /// of course.
+ ///
+ /// \note The algorithm calls only the \ref lemon::Path::addBack()
+ /// "addBack()" function of the given \ref Path "path" structure.
+ ///
+ /// \return <tt>(*this)</tt>
///
- /// \return \c (*this)
+ /// \sa cycle()
MinMeanCycle& cyclePath(Path &path) {
- if (local_path) {
- delete cycle_path;
- local_path = false;
+ if (_local_path) {
+ delete _cycle_path;
+ _local_path = false;
}
- cycle_path = &path;
+ _cycle_path = &path;
return *this;
}
More information about the Lemon-commits
mailing list