1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/lemon/gomory_hu.h Wed Mar 04 14:09:45 2009 +0000
1.3 @@ -0,0 +1,551 @@
1.4 +/* -*- C++ -*-
1.5 + *
1.6 + * This file is a part of LEMON, a generic C++ optimization library
1.7 + *
1.8 + * Copyright (C) 2003-2008
1.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
1.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
1.11 + *
1.12 + * Permission to use, modify and distribute this software is granted
1.13 + * provided that this copyright notice appears in all copies. For
1.14 + * precise terms see the accompanying LICENSE file.
1.15 + *
1.16 + * This software is provided "AS IS" with no warranty of any kind,
1.17 + * express or implied, and with no claim as to its suitability for any
1.18 + * purpose.
1.19 + *
1.20 + */
1.21 +
1.22 +#ifndef LEMON_GOMORY_HU_TREE_H
1.23 +#define LEMON_GOMORY_HU_TREE_H
1.24 +
1.25 +#include <limits>
1.26 +
1.27 +#include <lemon/core.h>
1.28 +#include <lemon/preflow.h>
1.29 +#include <lemon/concept_check.h>
1.30 +#include <lemon/concepts/maps.h>
1.31 +
1.32 +/// \ingroup min_cut
1.33 +/// \file
1.34 +/// \brief Gomory-Hu cut tree in graphs.
1.35 +
1.36 +namespace lemon {
1.37 +
1.38 + /// \ingroup min_cut
1.39 + ///
1.40 + /// \brief Gomory-Hu cut tree algorithm
1.41 + ///
1.42 + /// The Gomory-Hu tree is a tree on the node set of a given graph, but it
1.43 + /// may contain edges which are not in the original graph. It has the
1.44 + /// property that the minimum capacity edge of the path between two nodes
1.45 + /// in this tree has the same weight as the minimum cut in the graph
1.46 + /// between these nodes. Moreover the components obtained by removing
1.47 + /// this edge from the tree determine the corresponding minimum cut.
1.48 + ///
1.49 + /// Therefore once this tree is computed, the minimum cut between any pair
1.50 + /// of nodes can easily be obtained.
1.51 + ///
1.52 + /// The algorithm calculates \e n-1 distinct minimum cuts (currently with
1.53 + /// the \ref Preflow algorithm), therefore the algorithm has
1.54 + /// \f$(O(n^3\sqrt{e})\f$ overall time complexity. It calculates a
1.55 + /// rooted Gomory-Hu tree, its structure and the weights can be obtained
1.56 + /// by \c predNode(), \c predValue() and \c rootDist().
1.57 + ///
1.58 + /// The members \c minCutMap() and \c minCutValue() calculate
1.59 + /// the minimum cut and the minimum cut value between any two nodes
1.60 + /// in the graph. You can also list (iterate on) the nodes and the
1.61 + /// edges of the cuts using \c MinCutNodeIt and \c MinCutEdgeIt.
1.62 + ///
1.63 + /// \tparam GR The type of the undirected graph the algorithm runs on.
1.64 + /// \tparam CAP The type of the edge map describing the edge capacities.
1.65 + /// It is \ref concepts::Graph::EdgeMap "GR::EdgeMap<int>" by default.
1.66 +#ifdef DOXYGEN
1.67 + template <typename GR,
1.68 + typename CAP>
1.69 +#else
1.70 + template <typename GR,
1.71 + typename CAP = typename GR::template EdgeMap<int> >
1.72 +#endif
1.73 + class GomoryHu {
1.74 + public:
1.75 +
1.76 + /// The graph type
1.77 + typedef GR Graph;
1.78 + /// The type of the edge capacity map
1.79 + typedef CAP Capacity;
1.80 + /// The value type of capacities
1.81 + typedef typename Capacity::Value Value;
1.82 +
1.83 + private:
1.84 +
1.85 + TEMPLATE_GRAPH_TYPEDEFS(Graph);
1.86 +
1.87 + const Graph& _graph;
1.88 + const Capacity& _capacity;
1.89 +
1.90 + Node _root;
1.91 + typename Graph::template NodeMap<Node>* _pred;
1.92 + typename Graph::template NodeMap<Value>* _weight;
1.93 + typename Graph::template NodeMap<int>* _order;
1.94 +
1.95 + void createStructures() {
1.96 + if (!_pred) {
1.97 + _pred = new typename Graph::template NodeMap<Node>(_graph);
1.98 + }
1.99 + if (!_weight) {
1.100 + _weight = new typename Graph::template NodeMap<Value>(_graph);
1.101 + }
1.102 + if (!_order) {
1.103 + _order = new typename Graph::template NodeMap<int>(_graph);
1.104 + }
1.105 + }
1.106 +
1.107 + void destroyStructures() {
1.108 + if (_pred) {
1.109 + delete _pred;
1.110 + }
1.111 + if (_weight) {
1.112 + delete _weight;
1.113 + }
1.114 + if (_order) {
1.115 + delete _order;
1.116 + }
1.117 + }
1.118 +
1.119 + public:
1.120 +
1.121 + /// \brief Constructor
1.122 + ///
1.123 + /// Constructor
1.124 + /// \param graph The undirected graph the algorithm runs on.
1.125 + /// \param capacity The edge capacity map.
1.126 + GomoryHu(const Graph& graph, const Capacity& capacity)
1.127 + : _graph(graph), _capacity(capacity),
1.128 + _pred(0), _weight(0), _order(0)
1.129 + {
1.130 + checkConcept<concepts::ReadMap<Edge, Value>, Capacity>();
1.131 + }
1.132 +
1.133 +
1.134 + /// \brief Destructor
1.135 + ///
1.136 + /// Destructor
1.137 + ~GomoryHu() {
1.138 + destroyStructures();
1.139 + }
1.140 +
1.141 + private:
1.142 +
1.143 + // Initialize the internal data structures
1.144 + void init() {
1.145 + createStructures();
1.146 +
1.147 + _root = NodeIt(_graph);
1.148 + for (NodeIt n(_graph); n != INVALID; ++n) {
1.149 + _pred->set(n, _root);
1.150 + _order->set(n, -1);
1.151 + }
1.152 + _pred->set(_root, INVALID);
1.153 + _weight->set(_root, std::numeric_limits<Value>::max());
1.154 + }
1.155 +
1.156 +
1.157 + // Start the algorithm
1.158 + void start() {
1.159 + Preflow<Graph, Capacity> fa(_graph, _capacity, _root, INVALID);
1.160 +
1.161 + for (NodeIt n(_graph); n != INVALID; ++n) {
1.162 + if (n == _root) continue;
1.163 +
1.164 + Node pn = (*_pred)[n];
1.165 + fa.source(n);
1.166 + fa.target(pn);
1.167 +
1.168 + fa.runMinCut();
1.169 +
1.170 + _weight->set(n, fa.flowValue());
1.171 +
1.172 + for (NodeIt nn(_graph); nn != INVALID; ++nn) {
1.173 + if (nn != n && fa.minCut(nn) && (*_pred)[nn] == pn) {
1.174 + _pred->set(nn, n);
1.175 + }
1.176 + }
1.177 + if ((*_pred)[pn] != INVALID && fa.minCut((*_pred)[pn])) {
1.178 + _pred->set(n, (*_pred)[pn]);
1.179 + _pred->set(pn, n);
1.180 + _weight->set(n, (*_weight)[pn]);
1.181 + _weight->set(pn, fa.flowValue());
1.182 + }
1.183 + }
1.184 +
1.185 + _order->set(_root, 0);
1.186 + int index = 1;
1.187 +
1.188 + for (NodeIt n(_graph); n != INVALID; ++n) {
1.189 + std::vector<Node> st;
1.190 + Node nn = n;
1.191 + while ((*_order)[nn] == -1) {
1.192 + st.push_back(nn);
1.193 + nn = (*_pred)[nn];
1.194 + }
1.195 + while (!st.empty()) {
1.196 + _order->set(st.back(), index++);
1.197 + st.pop_back();
1.198 + }
1.199 + }
1.200 + }
1.201 +
1.202 + public:
1.203 +
1.204 + ///\name Execution Control
1.205 +
1.206 + ///@{
1.207 +
1.208 + /// \brief Run the Gomory-Hu algorithm.
1.209 + ///
1.210 + /// This function runs the Gomory-Hu algorithm.
1.211 + void run() {
1.212 + init();
1.213 + start();
1.214 + }
1.215 +
1.216 + /// @}
1.217 +
1.218 + ///\name Query Functions
1.219 + ///The results of the algorithm can be obtained using these
1.220 + ///functions.\n
1.221 + ///\ref run() "run()" should be called before using them.\n
1.222 + ///See also \ref MinCutNodeIt and \ref MinCutEdgeIt.
1.223 +
1.224 + ///@{
1.225 +
1.226 + /// \brief Return the predecessor node in the Gomory-Hu tree.
1.227 + ///
1.228 + /// This function returns the predecessor node in the Gomory-Hu tree.
1.229 + /// If the node is
1.230 + /// the root of the Gomory-Hu tree, then it returns \c INVALID.
1.231 + Node predNode(const Node& node) {
1.232 + return (*_pred)[node];
1.233 + }
1.234 +
1.235 + /// \brief Return the distance from the root node in the Gomory-Hu tree.
1.236 + ///
1.237 + /// This function returns the distance of \c node from the root node
1.238 + /// in the Gomory-Hu tree.
1.239 + int rootDist(const Node& node) {
1.240 + return (*_order)[node];
1.241 + }
1.242 +
1.243 + /// \brief Return the weight of the predecessor edge in the
1.244 + /// Gomory-Hu tree.
1.245 + ///
1.246 + /// This function returns the weight of the predecessor edge in the
1.247 + /// Gomory-Hu tree. If the node is the root, the result is undefined.
1.248 + Value predValue(const Node& node) {
1.249 + return (*_weight)[node];
1.250 + }
1.251 +
1.252 + /// \brief Return the minimum cut value between two nodes
1.253 + ///
1.254 + /// This function returns the minimum cut value between two nodes. The
1.255 + /// algorithm finds the nearest common ancestor in the Gomory-Hu
1.256 + /// tree and calculates the minimum weight edge on the paths to
1.257 + /// the ancestor.
1.258 + Value minCutValue(const Node& s, const Node& t) const {
1.259 + Node sn = s, tn = t;
1.260 + Value value = std::numeric_limits<Value>::max();
1.261 +
1.262 + while (sn != tn) {
1.263 + if ((*_order)[sn] < (*_order)[tn]) {
1.264 + if ((*_weight)[tn] <= value) value = (*_weight)[tn];
1.265 + tn = (*_pred)[tn];
1.266 + } else {
1.267 + if ((*_weight)[sn] <= value) value = (*_weight)[sn];
1.268 + sn = (*_pred)[sn];
1.269 + }
1.270 + }
1.271 + return value;
1.272 + }
1.273 +
1.274 + /// \brief Return the minimum cut between two nodes
1.275 + ///
1.276 + /// This function returns the minimum cut between the nodes \c s and \c t
1.277 + /// in the \c cutMap parameter by setting the nodes in the component of
1.278 + /// \c s to \c true and the other nodes to \c false.
1.279 + ///
1.280 + /// For higher level interfaces, see MinCutNodeIt and MinCutEdgeIt.
1.281 + template <typename CutMap>
1.282 + Value minCutMap(const Node& s, ///< The base node.
1.283 + const Node& t,
1.284 + ///< The node you want to separate from node \c s.
1.285 + CutMap& cutMap
1.286 + ///< The cut will be returned in this map.
1.287 + /// It must be a \c bool (or convertible)
1.288 + /// \ref concepts::ReadWriteMap "ReadWriteMap"
1.289 + /// on the graph nodes.
1.290 + ) const {
1.291 + Node sn = s, tn = t;
1.292 + bool s_root=false;
1.293 + Node rn = INVALID;
1.294 + Value value = std::numeric_limits<Value>::max();
1.295 +
1.296 + while (sn != tn) {
1.297 + if ((*_order)[sn] < (*_order)[tn]) {
1.298 + if ((*_weight)[tn] <= value) {
1.299 + rn = tn;
1.300 + s_root = false;
1.301 + value = (*_weight)[tn];
1.302 + }
1.303 + tn = (*_pred)[tn];
1.304 + } else {
1.305 + if ((*_weight)[sn] <= value) {
1.306 + rn = sn;
1.307 + s_root = true;
1.308 + value = (*_weight)[sn];
1.309 + }
1.310 + sn = (*_pred)[sn];
1.311 + }
1.312 + }
1.313 +
1.314 + typename Graph::template NodeMap<bool> reached(_graph, false);
1.315 + reached.set(_root, true);
1.316 + cutMap.set(_root, !s_root);
1.317 + reached.set(rn, true);
1.318 + cutMap.set(rn, s_root);
1.319 +
1.320 + std::vector<Node> st;
1.321 + for (NodeIt n(_graph); n != INVALID; ++n) {
1.322 + st.clear();
1.323 + Node nn = n;
1.324 + while (!reached[nn]) {
1.325 + st.push_back(nn);
1.326 + nn = (*_pred)[nn];
1.327 + }
1.328 + while (!st.empty()) {
1.329 + cutMap.set(st.back(), cutMap[nn]);
1.330 + st.pop_back();
1.331 + }
1.332 + }
1.333 +
1.334 + return value;
1.335 + }
1.336 +
1.337 + ///@}
1.338 +
1.339 + friend class MinCutNodeIt;
1.340 +
1.341 + /// Iterate on the nodes of a minimum cut
1.342 +
1.343 + /// This iterator class lists the nodes of a minimum cut found by
1.344 + /// GomoryHu. Before using it, you must allocate a GomoryHu class,
1.345 + /// and call its \ref GomoryHu::run() "run()" method.
1.346 + ///
1.347 + /// This example counts the nodes in the minimum cut separating \c s from
1.348 + /// \c t.
1.349 + /// \code
1.350 + /// GomoruHu<Graph> gom(g, capacities);
1.351 + /// gom.run();
1.352 + /// int cnt=0;
1.353 + /// for(GomoruHu<Graph>::MinCutNodeIt n(gom,s,t); n!=INVALID; ++n) ++cnt;
1.354 + /// \endcode
1.355 + class MinCutNodeIt
1.356 + {
1.357 + bool _side;
1.358 + typename Graph::NodeIt _node_it;
1.359 + typename Graph::template NodeMap<bool> _cut;
1.360 + public:
1.361 + /// Constructor
1.362 +
1.363 + /// Constructor.
1.364 + ///
1.365 + MinCutNodeIt(GomoryHu const &gomory,
1.366 + ///< The GomoryHu class. You must call its
1.367 + /// run() method
1.368 + /// before initializing this iterator.
1.369 + const Node& s, ///< The base node.
1.370 + const Node& t,
1.371 + ///< The node you want to separate from node \c s.
1.372 + bool side=true
1.373 + ///< If it is \c true (default) then the iterator lists
1.374 + /// the nodes of the component containing \c s,
1.375 + /// otherwise it lists the other component.
1.376 + /// \note As the minimum cut is not always unique,
1.377 + /// \code
1.378 + /// MinCutNodeIt(gomory, s, t, true);
1.379 + /// \endcode
1.380 + /// and
1.381 + /// \code
1.382 + /// MinCutNodeIt(gomory, t, s, false);
1.383 + /// \endcode
1.384 + /// does not necessarily give the same set of nodes.
1.385 + /// However it is ensured that
1.386 + /// \code
1.387 + /// MinCutNodeIt(gomory, s, t, true);
1.388 + /// \endcode
1.389 + /// and
1.390 + /// \code
1.391 + /// MinCutNodeIt(gomory, s, t, false);
1.392 + /// \endcode
1.393 + /// together list each node exactly once.
1.394 + )
1.395 + : _side(side), _cut(gomory._graph)
1.396 + {
1.397 + gomory.minCutMap(s,t,_cut);
1.398 + for(_node_it=typename Graph::NodeIt(gomory._graph);
1.399 + _node_it!=INVALID && _cut[_node_it]!=_side;
1.400 + ++_node_it) {}
1.401 + }
1.402 + /// Conversion to \c Node
1.403 +
1.404 + /// Conversion to \c Node.
1.405 + ///
1.406 + operator typename Graph::Node() const
1.407 + {
1.408 + return _node_it;
1.409 + }
1.410 + bool operator==(Invalid) { return _node_it==INVALID; }
1.411 + bool operator!=(Invalid) { return _node_it!=INVALID; }
1.412 + /// Next node
1.413 +
1.414 + /// Next node.
1.415 + ///
1.416 + MinCutNodeIt &operator++()
1.417 + {
1.418 + for(++_node_it;_node_it!=INVALID&&_cut[_node_it]!=_side;++_node_it) {}
1.419 + return *this;
1.420 + }
1.421 + /// Postfix incrementation
1.422 +
1.423 + /// Postfix incrementation.
1.424 + ///
1.425 + /// \warning This incrementation
1.426 + /// returns a \c Node, not a \c MinCutNodeIt, as one may
1.427 + /// expect.
1.428 + typename Graph::Node operator++(int)
1.429 + {
1.430 + typename Graph::Node n=*this;
1.431 + ++(*this);
1.432 + return n;
1.433 + }
1.434 + };
1.435 +
1.436 + friend class MinCutEdgeIt;
1.437 +
1.438 + /// Iterate on the edges of a minimum cut
1.439 +
1.440 + /// This iterator class lists the edges of a minimum cut found by
1.441 + /// GomoryHu. Before using it, you must allocate a GomoryHu class,
1.442 + /// and call its \ref GomoryHu::run() "run()" method.
1.443 + ///
1.444 + /// This example computes the value of the minimum cut separating \c s from
1.445 + /// \c t.
1.446 + /// \code
1.447 + /// GomoruHu<Graph> gom(g, capacities);
1.448 + /// gom.run();
1.449 + /// int value=0;
1.450 + /// for(GomoruHu<Graph>::MinCutEdgeIt e(gom,s,t); e!=INVALID; ++e)
1.451 + /// value+=capacities[e];
1.452 + /// \endcode
1.453 + /// the result will be the same as it is returned by
1.454 + /// \ref GomoryHu::minCutValue() "gom.minCutValue(s,t)"
1.455 + class MinCutEdgeIt
1.456 + {
1.457 + bool _side;
1.458 + const Graph &_graph;
1.459 + typename Graph::NodeIt _node_it;
1.460 + typename Graph::OutArcIt _arc_it;
1.461 + typename Graph::template NodeMap<bool> _cut;
1.462 + void step()
1.463 + {
1.464 + ++_arc_it;
1.465 + while(_node_it!=INVALID && _arc_it==INVALID)
1.466 + {
1.467 + for(++_node_it;_node_it!=INVALID&&!_cut[_node_it];++_node_it) {}
1.468 + if(_node_it!=INVALID)
1.469 + _arc_it=typename Graph::OutArcIt(_graph,_node_it);
1.470 + }
1.471 + }
1.472 +
1.473 + public:
1.474 + MinCutEdgeIt(GomoryHu const &gomory,
1.475 + ///< The GomoryHu class. You must call its
1.476 + /// run() method
1.477 + /// before initializing this iterator.
1.478 + const Node& s, ///< The base node.
1.479 + const Node& t,
1.480 + ///< The node you want to separate from node \c s.
1.481 + bool side=true
1.482 + ///< If it is \c true (default) then the listed arcs
1.483 + /// will be oriented from the
1.484 + /// the nodes of the component containing \c s,
1.485 + /// otherwise they will be oriented in the opposite
1.486 + /// direction.
1.487 + )
1.488 + : _graph(gomory._graph), _cut(_graph)
1.489 + {
1.490 + gomory.minCutMap(s,t,_cut);
1.491 + if(!side)
1.492 + for(typename Graph::NodeIt n(_graph);n!=INVALID;++n)
1.493 + _cut[n]=!_cut[n];
1.494 +
1.495 + for(_node_it=typename Graph::NodeIt(_graph);
1.496 + _node_it!=INVALID && !_cut[_node_it];
1.497 + ++_node_it) {}
1.498 + _arc_it = _node_it!=INVALID ?
1.499 + typename Graph::OutArcIt(_graph,_node_it) : INVALID;
1.500 + while(_node_it!=INVALID && _arc_it == INVALID)
1.501 + {
1.502 + for(++_node_it; _node_it!=INVALID&&!_cut[_node_it]; ++_node_it) {}
1.503 + if(_node_it!=INVALID)
1.504 + _arc_it= typename Graph::OutArcIt(_graph,_node_it);
1.505 + }
1.506 + while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
1.507 + }
1.508 + /// Conversion to \c Arc
1.509 +
1.510 + /// Conversion to \c Arc.
1.511 + ///
1.512 + operator typename Graph::Arc() const
1.513 + {
1.514 + return _arc_it;
1.515 + }
1.516 + /// Conversion to \c Edge
1.517 +
1.518 + /// Conversion to \c Edge.
1.519 + ///
1.520 + operator typename Graph::Edge() const
1.521 + {
1.522 + return _arc_it;
1.523 + }
1.524 + bool operator==(Invalid) { return _node_it==INVALID; }
1.525 + bool operator!=(Invalid) { return _node_it!=INVALID; }
1.526 + /// Next edge
1.527 +
1.528 + /// Next edge.
1.529 + ///
1.530 + MinCutEdgeIt &operator++()
1.531 + {
1.532 + step();
1.533 + while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
1.534 + return *this;
1.535 + }
1.536 + /// Postfix incrementation
1.537 +
1.538 + /// Postfix incrementation.
1.539 + ///
1.540 + /// \warning This incrementation
1.541 + /// returns an \c Arc, not a \c MinCutEdgeIt, as one may expect.
1.542 + typename Graph::Arc operator++(int)
1.543 + {
1.544 + typename Graph::Arc e=*this;
1.545 + ++(*this);
1.546 + return e;
1.547 + }
1.548 + };
1.549 +
1.550 + };
1.551 +
1.552 +}
1.553 +
1.554 +#endif