lemon/hypercube_graph.h
author Peter Kovacs <kpeter@inf.elte.hu>
Thu, 12 Nov 2009 23:26:13 +0100
changeset 806 fa6f37d7a25b
parent 787 c2230649a493
parent 786 e20173729589
permissions -rw-r--r--
Entirely rework CapacityScaling (#180)

- Use the new interface similarly to NetworkSimplex.
- Rework the implementation using an efficient internal structure
for handling the residual network. This improvement made the
code much faster (up to 2-5 times faster on large graphs).
- Handle GEQ supply type (LEQ is not supported).
- Handle negative costs for arcs of finite capacity.
(Note that this algorithm cannot handle arcs of negative cost
and infinite upper bound, thus it returns UNBOUNDED if such
an arc exists.)
- Extend the documentation.
     1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
     2  *
     3  * This file is a part of LEMON, a generic C++ optimization library.
     4  *
     5  * Copyright (C) 2003-2009
     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 HYPERCUBE_GRAPH_H
    20 #define HYPERCUBE_GRAPH_H
    21 
    22 #include <vector>
    23 #include <lemon/core.h>
    24 #include <lemon/assert.h>
    25 #include <lemon/bits/graph_extender.h>
    26 
    27 ///\ingroup graphs
    28 ///\file
    29 ///\brief HypercubeGraph class.
    30 
    31 namespace lemon {
    32 
    33   class HypercubeGraphBase {
    34 
    35   public:
    36 
    37     typedef HypercubeGraphBase Graph;
    38 
    39     class Node;
    40     class Edge;
    41     class Arc;
    42 
    43   public:
    44 
    45     HypercubeGraphBase() {}
    46 
    47   protected:
    48 
    49     void construct(int dim) {
    50       LEMON_ASSERT(dim >= 1, "The number of dimensions must be at least 1.");
    51       _dim = dim;
    52       _node_num = 1 << dim;
    53       _edge_num = dim * (1 << (dim-1));
    54     }
    55 
    56   public:
    57 
    58     typedef True NodeNumTag;
    59     typedef True EdgeNumTag;
    60     typedef True ArcNumTag;
    61 
    62     int nodeNum() const { return _node_num; }
    63     int edgeNum() const { return _edge_num; }
    64     int arcNum() const { return 2 * _edge_num; }
    65 
    66     int maxNodeId() const { return _node_num - 1; }
    67     int maxEdgeId() const { return _edge_num - 1; }
    68     int maxArcId() const { return 2 * _edge_num - 1; }
    69 
    70     static Node nodeFromId(int id) { return Node(id); }
    71     static Edge edgeFromId(int id) { return Edge(id); }
    72     static Arc arcFromId(int id) { return Arc(id); }
    73 
    74     static int id(Node node) { return node._id; }
    75     static int id(Edge edge) { return edge._id; }
    76     static int id(Arc arc) { return arc._id; }
    77 
    78     Node u(Edge edge) const {
    79       int base = edge._id & ((1 << (_dim-1)) - 1);
    80       int k = edge._id >> (_dim-1);
    81       return ((base >> k) << (k+1)) | (base & ((1 << k) - 1));
    82     }
    83 
    84     Node v(Edge edge) const {
    85       int base = edge._id & ((1 << (_dim-1)) - 1);
    86       int k = edge._id >> (_dim-1);
    87       return ((base >> k) << (k+1)) | (base & ((1 << k) - 1)) | (1 << k);
    88     }
    89 
    90     Node source(Arc arc) const {
    91       return (arc._id & 1) == 1 ? u(arc) : v(arc);
    92     }
    93 
    94     Node target(Arc arc) const {
    95       return (arc._id & 1) == 1 ? v(arc) : u(arc);
    96     }
    97 
    98     typedef True FindEdgeTag;
    99     typedef True FindArcTag;
   100 
   101     Edge findEdge(Node u, Node v, Edge prev = INVALID) const {
   102       if (prev != INVALID) return INVALID;
   103       int d = u._id ^ v._id;
   104       int k = 0;
   105       if (d == 0) return INVALID;
   106       for ( ; (d & 1) == 0; d >>= 1) ++k;
   107       if (d >> 1 != 0) return INVALID;
   108       return (k << (_dim-1)) | ((u._id >> (k+1)) << k) |
   109         (u._id & ((1 << k) - 1));
   110     }
   111 
   112     Arc findArc(Node u, Node v, Arc prev = INVALID) const {
   113       Edge edge = findEdge(u, v, prev);
   114       if (edge == INVALID) return INVALID;
   115       int k = edge._id >> (_dim-1);
   116       return ((u._id >> k) & 1) == 1 ? edge._id << 1 : (edge._id << 1) | 1;
   117     }
   118 
   119     class Node {
   120       friend class HypercubeGraphBase;
   121 
   122     protected:
   123       int _id;
   124       Node(int id) : _id(id) {}
   125     public:
   126       Node() {}
   127       Node (Invalid) : _id(-1) {}
   128       bool operator==(const Node node) const {return _id == node._id;}
   129       bool operator!=(const Node node) const {return _id != node._id;}
   130       bool operator<(const Node node) const {return _id < node._id;}
   131     };
   132 
   133     class Edge {
   134       friend class HypercubeGraphBase;
   135       friend class Arc;
   136 
   137     protected:
   138       int _id;
   139 
   140       Edge(int id) : _id(id) {}
   141 
   142     public:
   143       Edge() {}
   144       Edge (Invalid) : _id(-1) {}
   145       bool operator==(const Edge edge) const {return _id == edge._id;}
   146       bool operator!=(const Edge edge) const {return _id != edge._id;}
   147       bool operator<(const Edge edge) const {return _id < edge._id;}
   148     };
   149 
   150     class Arc {
   151       friend class HypercubeGraphBase;
   152 
   153     protected:
   154       int _id;
   155 
   156       Arc(int id) : _id(id) {}
   157 
   158     public:
   159       Arc() {}
   160       Arc (Invalid) : _id(-1) {}
   161       operator Edge() const { return _id != -1 ? Edge(_id >> 1) : INVALID; }
   162       bool operator==(const Arc arc) const {return _id == arc._id;}
   163       bool operator!=(const Arc arc) const {return _id != arc._id;}
   164       bool operator<(const Arc arc) const {return _id < arc._id;}
   165     };
   166 
   167     void first(Node& node) const {
   168       node._id = _node_num - 1;
   169     }
   170 
   171     static void next(Node& node) {
   172       --node._id;
   173     }
   174 
   175     void first(Edge& edge) const {
   176       edge._id = _edge_num - 1;
   177     }
   178 
   179     static void next(Edge& edge) {
   180       --edge._id;
   181     }
   182 
   183     void first(Arc& arc) const {
   184       arc._id = 2 * _edge_num - 1;
   185     }
   186 
   187     static void next(Arc& arc) {
   188       --arc._id;
   189     }
   190 
   191     void firstInc(Edge& edge, bool& dir, const Node& node) const {
   192       edge._id = node._id >> 1;
   193       dir = (node._id & 1) == 0;
   194     }
   195 
   196     void nextInc(Edge& edge, bool& dir) const {
   197       Node n = dir ? u(edge) : v(edge);
   198       int k = (edge._id >> (_dim-1)) + 1;
   199       if (k < _dim) {
   200         edge._id = (k << (_dim-1)) |
   201           ((n._id >> (k+1)) << k) | (n._id & ((1 << k) - 1));
   202         dir = ((n._id >> k) & 1) == 0;
   203       } else {
   204         edge._id = -1;
   205         dir = true;
   206       }
   207     }
   208 
   209     void firstOut(Arc& arc, const Node& node) const {
   210       arc._id = ((node._id >> 1) << 1) | (~node._id & 1);
   211     }
   212 
   213     void nextOut(Arc& arc) const {
   214       Node n = (arc._id & 1) == 1 ? u(arc) : v(arc);
   215       int k = (arc._id >> _dim) + 1;
   216       if (k < _dim) {
   217         arc._id = (k << (_dim-1)) |
   218           ((n._id >> (k+1)) << k) | (n._id & ((1 << k) - 1));
   219         arc._id = (arc._id << 1) | (~(n._id >> k) & 1);
   220       } else {
   221         arc._id = -1;
   222       }
   223     }
   224 
   225     void firstIn(Arc& arc, const Node& node) const {
   226       arc._id = ((node._id >> 1) << 1) | (node._id & 1);
   227     }
   228 
   229     void nextIn(Arc& arc) const {
   230       Node n = (arc._id & 1) == 1 ? v(arc) : u(arc);
   231       int k = (arc._id >> _dim) + 1;
   232       if (k < _dim) {
   233         arc._id = (k << (_dim-1)) |
   234           ((n._id >> (k+1)) << k) | (n._id & ((1 << k) - 1));
   235         arc._id = (arc._id << 1) | ((n._id >> k) & 1);
   236       } else {
   237         arc._id = -1;
   238       }
   239     }
   240 
   241     static bool direction(Arc arc) {
   242       return (arc._id & 1) == 1;
   243     }
   244 
   245     static Arc direct(Edge edge, bool dir) {
   246       return Arc((edge._id << 1) | (dir ? 1 : 0));
   247     }
   248 
   249     int dimension() const {
   250       return _dim;
   251     }
   252 
   253     bool projection(Node node, int n) const {
   254       return static_cast<bool>(node._id & (1 << n));
   255     }
   256 
   257     int dimension(Edge edge) const {
   258       return edge._id >> (_dim-1);
   259     }
   260 
   261     int dimension(Arc arc) const {
   262       return arc._id >> _dim;
   263     }
   264 
   265     static int index(Node node) {
   266       return node._id;
   267     }
   268 
   269     Node operator()(int ix) const {
   270       return Node(ix);
   271     }
   272 
   273   private:
   274     int _dim;
   275     int _node_num, _edge_num;
   276   };
   277 
   278 
   279   typedef GraphExtender<HypercubeGraphBase> ExtendedHypercubeGraphBase;
   280 
   281   /// \ingroup graphs
   282   ///
   283   /// \brief Hypercube graph class
   284   ///
   285   /// HypercubeGraph implements a special graph type. The nodes of the
   286   /// graph are indexed with integers having at most \c dim binary digits.
   287   /// Two nodes are connected in the graph if and only if their indices
   288   /// differ only on one position in the binary form.
   289   /// This class is completely static and it needs constant memory space.
   290   /// Thus you can neither add nor delete nodes or edges, however,
   291   /// the structure can be resized using resize().
   292   ///
   293   /// This type fully conforms to the \ref concepts::Graph "Graph concept".
   294   /// Most of its member functions and nested classes are documented
   295   /// only in the concept class.
   296   ///
   297   /// This class provides constant time counting for nodes, edges and arcs.
   298   ///
   299   /// \note The type of the indices is chosen to \c int for efficiency
   300   /// reasons. Thus the maximum dimension of this implementation is 26
   301   /// (assuming that the size of \c int is 32 bit).
   302   class HypercubeGraph : public ExtendedHypercubeGraphBase {
   303     typedef ExtendedHypercubeGraphBase Parent;
   304 
   305   public:
   306 
   307     /// \brief Constructs a hypercube graph with \c dim dimensions.
   308     ///
   309     /// Constructs a hypercube graph with \c dim dimensions.
   310     HypercubeGraph(int dim) { construct(dim); }
   311 
   312     /// \brief Resizes the graph
   313     ///
   314     /// This function resizes the graph. It fully destroys and
   315     /// rebuilds the structure, therefore the maps of the graph will be
   316     /// reallocated automatically and the previous values will be lost.
   317     void resize(int dim) {
   318       Parent::notifier(Arc()).clear();
   319       Parent::notifier(Edge()).clear();
   320       Parent::notifier(Node()).clear();
   321       construct(dim);
   322       Parent::notifier(Node()).build();
   323       Parent::notifier(Edge()).build();
   324       Parent::notifier(Arc()).build();
   325     }
   326 
   327     /// \brief The number of dimensions.
   328     ///
   329     /// Gives back the number of dimensions.
   330     int dimension() const {
   331       return Parent::dimension();
   332     }
   333 
   334     /// \brief Returns \c true if the n'th bit of the node is one.
   335     ///
   336     /// Returns \c true if the n'th bit of the node is one.
   337     bool projection(Node node, int n) const {
   338       return Parent::projection(node, n);
   339     }
   340 
   341     /// \brief The dimension id of an edge.
   342     ///
   343     /// Gives back the dimension id of the given edge.
   344     /// It is in the range <tt>[0..dim-1]</tt>.
   345     int dimension(Edge edge) const {
   346       return Parent::dimension(edge);
   347     }
   348 
   349     /// \brief The dimension id of an arc.
   350     ///
   351     /// Gives back the dimension id of the given arc.
   352     /// It is in the range <tt>[0..dim-1]</tt>.
   353     int dimension(Arc arc) const {
   354       return Parent::dimension(arc);
   355     }
   356 
   357     /// \brief The index of a node.
   358     ///
   359     /// Gives back the index of the given node.
   360     /// The lower bits of the integer describes the node.
   361     static int index(Node node) {
   362       return Parent::index(node);
   363     }
   364 
   365     /// \brief Gives back a node by its index.
   366     ///
   367     /// Gives back a node by its index.
   368     Node operator()(int ix) const {
   369       return Parent::operator()(ix);
   370     }
   371 
   372     /// \brief Number of nodes.
   373     int nodeNum() const { return Parent::nodeNum(); }
   374     /// \brief Number of edges.
   375     int edgeNum() const { return Parent::edgeNum(); }
   376     /// \brief Number of arcs.
   377     int arcNum() const { return Parent::arcNum(); }
   378 
   379     /// \brief Linear combination map.
   380     ///
   381     /// This map makes possible to give back a linear combination
   382     /// for each node. It works like the \c std::accumulate function,
   383     /// so it accumulates the \c bf binary function with the \c fv first
   384     /// value. The map accumulates only on that positions (dimensions)
   385     /// where the index of the node is one. The values that have to be
   386     /// accumulated should be given by the \c begin and \c end iterators
   387     /// and the length of this range should be equal to the dimension
   388     /// number of the graph.
   389     ///
   390     ///\code
   391     /// const int DIM = 3;
   392     /// HypercubeGraph graph(DIM);
   393     /// dim2::Point<double> base[DIM];
   394     /// for (int k = 0; k < DIM; ++k) {
   395     ///   base[k].x = rnd();
   396     ///   base[k].y = rnd();
   397     /// }
   398     /// HypercubeGraph::HyperMap<dim2::Point<double> >
   399     ///   pos(graph, base, base + DIM, dim2::Point<double>(0.0, 0.0));
   400     ///\endcode
   401     ///
   402     /// \see HypercubeGraph
   403     template <typename T, typename BF = std::plus<T> >
   404     class HyperMap {
   405     public:
   406 
   407       /// \brief The key type of the map
   408       typedef Node Key;
   409       /// \brief The value type of the map
   410       typedef T Value;
   411 
   412       /// \brief Constructor for HyperMap.
   413       ///
   414       /// Construct a HyperMap for the given graph. The values that have
   415       /// to be accumulated should be given by the \c begin and \c end
   416       /// iterators and the length of this range should be equal to the
   417       /// dimension number of the graph.
   418       ///
   419       /// This map accumulates the \c bf binary function with the \c fv
   420       /// first value on that positions (dimensions) where the index of
   421       /// the node is one.
   422       template <typename It>
   423       HyperMap(const Graph& graph, It begin, It end,
   424                T fv = 0, const BF& bf = BF())
   425         : _graph(graph), _values(begin, end), _first_value(fv), _bin_func(bf)
   426       {
   427         LEMON_ASSERT(_values.size() == graph.dimension(),
   428                      "Wrong size of range");
   429       }
   430 
   431       /// \brief The partial accumulated value.
   432       ///
   433       /// Gives back the partial accumulated value.
   434       Value operator[](const Key& k) const {
   435         Value val = _first_value;
   436         int id = _graph.index(k);
   437         int n = 0;
   438         while (id != 0) {
   439           if (id & 1) {
   440             val = _bin_func(val, _values[n]);
   441           }
   442           id >>= 1;
   443           ++n;
   444         }
   445         return val;
   446       }
   447 
   448     private:
   449       const Graph& _graph;
   450       std::vector<T> _values;
   451       T _first_value;
   452       BF _bin_func;
   453     };
   454 
   455   };
   456 
   457 }
   458 
   459 #endif