lemon/gomory_hu.h
changeset 1034 ef200e268af2
parent 786 e20173729589
child 1080 c5cd8960df74
equal deleted inserted replaced
5:00a797224124 6:3489607882e0
     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-2010
     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
   300                     ) const {
   300                     ) const {
   301       Node sn = s, tn = t;
   301       Node sn = s, tn = t;
   302       bool s_root=false;
   302       bool s_root=false;
   303       Node rn = INVALID;
   303       Node rn = INVALID;
   304       Value value = std::numeric_limits<Value>::max();
   304       Value value = std::numeric_limits<Value>::max();
   305       
   305 
   306       while (sn != tn) {
   306       while (sn != tn) {
   307 	if ((*_order)[sn] < (*_order)[tn]) {
   307         if ((*_order)[sn] < (*_order)[tn]) {
   308 	  if ((*_weight)[tn] <= value) {
   308           if ((*_weight)[tn] <= value) {
   309 	    rn = tn;
   309             rn = tn;
   310             s_root = false;
   310             s_root = false;
   311 	    value = (*_weight)[tn];
   311             value = (*_weight)[tn];
   312 	  }
   312           }
   313 	  tn = (*_pred)[tn];
   313           tn = (*_pred)[tn];
   314 	} else {
   314         } else {
   315 	  if ((*_weight)[sn] <= value) {
   315           if ((*_weight)[sn] <= value) {
   316 	    rn = sn;
   316             rn = sn;
   317             s_root = true;
   317             s_root = true;
   318 	    value = (*_weight)[sn];
   318             value = (*_weight)[sn];
   319 	  }
   319           }
   320 	  sn = (*_pred)[sn];
   320           sn = (*_pred)[sn];
   321 	}
   321         }
   322       }
   322       }
   323 
   323 
   324       typename Graph::template NodeMap<bool> reached(_graph, false);
   324       typename Graph::template NodeMap<bool> reached(_graph, false);
   325       reached[_root] = true;
   325       reached[_root] = true;
   326       cutMap.set(_root, !s_root);
   326       cutMap.set(_root, !s_root);
   327       reached[rn] = true;
   327       reached[rn] = true;
   328       cutMap.set(rn, s_root);
   328       cutMap.set(rn, s_root);
   329 
   329 
   330       std::vector<Node> st;
   330       std::vector<Node> st;
   331       for (NodeIt n(_graph); n != INVALID; ++n) {
   331       for (NodeIt n(_graph); n != INVALID; ++n) {
   332 	st.clear();
   332         st.clear();
   333         Node nn = n;
   333         Node nn = n;
   334 	while (!reached[nn]) {
   334         while (!reached[nn]) {
   335 	  st.push_back(nn);
   335           st.push_back(nn);
   336 	  nn = (*_pred)[nn];
   336           nn = (*_pred)[nn];
   337 	}
   337         }
   338 	while (!st.empty()) {
   338         while (!st.empty()) {
   339 	  cutMap.set(st.back(), cutMap[nn]);
   339           cutMap.set(st.back(), cutMap[nn]);
   340 	  st.pop_back();
   340           st.pop_back();
   341 	}
   341         }
   342       }
   342       }
   343       
   343 
   344       return value;
   344       return value;
   345     }
   345     }
   346 
   346 
   347     ///@}
   347     ///@}
   348 
   348 
   349     friend class MinCutNodeIt;
   349     friend class MinCutNodeIt;
   350 
   350 
   351     /// Iterate on the nodes of a minimum cut
   351     /// Iterate on the nodes of a minimum cut
   352     
   352 
   353     /// This iterator class lists the nodes of a minimum cut found by
   353     /// This iterator class lists the nodes of a minimum cut found by
   354     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   354     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   355     /// and call its \ref GomoryHu::run() "run()" method.
   355     /// and call its \ref GomoryHu::run() "run()" method.
   356     ///
   356     ///
   357     /// This example counts the nodes in the minimum cut separating \c s from
   357     /// This example counts the nodes in the minimum cut separating \c s from
   440         typename Graph::Node n=*this;
   440         typename Graph::Node n=*this;
   441         ++(*this);
   441         ++(*this);
   442         return n;
   442         return n;
   443       }
   443       }
   444     };
   444     };
   445     
   445 
   446     friend class MinCutEdgeIt;
   446     friend class MinCutEdgeIt;
   447     
   447 
   448     /// Iterate on the edges of a minimum cut
   448     /// Iterate on the edges of a minimum cut
   449     
   449 
   450     /// This iterator class lists the edges of a minimum cut found by
   450     /// This iterator class lists the edges of a minimum cut found by
   451     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   451     /// GomoryHu. Before using it, you must allocate a GomoryHu class
   452     /// and call its \ref GomoryHu::run() "run()" method.
   452     /// and call its \ref GomoryHu::run() "run()" method.
   453     ///
   453     ///
   454     /// This example computes the value of the minimum cut separating \c s from
   454     /// This example computes the value of the minimum cut separating \c s from
   477             for(++_node_it;_node_it!=INVALID&&!_cut[_node_it];++_node_it) {}
   477             for(++_node_it;_node_it!=INVALID&&!_cut[_node_it];++_node_it) {}
   478             if(_node_it!=INVALID)
   478             if(_node_it!=INVALID)
   479               _arc_it=typename Graph::OutArcIt(_graph,_node_it);
   479               _arc_it=typename Graph::OutArcIt(_graph,_node_it);
   480           }
   480           }
   481       }
   481       }
   482       
   482 
   483     public:
   483     public:
   484       /// Constructor
   484       /// Constructor
   485 
   485 
   486       /// Constructor.
   486       /// Constructor.
   487       ///
   487       ///
   546         step();
   546         step();
   547         while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
   547         while(_arc_it!=INVALID && _cut[_graph.target(_arc_it)]) step();
   548         return *this;
   548         return *this;
   549       }
   549       }
   550       /// Postfix incrementation
   550       /// Postfix incrementation
   551       
   551 
   552       /// Postfix incrementation.
   552       /// Postfix incrementation.
   553       ///
   553       ///
   554       /// \warning This incrementation
   554       /// \warning This incrementation
   555       /// returns an \c Arc, not a \c MinCutEdgeIt, as one may expect.
   555       /// returns an \c Arc, not a \c MinCutEdgeIt, as one may expect.
   556       typename Graph::Arc operator++(int)
   556       typename Graph::Arc operator++(int)