lemon/gomory_hu.h
branch1.1
changeset 782 1309a803a057
parent 588 293551ad254f
equal deleted inserted replaced
3:02c9547e1525 4:92759760dc0b
     1 /* -*- C++ -*-
     1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
     2  *
     2  *
     3  * This file is a part of LEMON, a generic C++ optimization library
     3  * This file is a part of LEMON, a generic C++ optimization library.
     4  *
     4  *
     5  * Copyright (C) 2003-2008
     5  * Copyright (C) 2003-2011
     6  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     6  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7  * (Egervary Research Group on Combinatorial Optimization, EGRES).
     7  * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8  *
     8  *
     9  * Permission to use, modify and distribute this software is granted
     9  * Permission to use, modify and distribute this software is granted
    10  * provided that this copyright notice appears in all copies. For
    10  * provided that this copyright notice appears in all copies. For
    25 #include <lemon/preflow.h>
    25 #include <lemon/preflow.h>
    26 #include <lemon/concept_check.h>
    26 #include <lemon/concept_check.h>
    27 #include <lemon/concepts/maps.h>
    27 #include <lemon/concepts/maps.h>
    28 
    28 
    29 /// \ingroup min_cut
    29 /// \ingroup min_cut
    30 /// \file 
    30 /// \file
    31 /// \brief Gomory-Hu cut tree in graphs.
    31 /// \brief Gomory-Hu cut tree in graphs.
    32 
    32 
    33 namespace lemon {
    33 namespace lemon {
    34 
    34 
    35   /// \ingroup min_cut
    35   /// \ingroup min_cut
    36   ///
    36   ///
    37   /// \brief Gomory-Hu cut tree algorithm
    37   /// \brief Gomory-Hu cut tree algorithm
    38   ///
    38   ///
    39   /// The Gomory-Hu tree is a tree on the node set of a given graph, but it
    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
    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 
    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
    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
    43   /// between these nodes. Moreover the components obtained by removing
    44   /// this edge from the tree determine the corresponding minimum cut.
    44   /// this edge from the tree determine the corresponding minimum cut.
    45   /// Therefore once this tree is computed, the minimum cut between any pair
    45   /// Therefore once this tree is computed, the minimum cut between any pair
    46   /// of nodes can easily be obtained.
    46   /// of nodes can easily be obtained.
    47   /// 
    47   ///
    48   /// The algorithm calculates \e n-1 distinct minimum cuts (currently with
    48   /// The algorithm calculates \e n-1 distinct minimum cuts (currently with
    49   /// the \ref Preflow algorithm), thus it has \f$O(n^3\sqrt{e})\f$ overall
    49   /// the \ref Preflow algorithm), thus it has \f$O(n^3\sqrt{e})\f$ overall
    50   /// time complexity. It calculates a rooted Gomory-Hu tree.
    50   /// time complexity. It calculates a rooted Gomory-Hu tree.
    51   /// The structure of the tree and the edge weights can be
    51   /// The structure of the tree and the edge weights can be
    52   /// obtained using \c predNode(), \c predValue() and \c rootDist().
    52   /// obtained using \c predNode(), \c predValue() and \c rootDist().
    58   /// \tparam GR The type of the undirected graph the algorithm runs on.
    58   /// \tparam GR The type of the undirected graph the algorithm runs on.
    59   /// \tparam CAP The type of the edge map containing the capacities.
    59   /// \tparam CAP The type of the edge map containing the capacities.
    60   /// The default map type is \ref concepts::Graph::EdgeMap "GR::EdgeMap<int>".
    60   /// The default map type is \ref concepts::Graph::EdgeMap "GR::EdgeMap<int>".
    61 #ifdef DOXYGEN
    61 #ifdef DOXYGEN
    62   template <typename GR,
    62   template <typename GR,
    63 	    typename CAP>
    63             typename CAP>
    64 #else
    64 #else
    65   template <typename GR,
    65   template <typename GR,
    66 	    typename CAP = typename GR::template EdgeMap<int> >
    66             typename CAP = typename GR::template EdgeMap<int> >
    67 #endif
    67 #endif
    68   class GomoryHu {
    68   class GomoryHu {
    69   public:
    69   public:
    70 
    70 
    71     /// The graph type of the algorithm
    71     /// The graph type of the algorithm
    72     typedef GR Graph;
    72     typedef GR Graph;
    73     /// The capacity map type of the algorithm
    73     /// The capacity map type of the algorithm
    74     typedef CAP Capacity;
    74     typedef CAP Capacity;
    75     /// The value type of capacities
    75     /// The value type of capacities
    76     typedef typename Capacity::Value Value;
    76     typedef typename Capacity::Value Value;
    77     
    77 
    78   private:
    78   private:
    79 
    79 
    80     TEMPLATE_GRAPH_TYPEDEFS(Graph);
    80     TEMPLATE_GRAPH_TYPEDEFS(Graph);
    81 
    81 
    82     const Graph& _graph;
    82     const Graph& _graph;
    87     typename Graph::template NodeMap<Value>* _weight;
    87     typename Graph::template NodeMap<Value>* _weight;
    88     typename Graph::template NodeMap<int>* _order;
    88     typename Graph::template NodeMap<int>* _order;
    89 
    89 
    90     void createStructures() {
    90     void createStructures() {
    91       if (!_pred) {
    91       if (!_pred) {
    92 	_pred = new typename Graph::template NodeMap<Node>(_graph);
    92         _pred = new typename Graph::template NodeMap<Node>(_graph);
    93       }
    93       }
    94       if (!_weight) {
    94       if (!_weight) {
    95 	_weight = new typename Graph::template NodeMap<Value>(_graph);
    95         _weight = new typename Graph::template NodeMap<Value>(_graph);
    96       }
    96       }
    97       if (!_order) {
    97       if (!_order) {
    98 	_order = new typename Graph::template NodeMap<int>(_graph);
    98         _order = new typename Graph::template NodeMap<int>(_graph);
    99       }
    99       }
   100     }
   100     }
   101 
   101 
   102     void destroyStructures() {
   102     void destroyStructures() {
   103       if (_pred) {
   103       if (_pred) {
   104 	delete _pred;
   104         delete _pred;
   105       }
   105       }
   106       if (_weight) {
   106       if (_weight) {
   107 	delete _weight;
   107         delete _weight;
   108       }
   108       }
   109       if (_order) {
   109       if (_order) {
   110 	delete _order;
   110         delete _order;
   111       }
   111       }
   112     }
   112     }
   113   
   113 
   114   public:
   114   public:
   115 
   115 
   116     /// \brief Constructor
   116     /// \brief Constructor
   117     ///
   117     ///
   118     /// Constructor.
   118     /// Constructor.
   119     /// \param graph The undirected graph the algorithm runs on.
   119     /// \param graph The undirected graph the algorithm runs on.
   120     /// \param capacity The edge capacity map.
   120     /// \param capacity The edge capacity map.
   121     GomoryHu(const Graph& graph, const Capacity& capacity) 
   121     GomoryHu(const Graph& graph, const Capacity& capacity)
   122       : _graph(graph), _capacity(capacity),
   122       : _graph(graph), _capacity(capacity),
   123 	_pred(0), _weight(0), _order(0) 
   123         _pred(0), _weight(0), _order(0)
   124     {
   124     {
   125       checkConcept<concepts::ReadMap<Edge, Value>, Capacity>();
   125       checkConcept<concepts::ReadMap<Edge, Value>, Capacity>();
   126     }
   126     }
   127 
   127 
   128 
   128 
   132     ~GomoryHu() {
   132     ~GomoryHu() {
   133       destroyStructures();
   133       destroyStructures();
   134     }
   134     }
   135 
   135 
   136   private:
   136   private:
   137   
   137 
   138     // Initialize the internal data structures
   138     // Initialize the internal data structures
   139     void init() {
   139     void init() {
   140       createStructures();
   140       createStructures();
   141 
   141 
   142       _root = NodeIt(_graph);
   142       _root = NodeIt(_graph);
   143       for (NodeIt n(_graph); n != INVALID; ++n) {
   143       for (NodeIt n(_graph); n != INVALID; ++n) {
   144         (*_pred)[n] = _root;
   144         (*_pred)[n] = _root;
   145         (*_order)[n] = -1;
   145         (*_order)[n] = -1;
   146       }
   146       }
   147       (*_pred)[_root] = INVALID;
   147       (*_pred)[_root] = INVALID;
   148       (*_weight)[_root] = std::numeric_limits<Value>::max(); 
   148       (*_weight)[_root] = std::numeric_limits<Value>::max();
   149     }
   149     }
   150 
   150 
   151 
   151 
   152     // Start the algorithm
   152     // Start the algorithm
   153     void start() {
   153     void start() {
   154       Preflow<Graph, Capacity> fa(_graph, _capacity, _root, INVALID);
   154       Preflow<Graph, Capacity> fa(_graph, _capacity, _root, INVALID);
   155 
   155 
   156       for (NodeIt n(_graph); n != INVALID; ++n) {
   156       for (NodeIt n(_graph); n != INVALID; ++n) {
   157 	if (n == _root) continue;
   157         if (n == _root) continue;
   158 
   158 
   159 	Node pn = (*_pred)[n];
   159         Node pn = (*_pred)[n];
   160 	fa.source(n);
   160         fa.source(n);
   161 	fa.target(pn);
   161         fa.target(pn);
   162 
   162 
   163 	fa.runMinCut();
   163         fa.runMinCut();
   164 
   164 
   165 	(*_weight)[n] = fa.flowValue();
   165         (*_weight)[n] = fa.flowValue();
   166 
   166 
   167 	for (NodeIt nn(_graph); nn != INVALID; ++nn) {
   167         for (NodeIt nn(_graph); nn != INVALID; ++nn) {
   168 	  if (nn != n && fa.minCut(nn) && (*_pred)[nn] == pn) {
   168           if (nn != n && fa.minCut(nn) && (*_pred)[nn] == pn) {
   169 	    (*_pred)[nn] = n;
   169             (*_pred)[nn] = n;
   170 	  }
   170           }
   171 	}
   171         }
   172 	if ((*_pred)[pn] != INVALID && fa.minCut((*_pred)[pn])) {
   172         if ((*_pred)[pn] != INVALID && fa.minCut((*_pred)[pn])) {
   173 	  (*_pred)[n] = (*_pred)[pn];
   173           (*_pred)[n] = (*_pred)[pn];
   174 	  (*_pred)[pn] = n;
   174           (*_pred)[pn] = n;
   175 	  (*_weight)[n] = (*_weight)[pn];
   175           (*_weight)[n] = (*_weight)[pn];
   176 	  (*_weight)[pn] = fa.flowValue();
   176           (*_weight)[pn] = fa.flowValue();
   177 	}
   177         }
   178       }
   178       }
   179 
   179 
   180       (*_order)[_root] = 0;
   180       (*_order)[_root] = 0;
   181       int index = 1;
   181       int index = 1;
   182 
   182 
   183       for (NodeIt n(_graph); n != INVALID; ++n) {
   183       for (NodeIt n(_graph); n != INVALID; ++n) {
   184 	std::vector<Node> st;
   184         std::vector<Node> st;
   185 	Node nn = n;
   185         Node nn = n;
   186 	while ((*_order)[nn] == -1) {
   186         while ((*_order)[nn] == -1) {
   187 	  st.push_back(nn);
   187           st.push_back(nn);
   188 	  nn = (*_pred)[nn];
   188           nn = (*_pred)[nn];
   189 	}
   189         }
   190 	while (!st.empty()) {
   190         while (!st.empty()) {
   191 	  (*_order)[st.back()] = index++;
   191           (*_order)[st.back()] = index++;
   192 	  st.pop_back();
   192           st.pop_back();
   193 	}
   193         }
   194       }
   194       }
   195     }
   195     }
   196 
   196 
   197   public:
   197   public:
   198 
   198 
   199     ///\name Execution Control
   199     ///\name Execution Control
   200  
   200 
   201     ///@{
   201     ///@{
   202 
   202 
   203     /// \brief Run the Gomory-Hu algorithm.
   203     /// \brief Run the Gomory-Hu algorithm.
   204     ///
   204     ///
   205     /// This function runs the Gomory-Hu algorithm.
   205     /// This function runs the Gomory-Hu algorithm.
   206     void run() {
   206     void run() {
   207       init();
   207       init();
   208       start();
   208       start();
   209     }
   209     }
   210     
   210 
   211     /// @}
   211     /// @}
   212 
   212 
   213     ///\name Query Functions
   213     ///\name Query Functions
   214     ///The results of the algorithm can be obtained using these
   214     ///The results of the algorithm can be obtained using these
   215     ///functions.\n
   215     ///functions.\n
   230     }
   230     }
   231 
   231 
   232     /// \brief Return the weight of the predecessor edge in the
   232     /// \brief Return the weight of the predecessor edge in the
   233     /// Gomory-Hu tree.
   233     /// Gomory-Hu tree.
   234     ///
   234     ///
   235     /// This function returns the weight of the predecessor edge of the 
   235     /// This function returns the weight of the predecessor edge of the
   236     /// given node in the Gomory-Hu tree.
   236     /// given node in the Gomory-Hu tree.
   237     /// If \c node is the root of the tree, the result is undefined.
   237     /// If \c node is the root of the tree, the result is undefined.
   238     ///
   238     ///
   239     /// \pre \ref run() must be called before using this function.
   239     /// \pre \ref run() must be called before using this function.
   240     Value predValue(const Node& node) const {
   240     Value predValue(const Node& node) const {
   252     }
   252     }
   253 
   253 
   254     /// \brief Return the minimum cut value between two nodes
   254     /// \brief Return the minimum cut value between two nodes
   255     ///
   255     ///
   256     /// This function returns the minimum cut value between the nodes
   256     /// This function returns the minimum cut value between the nodes
   257     /// \c s and \c t. 
   257     /// \c s and \c t.
   258     /// It finds the nearest common ancestor of the given nodes in the
   258     /// It finds the nearest common ancestor of the given nodes in the
   259     /// Gomory-Hu tree and calculates the minimum weight edge on the
   259     /// Gomory-Hu tree and calculates the minimum weight edge on the
   260     /// paths to the ancestor.
   260     /// paths to the ancestor.
   261     ///
   261     ///
   262     /// \pre \ref run() must be called before using this function.
   262     /// \pre \ref run() must be called before using this function.
   263     Value minCutValue(const Node& s, const Node& t) const {
   263     Value minCutValue(const Node& s, const Node& t) const {
   264       Node sn = s, tn = t;
   264       Node sn = s, tn = t;
   265       Value value = std::numeric_limits<Value>::max();
   265       Value value = std::numeric_limits<Value>::max();
   266       
   266 
   267       while (sn != tn) {
   267       while (sn != tn) {
   268 	if ((*_order)[sn] < (*_order)[tn]) {
   268         if ((*_order)[sn] < (*_order)[tn]) {
   269 	  if ((*_weight)[tn] <= value) value = (*_weight)[tn];
   269           if ((*_weight)[tn] <= value) value = (*_weight)[tn];
   270 	  tn = (*_pred)[tn];
   270           tn = (*_pred)[tn];
   271 	} else {
   271         } else {
   272 	  if ((*_weight)[sn] <= value) value = (*_weight)[sn];
   272           if ((*_weight)[sn] <= value) value = (*_weight)[sn];
   273 	  sn = (*_pred)[sn];
   273           sn = (*_pred)[sn];
   274 	}
   274         }
   275       }
   275       }
   276       return value;
   276       return value;
   277     }
   277     }
   278 
   278 
   279     /// \brief Return the minimum cut between two nodes
   279     /// \brief Return the minimum cut between two nodes
   292     ///
   292     ///
   293     /// \return The value of the minimum cut between \c s and \c t.
   293     /// \return The value of the minimum cut between \c s and \c t.
   294     ///
   294     ///
   295     /// \pre \ref run() must be called before using this function.
   295     /// \pre \ref run() must be called before using this function.
   296     template <typename CutMap>
   296     template <typename CutMap>
   297     Value minCutMap(const Node& s, ///< 
   297     Value minCutMap(const Node& s, ///<
   298                     const Node& t,
   298                     const Node& t,
   299                     ///< 
   299                     ///<
   300                     CutMap& cutMap
   300                     CutMap& cutMap
   301                     ///< 
   301                     ///<
   302                     ) const {
   302                     ) const {
   303       Node sn = s, tn = t;
   303       Node sn = s, tn = t;
   304       bool s_root=false;
   304       bool s_root=false;
   305       Node rn = INVALID;
   305       Node rn = INVALID;
   306       Value value = std::numeric_limits<Value>::max();
   306       Value value = std::numeric_limits<Value>::max();
   307       
   307 
   308       while (sn != tn) {
   308       while (sn != tn) {
   309 	if ((*_order)[sn] < (*_order)[tn]) {
   309         if ((*_order)[sn] < (*_order)[tn]) {
   310 	  if ((*_weight)[tn] <= value) {
   310           if ((*_weight)[tn] <= value) {
   311 	    rn = tn;
   311             rn = tn;
   312             s_root = false;
   312             s_root = false;
   313 	    value = (*_weight)[tn];
   313             value = (*_weight)[tn];
   314 	  }
   314           }
   315 	  tn = (*_pred)[tn];
   315           tn = (*_pred)[tn];
   316 	} else {
   316         } else {
   317 	  if ((*_weight)[sn] <= value) {
   317           if ((*_weight)[sn] <= value) {
   318 	    rn = sn;
   318             rn = sn;
   319             s_root = true;
   319             s_root = true;
   320 	    value = (*_weight)[sn];
   320             value = (*_weight)[sn];
   321 	  }
   321           }
   322 	  sn = (*_pred)[sn];
   322           sn = (*_pred)[sn];
   323 	}
   323         }
   324       }
   324       }
   325 
   325 
   326       typename Graph::template NodeMap<bool> reached(_graph, false);
   326       typename Graph::template NodeMap<bool> reached(_graph, false);
   327       reached[_root] = true;
   327       reached[_root] = true;
   328       cutMap.set(_root, !s_root);
   328       cutMap.set(_root, !s_root);
   329       reached[rn] = true;
   329       reached[rn] = true;
   330       cutMap.set(rn, s_root);
   330       cutMap.set(rn, s_root);
   331 
   331 
   332       std::vector<Node> st;
   332       std::vector<Node> st;
   333       for (NodeIt n(_graph); n != INVALID; ++n) {
   333       for (NodeIt n(_graph); n != INVALID; ++n) {
   334 	st.clear();
   334         st.clear();
   335         Node nn = n;
   335         Node nn = n;
   336 	while (!reached[nn]) {
   336         while (!reached[nn]) {
   337 	  st.push_back(nn);
   337           st.push_back(nn);
   338 	  nn = (*_pred)[nn];
   338           nn = (*_pred)[nn];
   339 	}
   339         }
   340 	while (!st.empty()) {
   340         while (!st.empty()) {
   341 	  cutMap.set(st.back(), cutMap[nn]);
   341           cutMap.set(st.back(), cutMap[nn]);
   342 	  st.pop_back();
   342           st.pop_back();
   343 	}
   343         }
   344       }
   344       }
   345       
   345 
   346       return value;
   346       return value;
   347     }
   347     }
   348 
   348 
   349     ///@}
   349     ///@}
   350 
   350 
   351     friend class MinCutNodeIt;
   351     friend class MinCutNodeIt;
   352 
   352 
   353     /// Iterate on the nodes of a minimum cut
   353     /// Iterate on the nodes of a minimum cut
   354     
   354 
   355     /// This iterator class lists the nodes of a minimum cut found by
   355     /// This iterator class lists the nodes of a minimum cut found by
   356     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   356     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   357     /// and call its \ref GomoryHu::run() "run()" method.
   357     /// and call its \ref GomoryHu::run() "run()" method.
   358     ///
   358     ///
   359     /// This example counts the nodes in the minimum cut separating \c s from
   359     /// This example counts the nodes in the minimum cut separating \c s from
   442         typename Graph::Node n=*this;
   442         typename Graph::Node n=*this;
   443         ++(*this);
   443         ++(*this);
   444         return n;
   444         return n;
   445       }
   445       }
   446     };
   446     };
   447     
   447 
   448     friend class MinCutEdgeIt;
   448     friend class MinCutEdgeIt;
   449     
   449 
   450     /// Iterate on the edges of a minimum cut
   450     /// Iterate on the edges of a minimum cut
   451     
   451 
   452     /// This iterator class lists the edges of a minimum cut found by
   452     /// This iterator class lists the edges of a minimum cut found by
   453     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   453     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   454     /// and call its \ref GomoryHu::run() "run()" method.
   454     /// and call its \ref GomoryHu::run() "run()" method.
   455     ///
   455     ///
   456     /// This example computes the value of the minimum cut separating \c s from
   456     /// This example computes the value of the minimum cut separating \c s from
   479             for(++_node_it;_node_it!=INVALID&&!_cut[_node_it];++_node_it) {}
   479             for(++_node_it;_node_it!=INVALID&&!_cut[_node_it];++_node_it) {}
   480             if(_node_it!=INVALID)
   480             if(_node_it!=INVALID)
   481               _arc_it=typename Graph::OutArcIt(_graph,_node_it);
   481               _arc_it=typename Graph::OutArcIt(_graph,_node_it);
   482           }
   482           }
   483       }
   483       }
   484       
   484 
   485     public:
   485     public:
   486       /// Constructor
   486       /// Constructor
   487 
   487 
   488       /// Constructor.
   488       /// Constructor.
   489       ///
   489       ///
   548         step();
   548         step();
   549         while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
   549         while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
   550         return *this;
   550         return *this;
   551       }
   551       }
   552       /// Postfix incrementation
   552       /// Postfix incrementation
   553       
   553 
   554       /// Postfix incrementation.
   554       /// Postfix incrementation.
   555       ///
   555       ///
   556       /// \warning This incrementation
   556       /// \warning This incrementation
   557       /// returns an \c Arc, not a \c MinCutEdgeIt, as one may expect.
   557       /// returns an \c Arc, not a \c MinCutEdgeIt, as one may expect.
   558       typename Graph::Arc operator++(int)
   558       typename Graph::Arc operator++(int)