COIN-OR::LEMON - Graph Library

Ticket #66: 33085d6328e4.patch

File 33085d6328e4.patch, 10.4 KB (added by Tapolcai János, 11 years ago)

ported from lemon 0.7

  • new file lemon/gomory_hu_tree.h

    # HG changeset patch
    # User Janos Tapolcai <tapolcai@tmit.bme.hu>
    # Date 1235146637 -3600
    # Node ID 33085d6328e470d6ad8cc11f1616bebb4f08066c
    # Parent  b9b3473327e344ba8647be9f8c674644892e25d1
    Proting Gomory-Hu algorithm
    
    diff --git a/lemon/gomory_hu_tree.h b/lemon/gomory_hu_tree.h
    new file mode 100644
    - +  
     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 <lemon/preflow.h>
     23#include <lemon/concept_check.h>
     24#include <lemon/concepts/maps.h>
     25
     26/// \ingroup min_cut
     27/// \file
     28/// \brief Gomory-Hu cut tree in graphs.
     29
     30namespace lemon {
     31
     32  /// \ingroup min_cut
     33  ///
     34  /// \brief Gomory-Hu cut tree algorithm
     35  ///
     36  /// The Gomory-Hu tree is a tree on the nodeset of the digraph, but it
     37  /// may contain arcs which are not in the original digraph. It helps
     38  /// to calculate the minimum cut between all pairs of nodes, because
     39  /// the minimum capacity arc on the tree path between two nodes has
     40  /// the same weight as the minimum cut in the digraph between these
     41  /// nodes. Moreover this arc separates the nodes to two parts which
     42  /// determine this minimum cut.
     43  ///
     44  /// The algorithm calculates \e n-1 distinict minimum cuts with
     45  /// preflow algorithm, therefore the algorithm has
     46  /// \f$(O(n^3\sqrt{e})\f$ overall time complexity. It calculates a
     47  /// rooted Gomory-Hu tree, the structure of the tree and the weights
     48  /// can be obtained with \c predNode() and \c predValue()
     49  /// functions. The \c minCutValue() and \c minCutMap() calculates
     50  /// the minimum cut and the minimum cut value between any two node
     51  /// in the digraph.
     52  template <typename _Graph,
     53            typename _Capacity = typename _Graph::template EdgeMap<int> >
     54  class GomoryHuTree {
     55  public:
     56
     57    /// The graph type
     58    typedef _Graph Graph;
     59    /// The capacity on edges
     60    typedef _Capacity Capacity;
     61    /// The value type of capacities
     62    typedef typename Capacity::Value Value;
     63   
     64  private:
     65
     66    TEMPLATE_GRAPH_TYPEDEFS(Graph);
     67
     68    const Graph& _graph;
     69    const Capacity& _capacity;
     70
     71    Node _root;
     72    typename Graph::template NodeMap<Node>* _pred;
     73    typename Graph::template NodeMap<Value>* _weight;
     74    typename Graph::template NodeMap<int>* _order;
     75
     76    void createStructures() {
     77      if (!_pred) {
     78        _pred = new typename Graph::template NodeMap<Node>(_graph);
     79      }
     80      if (!_weight) {
     81        _weight = new typename Graph::template NodeMap<Value>(_graph);
     82      }
     83      if (!_order) {
     84        _order = new typename Graph::template NodeMap<int>(_graph);
     85      }
     86    }
     87
     88    void destroyStructures() {
     89      if (_pred) {
     90        delete _pred;
     91      }
     92      if (_weight) {
     93        delete _weight;
     94      }
     95      if (_order) {
     96        delete _order;
     97      }
     98    }
     99 
     100  public:
     101
     102    /// \brief Constructor
     103    ///
     104    /// Constructor
     105    /// \param graph The graph type.
     106    /// \param capacity The capacity map.
     107    GomoryHuTree(const Graph& graph, const Capacity& capacity)
     108      : _graph(graph), _capacity(capacity),
     109        _pred(0), _weight(0), _order(0)
     110    {
     111      checkConcept<concepts::ReadMap<Edge, Value>, Capacity>();
     112    }
     113
     114
     115    /// \brief Destructor
     116    ///
     117    /// Destructor
     118    ~GomoryHuTree() {
     119      destroyStructures();
     120    }
     121
     122    /// \brief Initializes the internal data structures.
     123    ///
     124    /// Initializes the internal data structures.
     125    ///
     126    void init() {
     127      createStructures();
     128
     129      _root = NodeIt(_graph);
     130      for (NodeIt n(_graph); n != INVALID; ++n) {
     131        _pred->set(n, _root);
     132        _order->set(n, -1);
     133      }
     134      _pred->set(_root, INVALID);
     135      _weight->set(_root, std::numeric_limits<Value>::max());
     136    }
     137
     138
     139    /// \brief Starts the algorithm
     140    ///
     141    /// Starts the algorithm.
     142    void start() {
     143      Preflow<Graph, Capacity> fa(_graph, _capacity, _root, INVALID);
     144
     145      for (NodeIt n(_graph); n != INVALID; ++n) {
     146        if (n == _root) continue;
     147
     148        Node pn = (*_pred)[n];
     149        fa.source(n);
     150        fa.target(pn);
     151
     152        fa.runMinCut();
     153
     154        _weight->set(n, fa.flowValue());
     155
     156        for (NodeIt nn(_graph); nn != INVALID; ++nn) {
     157          if (nn != n && fa.minCut(nn) && (*_pred)[nn] == pn) {
     158            _pred->set(nn, n);
     159          }
     160        }
     161        if ((*_pred)[pn] != INVALID && fa.minCut((*_pred)[pn])) {
     162          _pred->set(n, (*_pred)[pn]);
     163          _pred->set(pn, n);
     164          _weight->set(n, (*_weight)[pn]);
     165          _weight->set(pn, fa.flowValue());     
     166        }
     167      }
     168
     169      _order->set(_root, 0);
     170      int index = 1;
     171
     172      for (NodeIt n(_graph); n != INVALID; ++n) {
     173        std::vector<Node> st;
     174        Node nn = n;
     175        while ((*_order)[nn] == -1) {
     176          st.push_back(nn);
     177          nn = (*_pred)[nn];
     178        }
     179        while (!st.empty()) {
     180          _order->set(st.back(), index++);
     181          st.pop_back();
     182        }
     183      }
     184    }
     185
     186    /// \brief Runs the Gomory-Hu algorithm. 
     187    ///
     188    /// Runs the Gomory-Hu algorithm.
     189    /// \note gh.run() is just a shortcut of the following code.
     190    /// \code
     191    ///   ght.init();
     192    ///   ght.start();
     193    /// \endcode
     194    void run() {
     195      init();
     196      start();
     197    }
     198
     199    /// \brief Returns the predecessor node in the Gomory-Hu tree.
     200    ///
     201    /// Returns the predecessor node in the Gomory-Hu tree. If the node is
     202    /// the root of the Gomory-Hu tree, then it returns \c INVALID.
     203    Node predNode(const Node& node) {
     204      return (*_pred)[node];
     205    }
     206
     207    /// \brief Returns the weight of the predecessor arc in the
     208    /// Gomory-Hu tree.
     209    ///
     210    /// Returns the weight of the predecessor arc in the Gomory-Hu
     211    /// tree.  If the node is the root of the Gomory-Hu tree, the
     212    /// result is undefined.
     213    Value predValue(const Node& node) {
     214      return (*_weight)[node];
     215    }
     216
     217    /// \brief Returns the minimum cut value between two nodes
     218    ///
     219    /// Returns the minimum cut value between two nodes. The
     220    /// algorithm finds the nearest common ancestor in the Gomory-Hu
     221    /// tree and calculates the minimum weight arc on the paths to
     222    /// the ancestor.
     223    Value minCutValue(const Node& s, const Node& t) const {
     224      Node sn = s, tn = t;
     225      Value value = std::numeric_limits<Value>::max();
     226     
     227      while (sn != tn) {
     228        if ((*_order)[sn] < (*_order)[tn]) {
     229          if ((*_weight)[tn] < value) value = (*_weight)[tn];
     230          tn = (*_pred)[tn];
     231        } else {
     232          if ((*_weight)[sn] < value) value = (*_weight)[sn];
     233          sn = (*_pred)[sn];
     234        }
     235      }
     236      return value;
     237    }
     238
     239    /// \brief Returns the minimum cut between two nodes
     240    ///
     241    /// Returns the minimum cut value between two nodes. The
     242    /// algorithm finds the nearest common ancestor in the Gomory-Hu
     243    /// tree and calculates the minimum weight arc on the paths to
     244    /// the ancestor. Then it sets all nodes to the cut determined by
     245    /// this arc. The \c cutMap should be \ref concepts::ReadWriteMap
     246    /// "ReadWriteMap".
     247    template <typename CutMap>
     248    Value minCutMap(const Node& s, const Node& t, CutMap& cutMap) const {
     249      Node sn = s, tn = t;
     250
     251      Node rn = INVALID;
     252      Value value = std::numeric_limits<Value>::max();
     253     
     254      while (sn != tn) {
     255        if ((*_order)[sn] < (*_order)[tn]) {
     256          if ((*_weight)[tn] < value) {
     257            rn = tn;
     258            value = (*_weight)[tn];
     259          }
     260          tn = (*_pred)[tn];
     261        } else {
     262          if ((*_weight)[sn] < value) {
     263            rn = sn;
     264            value = (*_weight)[sn];
     265          }
     266          sn = (*_pred)[sn];
     267        }
     268      }
     269
     270      typename Graph::template NodeMap<bool> reached(_graph, false);
     271      reached.set(_root, true);
     272      cutMap.set(_root, false);
     273      reached.set(rn, true);
     274      cutMap.set(rn, true);
     275
     276      for (NodeIt n(_graph); n != INVALID; ++n) {
     277        std::vector<Node> st;
     278        Node nn = n;
     279        while (!reached[nn]) {
     280          st.push_back(nn);
     281          nn = (*_pred)[nn];
     282        }
     283        while (!st.empty()) {
     284          cutMap.set(st.back(), cutMap[nn]);
     285          st.pop_back();
     286        }
     287      }
     288     
     289      return value;
     290    }
     291
     292  };
     293
     294}
     295
     296#endif
  • new file test/gomory_hu_test.cc

    diff --git a/test/gomory_hu_test.cc b/test/gomory_hu_test.cc
    new file mode 100644
    - +  
     1#include <iostream>
     2
     3#include "test_tools.h"
     4#include <lemon/smart_graph.h>
     5#include <lemon/adaptors.h>
     6#include <lemon/lgf_reader.h>
     7#include <lemon/lgf_writer.h>
     8#include <lemon/dimacs.h>
     9#include <lemon/time_measure.h>
     10#include <lemon/gomory_hu_tree.h>
     11#include <cstdlib>
     12
     13using namespace std;
     14using namespace lemon;
     15
     16typedef SmartGraph Graph;
     17
     18char test_lgf[] =
     19  "@nodes\n"
     20  "label\n"
     21  "0\n"
     22  "1\n"
     23  "2\n"
     24  "3\n"
     25  "4\n"
     26  "@arcs\n"
     27  "     label length\n"
     28  "0 1  0     1\n"
     29  "1 2  1     1\n"
     30  "2 3  2     1\n"
     31  "0 3  4     5\n"
     32  "0 3  5     10\n"
     33  "0 3  6     7\n"
     34  "4 2  7     1\n"
     35  "@attributes\n"
     36  "source 0\n"
     37  "target 3\n";
     38 
     39GRAPH_TYPEDEFS(Graph);
     40typedef Graph::EdgeMap<int> IntEdgeMap;
     41typedef Graph::NodeMap<bool> BoolNodeMap;
     42
     43int cutValue(const Graph& graph, const BoolNodeMap& cut,
     44             const IntEdgeMap& capacity) {
     45
     46  int sum = 0;
     47  for (EdgeIt e(graph); e != INVALID; ++e) {
     48    Node s = graph.u(e);
     49    Node t = graph.v(e);
     50
     51    if (cut[s] && !cut[t]) {
     52      sum += capacity[e];
     53    }
     54  }
     55  return sum;
     56}
     57
     58
     59int main(int argc, const char *argv[]) {
     60  Graph graph;
     61  IntEdgeMap capacity(graph);
     62
     63  std::istringstream input(test_lgf);
     64  GraphReader<Graph>(graph, input).
     65    edgeMap("capacity", capacity).run();
     66
     67  GomoryHuTree<Graph> ght(graph, capacity);
     68  ght.init();
     69  ght.run();
     70
     71
     72  int cnt = 0;
     73  std::cerr << countNodes(graph) << std::endl;
     74  for (NodeIt u(graph); u != INVALID; ++u) {
     75    std::cerr << ++cnt << std::endl;
     76    int bcnt = 0;
     77    for (NodeIt v(graph); v != u; ++v) {
     78      Preflow<Graph, IntEdgeMap> pf(graph, capacity, u, v);
     79      pf.runMinCut();
     80      BoolNodeMap cm(graph);
     81      ght.minCutMap(u, v, cm);
     82      check(pf.flowValue() == ght.minCutValue(u, v), "Wrong cut 1");
     83      check(cm[u] != cm[v], "Wrong cut 3");
     84      check(pf.flowValue() == cutValue(graph, cm, capacity), "Wrong cut 2");
     85     
     86    }
     87  }
     88 
     89  return 0;
     90}