COIN-OR::LEMON - Graph Library

source: lemon-1.2/lemon/gomory_hu_tree.h @ 544:ccd2d3a3001e

Last change on this file since 544:ccd2d3a3001e was 544:ccd2d3a3001e, checked in by Alpar Juttner <alpar@…>, 11 years ago

Cut iterators for GomoryHuTree? + doc cleanup + bug fixes (#66)

File size: 16.0 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_GOMORY_HU_TREE_H
20#define LEMON_GOMORY_HU_TREE_H
21
22#include <limits>
23
24#include <lemon/core.h>
25#include <lemon/preflow.h>
26#include <lemon/concept_check.h>
27#include <lemon/concepts/maps.h>
28
29/// \ingroup min_cut
30/// \file
31/// \brief Gomory-Hu cut tree in graphs.
32
33namespace lemon {
34
35  /// \ingroup min_cut
36  ///
37  /// \brief Gomory-Hu cut tree algorithm
38  ///
39  /// The Gomory-Hu tree is a tree on the node set of the graph, but it
40  /// may contain edges which are not in the original digraph. It has the
41  /// property that the minimum capacity edge of the path between two nodes
42  /// in this tree has the same weight as the minimum cut in the digraph
43  /// between these nodes. Moreover the components obtained by removing
44  /// this edge from the tree determine the corresponding minimum cut.
45  ///
46  /// Therefore once this tree is computed, the minimum cut between any pair
47  /// of nodes can easily be obtained.
48  ///
49  /// The algorithm calculates \e n-1 distinct minimum cuts (currently with
50  /// the \ref Preflow algorithm), therefore the algorithm has
51  /// \f$(O(n^3\sqrt{e})\f$ overall time complexity. It calculates a
52  /// rooted Gomory-Hu tree, its structure and the weights can be obtained
53  /// by \c predNode(), \c predValue() and \c rootDist().
54  ///
55  /// The members \c minCutMap() and \c minCutValue() calculate
56  /// the minimum cut and the minimum cut value between any two node
57  /// in the digraph. You can also list (iterate on) the nodes and the
58  /// edges of the cuts using MinCutNodeIt and MinCutEdgeIt.
59  ///
60  /// \tparam GR The undirected graph data structure the algorithm will run on
61  /// \tparam CAP type of the EdgeMap describing the Edge capacities.
62  /// it is typename GR::template EdgeMap<int> by default.
63  template <typename GR,
64            typename CAP = typename GR::template EdgeMap<int>
65            >
66  class GomoryHuTree {
67  public:
68
69    /// The graph type
70    typedef GR Graph;
71    /// The type if the edge capacity map
72    typedef CAP Capacity;
73    /// The value type of capacities
74    typedef typename Capacity::Value Value;
75   
76  private:
77
78    TEMPLATE_GRAPH_TYPEDEFS(Graph);
79
80    const Graph& _graph;
81    const Capacity& _capacity;
82
83    Node _root;
84    typename Graph::template NodeMap<Node>* _pred;
85    typename Graph::template NodeMap<Value>* _weight;
86    typename Graph::template NodeMap<int>* _order;
87
88    void createStructures() {
89      if (!_pred) {
90        _pred = new typename Graph::template NodeMap<Node>(_graph);
91      }
92      if (!_weight) {
93        _weight = new typename Graph::template NodeMap<Value>(_graph);
94      }
95      if (!_order) {
96        _order = new typename Graph::template NodeMap<int>(_graph);
97      }
98    }
99
100    void destroyStructures() {
101      if (_pred) {
102        delete _pred;
103      }
104      if (_weight) {
105        delete _weight;
106      }
107      if (_order) {
108        delete _order;
109      }
110    }
111 
112  public:
113
114    /// \brief Constructor
115    ///
116    /// Constructor
117    /// \param graph The graph the algorithm will run on.
118    /// \param capacity The capacity map.
119    GomoryHuTree(const Graph& graph, const Capacity& capacity)
120      : _graph(graph), _capacity(capacity),
121        _pred(0), _weight(0), _order(0)
122    {
123      checkConcept<concepts::ReadMap<Edge, Value>, Capacity>();
124    }
125
126
127    /// \brief Destructor
128    ///
129    /// Destructor
130    ~GomoryHuTree() {
131      destroyStructures();
132    }
133
134    // \brief Initialize the internal data structures.
135    //
136    // This function initializes the internal data structures.
137    //
138    void init() {
139      createStructures();
140
141      _root = NodeIt(_graph);
142      for (NodeIt n(_graph); n != INVALID; ++n) {
143        _pred->set(n, _root);
144        _order->set(n, -1);
145      }
146      _pred->set(_root, INVALID);
147      _weight->set(_root, std::numeric_limits<Value>::max());
148    }
149
150
151    // \brief Start the algorithm
152    //
153    // This function starts the algorithm.
154    //
155    // \pre \ref init() must be called before using this function.
156    //
157    void start() {
158      Preflow<Graph, Capacity> fa(_graph, _capacity, _root, INVALID);
159
160      for (NodeIt n(_graph); n != INVALID; ++n) {
161        if (n == _root) continue;
162
163        Node pn = (*_pred)[n];
164        fa.source(n);
165        fa.target(pn);
166
167        fa.runMinCut();
168
169        _weight->set(n, fa.flowValue());
170
171        for (NodeIt nn(_graph); nn != INVALID; ++nn) {
172          if (nn != n && fa.minCut(nn) && (*_pred)[nn] == pn) {
173            _pred->set(nn, n);
174          }
175        }
176        if ((*_pred)[pn] != INVALID && fa.minCut((*_pred)[pn])) {
177          _pred->set(n, (*_pred)[pn]);
178          _pred->set(pn, n);
179          _weight->set(n, (*_weight)[pn]);
180          _weight->set(pn, fa.flowValue());     
181        }
182      }
183
184      _order->set(_root, 0);
185      int index = 1;
186
187      for (NodeIt n(_graph); n != INVALID; ++n) {
188        std::vector<Node> st;
189        Node nn = n;
190        while ((*_order)[nn] == -1) {
191          st.push_back(nn);
192          nn = (*_pred)[nn];
193        }
194        while (!st.empty()) {
195          _order->set(st.back(), index++);
196          st.pop_back();
197        }
198      }
199    }
200
201    ///\name Execution Control
202 
203    ///@{
204
205    /// \brief Run the Gomory-Hu algorithm.
206    ///
207    /// This function runs the Gomory-Hu algorithm.
208    void run() {
209      init();
210      start();
211    }
212   
213    /// @}
214
215    ///\name Query Functions
216    ///The results of the algorithm can be obtained using these
217    ///functions.\n
218    ///The \ref run() "run()" should be called before using them.\n
219    ///See also MinCutNodeIt and MinCutEdgeIt
220
221    ///@{
222
223    /// \brief Return the predecessor node in the Gomory-Hu tree.
224    ///
225    /// This function returns the predecessor node in the Gomory-Hu tree.
226    /// If the node is
227    /// the root of the Gomory-Hu tree, then it returns \c INVALID.
228    Node predNode(const Node& node) {
229      return (*_pred)[node];
230    }
231
232    /// \brief Return the distance from the root node in the Gomory-Hu tree.
233    ///
234    /// This function returns the distance of \c node from the root node
235    /// in the Gomory-Hu tree.
236    int rootDist(const Node& node) {
237      return (*_order)[node];
238    }
239
240    /// \brief Return the weight of the predecessor edge in the
241    /// Gomory-Hu tree.
242    ///
243    /// This function returns the weight of the predecessor edge in the
244    /// Gomory-Hu tree.  If the node is the root, the result is undefined.
245    Value predValue(const Node& node) {
246      return (*_weight)[node];
247    }
248
249    /// \brief Return the minimum cut value between two nodes
250    ///
251    /// This function returns the minimum cut value between two nodes. The
252    /// algorithm finds the nearest common ancestor in the Gomory-Hu
253    /// tree and calculates the minimum weight arc on the paths to
254    /// the ancestor.
255    Value minCutValue(const Node& s, const Node& t) const {
256      Node sn = s, tn = t;
257      Value value = std::numeric_limits<Value>::max();
258     
259      while (sn != tn) {
260        if ((*_order)[sn] < (*_order)[tn]) {
261          if ((*_weight)[tn] <= value) value = (*_weight)[tn];
262          tn = (*_pred)[tn];
263        } else {
264          if ((*_weight)[sn] <= value) value = (*_weight)[sn];
265          sn = (*_pred)[sn];
266        }
267      }
268      return value;
269    }
270
271    /// \brief Return the minimum cut between two nodes
272    ///
273    /// This function returns the minimum cut between the nodes \c s and \c t
274    /// the \r cutMap parameter by setting the nodes in the component of
275    /// \c \s to true and the other nodes to false.
276    ///
277    /// The \c cutMap should be \ref concepts::ReadWriteMap
278    /// "ReadWriteMap".
279    ///
280    /// For higher level interfaces, see MinCutNodeIt and MinCutEdgeIt
281    template <typename CutMap>
282    Value minCutMap(const Node& s, ///< Base node
283                    const Node& t,
284                    ///< The node you want to separate from Node s.
285                    CutMap& cutMap
286                    ///< The cut will be return in this map.
287                    /// It must be a \c bool \ref concepts::ReadWriteMap
288                    /// "ReadWriteMap" on the graph nodes.
289                    ) const {
290      Node sn = s, tn = t;
291      bool s_root=false;
292      Node rn = INVALID;
293      Value value = std::numeric_limits<Value>::max();
294     
295      while (sn != tn) {
296        if ((*_order)[sn] < (*_order)[tn]) {
297          if ((*_weight)[tn] <= value) {
298            rn = tn;
299            s_root = false;
300            value = (*_weight)[tn];
301          }
302          tn = (*_pred)[tn];
303        } else {
304          if ((*_weight)[sn] <= value) {
305            rn = sn;
306            s_root = true;
307            value = (*_weight)[sn];
308          }
309          sn = (*_pred)[sn];
310        }
311      }
312
313      typename Graph::template NodeMap<bool> reached(_graph, false);
314      reached.set(_root, true);
315      cutMap.set(_root, !s_root);
316      reached.set(rn, true);
317      cutMap.set(rn, s_root);
318
319      std::vector<Node> st;
320      for (NodeIt n(_graph); n != INVALID; ++n) {
321        st.clear();
322        Node nn = n;
323        while (!reached[nn]) {
324          st.push_back(nn);
325          nn = (*_pred)[nn];
326        }
327        while (!st.empty()) {
328          cutMap.set(st.back(), cutMap[nn]);
329          st.pop_back();
330        }
331      }
332     
333      return value;
334    }
335
336    ///@}
337
338    friend class MinCutNodeIt;
339
340    /// Iterate on the nodes of a minimum cut
341   
342    /// This iterator class lists the nodes of a minimum cut found by
343    /// GomoryHuTree. Before using it, you must allocate a GomoryHuTree class,
344    /// and call its \ref GomoryHuTree::run() "run()" method.
345    ///
346    /// This example counts the nodes in the minimum cut separating \c s from
347    /// \c t.
348    /// \code
349    /// GomoruHuTree<Graph> gom(g, capacities);
350    /// gom.run();
351    /// int sum=0;
352    /// for(GomoruHuTree<Graph>::MinCutNodeIt n(gom,s,t);n!=INVALID;++n) ++sum;
353    /// \endcode
354    class MinCutNodeIt
355    {
356      bool _side;
357      typename Graph::NodeIt _node_it;
358      typename Graph::template NodeMap<bool> _cut;
359    public:
360      /// Constructor
361
362      /// Constructor
363      ///
364      MinCutNodeIt(GomoryHuTree const &gomory,
365                   ///< The GomoryHuTree class. You must call its
366                   ///  run() method
367                   ///  before initializing this iterator
368                   const Node& s, ///< Base node
369                   const Node& t,
370                   ///< The node you want to separate from Node s.
371                   bool side=true
372                   ///< If it is \c true (default) then the iterator lists
373                   ///  the nodes of the component containing \c s,
374                   ///  otherwise it lists the other component.
375                   /// \note As the minimum cut is not always unique,
376                   /// \code
377                   /// MinCutNodeIt(gomory, s, t, true);
378                   /// \endcode
379                   /// and
380                   /// \code
381                   /// MinCutNodeIt(gomory, t, s, false);
382                   /// \endcode
383                   /// does not necessarily give the same set of nodes.
384                   /// However it is ensured that
385                   /// \code
386                   /// MinCutNodeIt(gomory, s, t, true);
387                   /// \endcode
388                   /// and
389                   /// \code
390                   /// MinCutNodeIt(gomory, s, t, false);
391                   /// \endcode
392                   /// together list each node exactly once.
393                   )
394        : _side(side), _cut(gomory._graph)
395      {
396        gomory.minCutMap(s,t,_cut);
397        for(_node_it=typename Graph::NodeIt(gomory._graph);
398            _node_it!=INVALID && _cut[_node_it]!=_side;
399            ++_node_it) {}
400      }
401      /// Conversion to Node
402
403      /// Conversion to Node
404      ///
405      operator typename Graph::Node() const
406      {
407        return _node_it;
408      }
409      bool operator==(Invalid) { return _node_it==INVALID; }
410      bool operator!=(Invalid) { return _node_it!=INVALID; }
411      /// Next node
412
413      /// Next node
414      ///
415      MinCutNodeIt &operator++()
416      {
417        for(++_node_it;_node_it!=INVALID&&_cut[_node_it]!=_side;++_node_it) {}
418        return *this;
419      }
420      /// Postfix incrementation
421
422      /// Postfix incrementation
423      ///
424      /// \warning This incrementation
425      /// returns a \c Node, not a \ref MinCutNodeIt, as one may
426      /// expect.
427      typename Graph::Node operator++(int)
428      {
429        typename Graph::Node n=*this;
430        ++(*this);
431        return n;
432      }
433    };
434   
435    friend class MinCutEdgeIt;
436   
437    /// Iterate on the edges of a minimum cut
438   
439    /// This iterator class lists the edges of a minimum cut found by
440    /// GomoryHuTree. Before using it, you must allocate a GomoryHuTree class,
441    /// and call its \ref GomoryHuTree::run() "run()" method.
442    ///
443    /// This example computes the value of the minimum cut separating \c s from
444    /// \c t.
445    /// \code
446    /// GomoruHuTree<Graph> gom(g, capacities);
447    /// gom.run();
448    /// int value=0;
449    /// for(GomoruHuTree<Graph>::MinCutEdgeIt e(gom,s,t);e!=INVALID;++e)
450    ///   value+=capacities[e];
451    /// \endcode
452    /// the result will be the same as it is returned by
453    /// \ref GomoryHuTree::minCostValue() "gom.minCostValue(s,t)"
454    class MinCutEdgeIt
455    {
456      bool _side;
457      const Graph &_graph;
458      typename Graph::NodeIt _node_it;
459      typename Graph::OutArcIt _arc_it;
460      typename Graph::template NodeMap<bool> _cut;
461      void step()
462      {
463        ++_arc_it;
464        while(_node_it!=INVALID && _arc_it==INVALID)
465          {
466            for(++_node_it;_node_it!=INVALID&&!_cut[_node_it];++_node_it) {}
467            if(_node_it!=INVALID)
468              _arc_it=typename Graph::OutArcIt(_graph,_node_it);
469          }
470      }
471     
472    public:
473      MinCutEdgeIt(GomoryHuTree const &gomory,
474                   ///< The GomoryHuTree class. You must call its
475                   ///  run() method
476                   ///  before initializing this iterator
477                   const Node& s,  ///< Base node
478                   const Node& t,
479                   ///< The node you want to separate from Node s.
480                   bool side=true
481                   ///< If it is \c true (default) then the listed arcs
482                   ///  will be oriented from the
483                   ///  the nodes of the component containing \c s,
484                   ///  otherwise they will be oriented in the opposite
485                   ///  direction.
486                   )
487        : _graph(gomory._graph), _cut(_graph)
488      {
489        gomory.minCutMap(s,t,_cut);
490        if(!side)
491          for(typename Graph::NodeIt n(_graph);n!=INVALID;++n)
492            _cut[n]=!_cut[n];
493
494        for(_node_it=typename Graph::NodeIt(_graph);
495            _node_it!=INVALID && !_cut[_node_it];
496            ++_node_it) {}
497        _arc_it = _node_it!=INVALID ?
498          typename Graph::OutArcIt(_graph,_node_it) : INVALID;
499        while(_node_it!=INVALID && _arc_it == INVALID)
500          {
501            for(++_node_it; _node_it!=INVALID&&!_cut[_node_it]; ++_node_it) {}
502            if(_node_it!=INVALID)
503              _arc_it= typename Graph::OutArcIt(_graph,_node_it);
504          }
505        while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
506      }
507      /// Conversion to Arc
508
509      /// Conversion to Arc
510      ///
511      operator typename Graph::Arc() const
512      {
513        return _arc_it;
514      }
515      /// Conversion to Edge
516
517      /// Conversion to Edge
518      ///
519      operator typename Graph::Edge() const
520      {
521        return _arc_it;
522      }
523      bool operator==(Invalid) { return _node_it==INVALID; }
524      bool operator!=(Invalid) { return _node_it!=INVALID; }
525      /// Next edge
526
527      /// Next edge
528      ///
529      MinCutEdgeIt &operator++()
530      {
531        step();
532        while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
533        return *this;
534      }
535      /// Postfix incrementation
536     
537      /// Postfix incrementation
538      ///
539      /// \warning This incrementation
540      /// returns a \c Arc, not a \ref MinCutEdgeIt, as one may
541      /// expect.
542      typename Graph::Arc operator++(int)
543      {
544        typename Graph::Arc e=*this;
545        ++(*this);
546        return e;
547      }
548    };
549
550  };
551
552}
553
554#endif
Note: See TracBrowser for help on using the repository browser.