COIN-OR::LEMON - Graph Library

source: lemon-1.2/lemon/gomory_hu.h @ 571:d5c39e9d1a4e

Last change on this file since 571:d5c39e9d1a4e was 546:d6b40ebb2617, checked in by Peter Kovacs <kpeter@…>, 15 years ago

Doc improvements in GomoryHu? (#66)
And make init() and start() private + bug fix in the test file.

File size: 15.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_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 a given graph, but it
40  /// may contain edges which are not in the original graph. 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 graph
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 nodes
57  /// in the graph. You can also list (iterate on) the nodes and the
58  /// edges of the cuts using \c MinCutNodeIt and \c MinCutEdgeIt.
59  ///
60  /// \tparam GR The type of the undirected graph the algorithm runs on.
61  /// \tparam CAP The type of the edge map describing the edge capacities.
62  /// It is \ref concepts::Graph::EdgeMap "GR::EdgeMap<int>" by default.
63#ifdef DOXYGEN
64  template <typename GR,
65            typename CAP>
66#else
67  template <typename GR,
68            typename CAP = typename GR::template EdgeMap<int> >
69#endif
70  class GomoryHu {
71  public:
72
73    /// The graph type
74    typedef GR Graph;
75    /// The type of the edge capacity map
76    typedef CAP Capacity;
77    /// The value type of capacities
78    typedef typename Capacity::Value Value;
79   
80  private:
81
82    TEMPLATE_GRAPH_TYPEDEFS(Graph);
83
84    const Graph& _graph;
85    const Capacity& _capacity;
86
87    Node _root;
88    typename Graph::template NodeMap<Node>* _pred;
89    typename Graph::template NodeMap<Value>* _weight;
90    typename Graph::template NodeMap<int>* _order;
91
92    void createStructures() {
93      if (!_pred) {
94        _pred = new typename Graph::template NodeMap<Node>(_graph);
95      }
96      if (!_weight) {
97        _weight = new typename Graph::template NodeMap<Value>(_graph);
98      }
99      if (!_order) {
100        _order = new typename Graph::template NodeMap<int>(_graph);
101      }
102    }
103
104    void destroyStructures() {
105      if (_pred) {
106        delete _pred;
107      }
108      if (_weight) {
109        delete _weight;
110      }
111      if (_order) {
112        delete _order;
113      }
114    }
115 
116  public:
117
118    /// \brief Constructor
119    ///
120    /// Constructor
121    /// \param graph The undirected graph the algorithm runs on.
122    /// \param capacity The edge capacity map.
123    GomoryHu(const Graph& graph, const Capacity& capacity)
124      : _graph(graph), _capacity(capacity),
125        _pred(0), _weight(0), _order(0)
126    {
127      checkConcept<concepts::ReadMap<Edge, Value>, Capacity>();
128    }
129
130
131    /// \brief Destructor
132    ///
133    /// Destructor
134    ~GomoryHu() {
135      destroyStructures();
136    }
137
138  private:
139 
140    // Initialize the internal data structures
141    void init() {
142      createStructures();
143
144      _root = NodeIt(_graph);
145      for (NodeIt n(_graph); n != INVALID; ++n) {
146        _pred->set(n, _root);
147        _order->set(n, -1);
148      }
149      _pred->set(_root, INVALID);
150      _weight->set(_root, std::numeric_limits<Value>::max());
151    }
152
153
154    // Start the algorithm
155    void start() {
156      Preflow<Graph, Capacity> fa(_graph, _capacity, _root, INVALID);
157
158      for (NodeIt n(_graph); n != INVALID; ++n) {
159        if (n == _root) continue;
160
161        Node pn = (*_pred)[n];
162        fa.source(n);
163        fa.target(pn);
164
165        fa.runMinCut();
166
167        _weight->set(n, fa.flowValue());
168
169        for (NodeIt nn(_graph); nn != INVALID; ++nn) {
170          if (nn != n && fa.minCut(nn) && (*_pred)[nn] == pn) {
171            _pred->set(nn, n);
172          }
173        }
174        if ((*_pred)[pn] != INVALID && fa.minCut((*_pred)[pn])) {
175          _pred->set(n, (*_pred)[pn]);
176          _pred->set(pn, n);
177          _weight->set(n, (*_weight)[pn]);
178          _weight->set(pn, fa.flowValue());     
179        }
180      }
181
182      _order->set(_root, 0);
183      int index = 1;
184
185      for (NodeIt n(_graph); n != INVALID; ++n) {
186        std::vector<Node> st;
187        Node nn = n;
188        while ((*_order)[nn] == -1) {
189          st.push_back(nn);
190          nn = (*_pred)[nn];
191        }
192        while (!st.empty()) {
193          _order->set(st.back(), index++);
194          st.pop_back();
195        }
196      }
197    }
198
199  public:
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    ///\ref run() "run()" should be called before using them.\n
219    ///See also \ref MinCutNodeIt and \ref 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 edge 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    /// in the \c cutMap parameter by setting the nodes in the component of
275    /// \c s to \c true and the other nodes to \c false.
276    ///
277    /// For higher level interfaces, see MinCutNodeIt and MinCutEdgeIt.
278    template <typename CutMap>
279    Value minCutMap(const Node& s, ///< The base node.
280                    const Node& t,
281                    ///< The node you want to separate from node \c s.
282                    CutMap& cutMap
283                    ///< The cut will be returned in this map.
284                    /// It must be a \c bool (or convertible)
285                    /// \ref concepts::ReadWriteMap "ReadWriteMap"
286                    /// on the graph nodes.
287                    ) const {
288      Node sn = s, tn = t;
289      bool s_root=false;
290      Node rn = INVALID;
291      Value value = std::numeric_limits<Value>::max();
292     
293      while (sn != tn) {
294        if ((*_order)[sn] < (*_order)[tn]) {
295          if ((*_weight)[tn] <= value) {
296            rn = tn;
297            s_root = false;
298            value = (*_weight)[tn];
299          }
300          tn = (*_pred)[tn];
301        } else {
302          if ((*_weight)[sn] <= value) {
303            rn = sn;
304            s_root = true;
305            value = (*_weight)[sn];
306          }
307          sn = (*_pred)[sn];
308        }
309      }
310
311      typename Graph::template NodeMap<bool> reached(_graph, false);
312      reached.set(_root, true);
313      cutMap.set(_root, !s_root);
314      reached.set(rn, true);
315      cutMap.set(rn, s_root);
316
317      std::vector<Node> st;
318      for (NodeIt n(_graph); n != INVALID; ++n) {
319        st.clear();
320        Node nn = n;
321        while (!reached[nn]) {
322          st.push_back(nn);
323          nn = (*_pred)[nn];
324        }
325        while (!st.empty()) {
326          cutMap.set(st.back(), cutMap[nn]);
327          st.pop_back();
328        }
329      }
330     
331      return value;
332    }
333
334    ///@}
335
336    friend class MinCutNodeIt;
337
338    /// Iterate on the nodes of a minimum cut
339   
340    /// This iterator class lists the nodes of a minimum cut found by
341    /// GomoryHu. Before using it, you must allocate a GomoryHu class,
342    /// and call its \ref GomoryHu::run() "run()" method.
343    ///
344    /// This example counts the nodes in the minimum cut separating \c s from
345    /// \c t.
346    /// \code
347    /// GomoruHu<Graph> gom(g, capacities);
348    /// gom.run();
349    /// int cnt=0;
350    /// for(GomoruHu<Graph>::MinCutNodeIt n(gom,s,t); n!=INVALID; ++n) ++cnt;
351    /// \endcode
352    class MinCutNodeIt
353    {
354      bool _side;
355      typename Graph::NodeIt _node_it;
356      typename Graph::template NodeMap<bool> _cut;
357    public:
358      /// Constructor
359
360      /// Constructor.
361      ///
362      MinCutNodeIt(GomoryHu const &gomory,
363                   ///< The GomoryHu class. You must call its
364                   ///  run() method
365                   ///  before initializing this iterator.
366                   const Node& s, ///< The base node.
367                   const Node& t,
368                   ///< The node you want to separate from node \c s.
369                   bool side=true
370                   ///< If it is \c true (default) then the iterator lists
371                   ///  the nodes of the component containing \c s,
372                   ///  otherwise it lists the other component.
373                   /// \note As the minimum cut is not always unique,
374                   /// \code
375                   /// MinCutNodeIt(gomory, s, t, true);
376                   /// \endcode
377                   /// and
378                   /// \code
379                   /// MinCutNodeIt(gomory, t, s, false);
380                   /// \endcode
381                   /// does not necessarily give the same set of nodes.
382                   /// However it is ensured that
383                   /// \code
384                   /// MinCutNodeIt(gomory, s, t, true);
385                   /// \endcode
386                   /// and
387                   /// \code
388                   /// MinCutNodeIt(gomory, s, t, false);
389                   /// \endcode
390                   /// together list each node exactly once.
391                   )
392        : _side(side), _cut(gomory._graph)
393      {
394        gomory.minCutMap(s,t,_cut);
395        for(_node_it=typename Graph::NodeIt(gomory._graph);
396            _node_it!=INVALID && _cut[_node_it]!=_side;
397            ++_node_it) {}
398      }
399      /// Conversion to \c Node
400
401      /// Conversion to \c Node.
402      ///
403      operator typename Graph::Node() const
404      {
405        return _node_it;
406      }
407      bool operator==(Invalid) { return _node_it==INVALID; }
408      bool operator!=(Invalid) { return _node_it!=INVALID; }
409      /// Next node
410
411      /// Next node.
412      ///
413      MinCutNodeIt &operator++()
414      {
415        for(++_node_it;_node_it!=INVALID&&_cut[_node_it]!=_side;++_node_it) {}
416        return *this;
417      }
418      /// Postfix incrementation
419
420      /// Postfix incrementation.
421      ///
422      /// \warning This incrementation
423      /// returns a \c Node, not a \c MinCutNodeIt, as one may
424      /// expect.
425      typename Graph::Node operator++(int)
426      {
427        typename Graph::Node n=*this;
428        ++(*this);
429        return n;
430      }
431    };
432   
433    friend class MinCutEdgeIt;
434   
435    /// Iterate on the edges of a minimum cut
436   
437    /// This iterator class lists the edges of a minimum cut found by
438    /// GomoryHu. Before using it, you must allocate a GomoryHu class,
439    /// and call its \ref GomoryHu::run() "run()" method.
440    ///
441    /// This example computes the value of the minimum cut separating \c s from
442    /// \c t.
443    /// \code
444    /// GomoruHu<Graph> gom(g, capacities);
445    /// gom.run();
446    /// int value=0;
447    /// for(GomoruHu<Graph>::MinCutEdgeIt e(gom,s,t); e!=INVALID; ++e)
448    ///   value+=capacities[e];
449    /// \endcode
450    /// the result will be the same as it is returned by
451    /// \ref GomoryHu::minCutValue() "gom.minCutValue(s,t)"
452    class MinCutEdgeIt
453    {
454      bool _side;
455      const Graph &_graph;
456      typename Graph::NodeIt _node_it;
457      typename Graph::OutArcIt _arc_it;
458      typename Graph::template NodeMap<bool> _cut;
459      void step()
460      {
461        ++_arc_it;
462        while(_node_it!=INVALID && _arc_it==INVALID)
463          {
464            for(++_node_it;_node_it!=INVALID&&!_cut[_node_it];++_node_it) {}
465            if(_node_it!=INVALID)
466              _arc_it=typename Graph::OutArcIt(_graph,_node_it);
467          }
468      }
469     
470    public:
471      MinCutEdgeIt(GomoryHu const &gomory,
472                   ///< The GomoryHu class. You must call its
473                   ///  run() method
474                   ///  before initializing this iterator.
475                   const Node& s,  ///< The base node.
476                   const Node& t,
477                   ///< The node you want to separate from node \c s.
478                   bool side=true
479                   ///< If it is \c true (default) then the listed arcs
480                   ///  will be oriented from the
481                   ///  the nodes of the component containing \c s,
482                   ///  otherwise they will be oriented in the opposite
483                   ///  direction.
484                   )
485        : _graph(gomory._graph), _cut(_graph)
486      {
487        gomory.minCutMap(s,t,_cut);
488        if(!side)
489          for(typename Graph::NodeIt n(_graph);n!=INVALID;++n)
490            _cut[n]=!_cut[n];
491
492        for(_node_it=typename Graph::NodeIt(_graph);
493            _node_it!=INVALID && !_cut[_node_it];
494            ++_node_it) {}
495        _arc_it = _node_it!=INVALID ?
496          typename Graph::OutArcIt(_graph,_node_it) : INVALID;
497        while(_node_it!=INVALID && _arc_it == INVALID)
498          {
499            for(++_node_it; _node_it!=INVALID&&!_cut[_node_it]; ++_node_it) {}
500            if(_node_it!=INVALID)
501              _arc_it= typename Graph::OutArcIt(_graph,_node_it);
502          }
503        while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
504      }
505      /// Conversion to \c Arc
506
507      /// Conversion to \c Arc.
508      ///
509      operator typename Graph::Arc() const
510      {
511        return _arc_it;
512      }
513      /// Conversion to \c Edge
514
515      /// Conversion to \c Edge.
516      ///
517      operator typename Graph::Edge() const
518      {
519        return _arc_it;
520      }
521      bool operator==(Invalid) { return _node_it==INVALID; }
522      bool operator!=(Invalid) { return _node_it!=INVALID; }
523      /// Next edge
524
525      /// Next edge.
526      ///
527      MinCutEdgeIt &operator++()
528      {
529        step();
530        while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
531        return *this;
532      }
533      /// Postfix incrementation
534     
535      /// Postfix incrementation.
536      ///
537      /// \warning This incrementation
538      /// returns an \c Arc, not a \c MinCutEdgeIt, as one may expect.
539      typename Graph::Arc operator++(int)
540      {
541        typename Graph::Arc e=*this;
542        ++(*this);
543        return e;
544      }
545    };
546
547  };
548
549}
550
551#endif
Note: See TracBrowser for help on using the repository browser.