lemon/hypercube_graph.h
author Balazs Dezso <deba@inf.elte.hu>
Thu, 24 Jun 2010 09:27:53 +0200
changeset 891 bb70ad62c95f
parent 582 7a28e215f715
child 735 853fcddcf282
child 778 a143f19f465b
permissions -rw-r--r--
Fix critical bug in preflow (#372)

The wrong transition between the bound decrease and highest active
heuristics caused the bug. The last node chosen in bound decrease mode
is used in the first iteration in highest active mode.
kpeter@364
     1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
kpeter@364
     2
 *
kpeter@364
     3
 * This file is a part of LEMON, a generic C++ optimization library.
kpeter@364
     4
 *
alpar@440
     5
 * Copyright (C) 2003-2009
kpeter@364
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
kpeter@364
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
kpeter@364
     8
 *
kpeter@364
     9
 * Permission to use, modify and distribute this software is granted
kpeter@364
    10
 * provided that this copyright notice appears in all copies. For
kpeter@364
    11
 * precise terms see the accompanying LICENSE file.
kpeter@364
    12
 *
kpeter@364
    13
 * This software is provided "AS IS" with no warranty of any kind,
kpeter@364
    14
 * express or implied, and with no claim as to its suitability for any
kpeter@364
    15
 * purpose.
kpeter@364
    16
 *
kpeter@364
    17
 */
kpeter@364
    18
kpeter@364
    19
#ifndef HYPERCUBE_GRAPH_H
kpeter@364
    20
#define HYPERCUBE_GRAPH_H
kpeter@364
    21
kpeter@364
    22
#include <vector>
kpeter@364
    23
#include <lemon/core.h>
kpeter@365
    24
#include <lemon/assert.h>
kpeter@364
    25
#include <lemon/bits/graph_extender.h>
kpeter@364
    26
kpeter@364
    27
///\ingroup graphs
kpeter@364
    28
///\file
kpeter@365
    29
///\brief HypercubeGraph class.
kpeter@364
    30
kpeter@364
    31
namespace lemon {
kpeter@364
    32
kpeter@365
    33
  class HypercubeGraphBase {
kpeter@364
    34
kpeter@364
    35
  public:
kpeter@364
    36
kpeter@365
    37
    typedef HypercubeGraphBase Graph;
kpeter@364
    38
kpeter@364
    39
    class Node;
kpeter@365
    40
    class Edge;
kpeter@364
    41
    class Arc;
kpeter@364
    42
kpeter@364
    43
  public:
kpeter@364
    44
kpeter@365
    45
    HypercubeGraphBase() {}
kpeter@364
    46
kpeter@364
    47
  protected:
kpeter@364
    48
kpeter@364
    49
    void construct(int dim) {
kpeter@365
    50
      LEMON_ASSERT(dim >= 1, "The number of dimensions must be at least 1.");
kpeter@364
    51
      _dim = dim;
kpeter@365
    52
      _node_num = 1 << dim;
alpar@372
    53
      _edge_num = dim * (1 << (dim-1));
kpeter@364
    54
    }
kpeter@364
    55
kpeter@364
    56
  public:
kpeter@364
    57
kpeter@364
    58
    typedef True NodeNumTag;
kpeter@365
    59
    typedef True EdgeNumTag;
kpeter@364
    60
    typedef True ArcNumTag;
kpeter@364
    61
kpeter@365
    62
    int nodeNum() const { return _node_num; }
kpeter@365
    63
    int edgeNum() const { return _edge_num; }
kpeter@365
    64
    int arcNum() const { return 2 * _edge_num; }
kpeter@364
    65
kpeter@365
    66
    int maxNodeId() const { return _node_num - 1; }
kpeter@365
    67
    int maxEdgeId() const { return _edge_num - 1; }
kpeter@365
    68
    int maxArcId() const { return 2 * _edge_num - 1; }
kpeter@364
    69
kpeter@365
    70
    static Node nodeFromId(int id) { return Node(id); }
kpeter@365
    71
    static Edge edgeFromId(int id) { return Edge(id); }
kpeter@365
    72
    static Arc arcFromId(int id) { return Arc(id); }
kpeter@365
    73
kpeter@365
    74
    static int id(Node node) { return node._id; }
kpeter@365
    75
    static int id(Edge edge) { return edge._id; }
kpeter@365
    76
    static int id(Arc arc) { return arc._id; }
kpeter@365
    77
kpeter@365
    78
    Node u(Edge edge) const {
alpar@372
    79
      int base = edge._id & ((1 << (_dim-1)) - 1);
alpar@372
    80
      int k = edge._id >> (_dim-1);
alpar@372
    81
      return ((base >> k) << (k+1)) | (base & ((1 << k) - 1));
kpeter@364
    82
    }
kpeter@364
    83
kpeter@365
    84
    Node v(Edge edge) const {
alpar@372
    85
      int base = edge._id & ((1 << (_dim-1)) - 1);
alpar@372
    86
      int k = edge._id >> (_dim-1);
alpar@372
    87
      return ((base >> k) << (k+1)) | (base & ((1 << k) - 1)) | (1 << k);
kpeter@364
    88
    }
kpeter@364
    89
kpeter@365
    90
    Node source(Arc arc) const {
kpeter@365
    91
      return (arc._id & 1) == 1 ? u(arc) : v(arc);
kpeter@365
    92
    }
kpeter@364
    93
kpeter@365
    94
    Node target(Arc arc) const {
kpeter@365
    95
      return (arc._id & 1) == 1 ? v(arc) : u(arc);
kpeter@365
    96
    }
kpeter@364
    97
kpeter@365
    98
    typedef True FindEdgeTag;
kpeter@365
    99
    typedef True FindArcTag;
kpeter@365
   100
kpeter@365
   101
    Edge findEdge(Node u, Node v, Edge prev = INVALID) const {
kpeter@365
   102
      if (prev != INVALID) return INVALID;
kpeter@365
   103
      int d = u._id ^ v._id;
kpeter@365
   104
      int k = 0;
kpeter@365
   105
      if (d == 0) return INVALID;
kpeter@365
   106
      for ( ; (d & 1) == 0; d >>= 1) ++k;
kpeter@365
   107
      if (d >> 1 != 0) return INVALID;
alpar@372
   108
      return (k << (_dim-1)) | ((u._id >> (k+1)) << k) |
alpar@372
   109
        (u._id & ((1 << k) - 1));
kpeter@365
   110
    }
kpeter@365
   111
kpeter@365
   112
    Arc findArc(Node u, Node v, Arc prev = INVALID) const {
kpeter@365
   113
      Edge edge = findEdge(u, v, prev);
kpeter@365
   114
      if (edge == INVALID) return INVALID;
alpar@372
   115
      int k = edge._id >> (_dim-1);
kpeter@365
   116
      return ((u._id >> k) & 1) == 1 ? edge._id << 1 : (edge._id << 1) | 1;
kpeter@365
   117
    }
kpeter@364
   118
kpeter@364
   119
    class Node {
kpeter@365
   120
      friend class HypercubeGraphBase;
kpeter@365
   121
kpeter@364
   122
    protected:
kpeter@365
   123
      int _id;
kpeter@365
   124
      Node(int id) : _id(id) {}
kpeter@364
   125
    public:
kpeter@364
   126
      Node() {}
kpeter@365
   127
      Node (Invalid) : _id(-1) {}
kpeter@365
   128
      bool operator==(const Node node) const {return _id == node._id;}
kpeter@365
   129
      bool operator!=(const Node node) const {return _id != node._id;}
kpeter@365
   130
      bool operator<(const Node node) const {return _id < node._id;}
kpeter@365
   131
    };
kpeter@365
   132
kpeter@365
   133
    class Edge {
kpeter@365
   134
      friend class HypercubeGraphBase;
kpeter@365
   135
      friend class Arc;
kpeter@365
   136
kpeter@365
   137
    protected:
kpeter@365
   138
      int _id;
kpeter@365
   139
kpeter@365
   140
      Edge(int id) : _id(id) {}
kpeter@365
   141
kpeter@365
   142
    public:
kpeter@365
   143
      Edge() {}
kpeter@365
   144
      Edge (Invalid) : _id(-1) {}
kpeter@365
   145
      bool operator==(const Edge edge) const {return _id == edge._id;}
kpeter@365
   146
      bool operator!=(const Edge edge) const {return _id != edge._id;}
kpeter@365
   147
      bool operator<(const Edge edge) const {return _id < edge._id;}
kpeter@364
   148
    };
kpeter@364
   149
kpeter@364
   150
    class Arc {
kpeter@365
   151
      friend class HypercubeGraphBase;
kpeter@365
   152
kpeter@364
   153
    protected:
kpeter@365
   154
      int _id;
kpeter@365
   155
kpeter@365
   156
      Arc(int id) : _id(id) {}
kpeter@365
   157
kpeter@364
   158
    public:
kpeter@365
   159
      Arc() {}
kpeter@365
   160
      Arc (Invalid) : _id(-1) {}
kpeter@365
   161
      operator Edge() const { return _id != -1 ? Edge(_id >> 1) : INVALID; }
kpeter@365
   162
      bool operator==(const Arc arc) const {return _id == arc._id;}
kpeter@365
   163
      bool operator!=(const Arc arc) const {return _id != arc._id;}
kpeter@365
   164
      bool operator<(const Arc arc) const {return _id < arc._id;}
kpeter@364
   165
    };
kpeter@364
   166
kpeter@364
   167
    void first(Node& node) const {
kpeter@365
   168
      node._id = _node_num - 1;
kpeter@364
   169
    }
kpeter@364
   170
kpeter@364
   171
    static void next(Node& node) {
kpeter@365
   172
      --node._id;
kpeter@365
   173
    }
kpeter@365
   174
kpeter@365
   175
    void first(Edge& edge) const {
kpeter@365
   176
      edge._id = _edge_num - 1;
kpeter@365
   177
    }
kpeter@365
   178
kpeter@365
   179
    static void next(Edge& edge) {
kpeter@365
   180
      --edge._id;
kpeter@364
   181
    }
kpeter@364
   182
kpeter@364
   183
    void first(Arc& arc) const {
kpeter@365
   184
      arc._id = 2 * _edge_num - 1;
kpeter@364
   185
    }
kpeter@364
   186
kpeter@364
   187
    static void next(Arc& arc) {
kpeter@365
   188
      --arc._id;
kpeter@365
   189
    }
kpeter@365
   190
kpeter@365
   191
    void firstInc(Edge& edge, bool& dir, const Node& node) const {
kpeter@365
   192
      edge._id = node._id >> 1;
kpeter@365
   193
      dir = (node._id & 1) == 0;
kpeter@365
   194
    }
kpeter@365
   195
kpeter@365
   196
    void nextInc(Edge& edge, bool& dir) const {
kpeter@365
   197
      Node n = dir ? u(edge) : v(edge);
alpar@372
   198
      int k = (edge._id >> (_dim-1)) + 1;
kpeter@365
   199
      if (k < _dim) {
alpar@372
   200
        edge._id = (k << (_dim-1)) |
alpar@372
   201
          ((n._id >> (k+1)) << k) | (n._id & ((1 << k) - 1));
kpeter@365
   202
        dir = ((n._id >> k) & 1) == 0;
kpeter@365
   203
      } else {
kpeter@365
   204
        edge._id = -1;
kpeter@365
   205
        dir = true;
kpeter@365
   206
      }
kpeter@364
   207
    }
kpeter@364
   208
kpeter@364
   209
    void firstOut(Arc& arc, const Node& node) const {
kpeter@365
   210
      arc._id = ((node._id >> 1) << 1) | (~node._id & 1);
kpeter@364
   211
    }
kpeter@364
   212
kpeter@364
   213
    void nextOut(Arc& arc) const {
kpeter@365
   214
      Node n = (arc._id & 1) == 1 ? u(arc) : v(arc);
kpeter@365
   215
      int k = (arc._id >> _dim) + 1;
kpeter@365
   216
      if (k < _dim) {
alpar@372
   217
        arc._id = (k << (_dim-1)) |
alpar@372
   218
          ((n._id >> (k+1)) << k) | (n._id & ((1 << k) - 1));
kpeter@365
   219
        arc._id = (arc._id << 1) | (~(n._id >> k) & 1);
kpeter@365
   220
      } else {
kpeter@365
   221
        arc._id = -1;
kpeter@365
   222
      }
kpeter@364
   223
    }
kpeter@364
   224
kpeter@364
   225
    void firstIn(Arc& arc, const Node& node) const {
kpeter@365
   226
      arc._id = ((node._id >> 1) << 1) | (node._id & 1);
kpeter@364
   227
    }
kpeter@364
   228
kpeter@364
   229
    void nextIn(Arc& arc) const {
kpeter@365
   230
      Node n = (arc._id & 1) == 1 ? v(arc) : u(arc);
kpeter@365
   231
      int k = (arc._id >> _dim) + 1;
kpeter@365
   232
      if (k < _dim) {
alpar@372
   233
        arc._id = (k << (_dim-1)) |
alpar@372
   234
          ((n._id >> (k+1)) << k) | (n._id & ((1 << k) - 1));
kpeter@365
   235
        arc._id = (arc._id << 1) | ((n._id >> k) & 1);
kpeter@364
   236
      } else {
kpeter@365
   237
        arc._id = -1;
kpeter@364
   238
      }
kpeter@364
   239
    }
kpeter@364
   240
kpeter@365
   241
    static bool direction(Arc arc) {
kpeter@365
   242
      return (arc._id & 1) == 1;
kpeter@365
   243
    }
kpeter@365
   244
kpeter@365
   245
    static Arc direct(Edge edge, bool dir) {
kpeter@365
   246
      return Arc((edge._id << 1) | (dir ? 1 : 0));
kpeter@365
   247
    }
kpeter@365
   248
kpeter@364
   249
    int dimension() const {
kpeter@364
   250
      return _dim;
kpeter@364
   251
    }
kpeter@364
   252
kpeter@364
   253
    bool projection(Node node, int n) const {
kpeter@365
   254
      return static_cast<bool>(node._id & (1 << n));
kpeter@365
   255
    }
kpeter@365
   256
kpeter@365
   257
    int dimension(Edge edge) const {
alpar@372
   258
      return edge._id >> (_dim-1);
kpeter@364
   259
    }
kpeter@364
   260
kpeter@364
   261
    int dimension(Arc arc) const {
kpeter@365
   262
      return arc._id >> _dim;
kpeter@364
   263
    }
kpeter@364
   264
kpeter@364
   265
    int index(Node node) const {
kpeter@365
   266
      return node._id;
kpeter@364
   267
    }
kpeter@364
   268
kpeter@364
   269
    Node operator()(int ix) const {
kpeter@364
   270
      return Node(ix);
kpeter@364
   271
    }
kpeter@364
   272
kpeter@364
   273
  private:
kpeter@365
   274
    int _dim;
kpeter@365
   275
    int _node_num, _edge_num;
kpeter@364
   276
  };
kpeter@364
   277
kpeter@364
   278
kpeter@365
   279
  typedef GraphExtender<HypercubeGraphBase> ExtendedHypercubeGraphBase;
kpeter@364
   280
kpeter@365
   281
  /// \ingroup graphs
kpeter@364
   282
  ///
kpeter@365
   283
  /// \brief Hypercube graph class
kpeter@364
   284
  ///
kpeter@365
   285
  /// This class implements a special graph type. The nodes of the graph
kpeter@365
   286
  /// are indiced with integers with at most \c dim binary digits.
kpeter@365
   287
  /// Two nodes are connected in the graph if and only if their indices
kpeter@365
   288
  /// differ only on one position in the binary form.
kpeter@364
   289
  ///
kpeter@365
   290
  /// \note The type of the indices is chosen to \c int for efficiency
kpeter@365
   291
  /// reasons. Thus the maximum dimension of this implementation is 26
kpeter@365
   292
  /// (assuming that the size of \c int is 32 bit).
kpeter@364
   293
  ///
kpeter@559
   294
  /// This graph type fully conforms to the \ref concepts::Graph
kpeter@582
   295
  /// "Graph concept".
kpeter@365
   296
  class HypercubeGraph : public ExtendedHypercubeGraphBase {
kpeter@617
   297
    typedef ExtendedHypercubeGraphBase Parent;
kpeter@617
   298
kpeter@364
   299
  public:
kpeter@364
   300
kpeter@365
   301
    /// \brief Constructs a hypercube graph with \c dim dimensions.
kpeter@364
   302
    ///
kpeter@365
   303
    /// Constructs a hypercube graph with \c dim dimensions.
kpeter@365
   304
    HypercubeGraph(int dim) { construct(dim); }
kpeter@364
   305
kpeter@365
   306
    /// \brief The number of dimensions.
kpeter@364
   307
    ///
kpeter@365
   308
    /// Gives back the number of dimensions.
kpeter@364
   309
    int dimension() const {
kpeter@364
   310
      return Parent::dimension();
kpeter@364
   311
    }
kpeter@364
   312
kpeter@365
   313
    /// \brief Returns \c true if the n'th bit of the node is one.
kpeter@364
   314
    ///
kpeter@365
   315
    /// Returns \c true if the n'th bit of the node is one.
kpeter@364
   316
    bool projection(Node node, int n) const {
kpeter@364
   317
      return Parent::projection(node, n);
kpeter@364
   318
    }
kpeter@364
   319
kpeter@365
   320
    /// \brief The dimension id of an edge.
kpeter@364
   321
    ///
kpeter@365
   322
    /// Gives back the dimension id of the given edge.
kpeter@365
   323
    /// It is in the [0..dim-1] range.
kpeter@365
   324
    int dimension(Edge edge) const {
kpeter@365
   325
      return Parent::dimension(edge);
kpeter@365
   326
    }
kpeter@365
   327
kpeter@365
   328
    /// \brief The dimension id of an arc.
kpeter@365
   329
    ///
kpeter@365
   330
    /// Gives back the dimension id of the given arc.
kpeter@365
   331
    /// It is in the [0..dim-1] range.
kpeter@364
   332
    int dimension(Arc arc) const {
kpeter@364
   333
      return Parent::dimension(arc);
kpeter@364
   334
    }
kpeter@364
   335
kpeter@365
   336
    /// \brief The index of a node.
kpeter@364
   337
    ///
kpeter@365
   338
    /// Gives back the index of the given node.
kpeter@365
   339
    /// The lower bits of the integer describes the node.
kpeter@364
   340
    int index(Node node) const {
kpeter@364
   341
      return Parent::index(node);
kpeter@364
   342
    }
kpeter@364
   343
kpeter@365
   344
    /// \brief Gives back a node by its index.
kpeter@364
   345
    ///
kpeter@365
   346
    /// Gives back a node by its index.
kpeter@364
   347
    Node operator()(int ix) const {
kpeter@364
   348
      return Parent::operator()(ix);
kpeter@364
   349
    }
kpeter@364
   350
kpeter@364
   351
    /// \brief Number of nodes.
kpeter@364
   352
    int nodeNum() const { return Parent::nodeNum(); }
kpeter@365
   353
    /// \brief Number of edges.
kpeter@365
   354
    int edgeNum() const { return Parent::edgeNum(); }
kpeter@364
   355
    /// \brief Number of arcs.
kpeter@364
   356
    int arcNum() const { return Parent::arcNum(); }
kpeter@364
   357
kpeter@364
   358
    /// \brief Linear combination map.
kpeter@364
   359
    ///
kpeter@365
   360
    /// This map makes possible to give back a linear combination
kpeter@365
   361
    /// for each node. It works like the \c std::accumulate function,
kpeter@365
   362
    /// so it accumulates the \c bf binary function with the \c fv first
kpeter@365
   363
    /// value. The map accumulates only on that positions (dimensions)
kpeter@365
   364
    /// where the index of the node is one. The values that have to be
kpeter@365
   365
    /// accumulated should be given by the \c begin and \c end iterators
kpeter@365
   366
    /// and the length of this range should be equal to the dimension
kpeter@365
   367
    /// number of the graph.
kpeter@364
   368
    ///
kpeter@364
   369
    ///\code
kpeter@364
   370
    /// const int DIM = 3;
kpeter@365
   371
    /// HypercubeGraph graph(DIM);
kpeter@364
   372
    /// dim2::Point<double> base[DIM];
kpeter@364
   373
    /// for (int k = 0; k < DIM; ++k) {
kpeter@364
   374
    ///   base[k].x = rnd();
kpeter@364
   375
    ///   base[k].y = rnd();
kpeter@364
   376
    /// }
kpeter@365
   377
    /// HypercubeGraph::HyperMap<dim2::Point<double> >
kpeter@365
   378
    ///   pos(graph, base, base + DIM, dim2::Point<double>(0.0, 0.0));
kpeter@364
   379
    ///\endcode
kpeter@364
   380
    ///
kpeter@365
   381
    /// \see HypercubeGraph
kpeter@364
   382
    template <typename T, typename BF = std::plus<T> >
kpeter@364
   383
    class HyperMap {
kpeter@364
   384
    public:
kpeter@364
   385
kpeter@365
   386
      /// \brief The key type of the map
kpeter@364
   387
      typedef Node Key;
kpeter@365
   388
      /// \brief The value type of the map
kpeter@364
   389
      typedef T Value;
kpeter@364
   390
kpeter@364
   391
      /// \brief Constructor for HyperMap.
kpeter@364
   392
      ///
kpeter@365
   393
      /// Construct a HyperMap for the given graph. The values that have
kpeter@365
   394
      /// to be accumulated should be given by the \c begin and \c end
kpeter@365
   395
      /// iterators and the length of this range should be equal to the
kpeter@365
   396
      /// dimension number of the graph.
kpeter@364
   397
      ///
kpeter@365
   398
      /// This map accumulates the \c bf binary function with the \c fv
kpeter@365
   399
      /// first value on that positions (dimensions) where the index of
kpeter@365
   400
      /// the node is one.
kpeter@364
   401
      template <typename It>
kpeter@365
   402
      HyperMap(const Graph& graph, It begin, It end,
kpeter@365
   403
               T fv = 0, const BF& bf = BF())
kpeter@365
   404
        : _graph(graph), _values(begin, end), _first_value(fv), _bin_func(bf)
kpeter@364
   405
      {
kpeter@365
   406
        LEMON_ASSERT(_values.size() == graph.dimension(),
kpeter@365
   407
                     "Wrong size of range");
kpeter@364
   408
      }
kpeter@364
   409
kpeter@365
   410
      /// \brief The partial accumulated value.
kpeter@364
   411
      ///
kpeter@364
   412
      /// Gives back the partial accumulated value.
kpeter@365
   413
      Value operator[](const Key& k) const {
kpeter@364
   414
        Value val = _first_value;
kpeter@364
   415
        int id = _graph.index(k);
kpeter@364
   416
        int n = 0;
kpeter@364
   417
        while (id != 0) {
kpeter@364
   418
          if (id & 1) {
kpeter@364
   419
            val = _bin_func(val, _values[n]);
kpeter@364
   420
          }
kpeter@364
   421
          id >>= 1;
kpeter@364
   422
          ++n;
kpeter@364
   423
        }
kpeter@364
   424
        return val;
kpeter@364
   425
      }
kpeter@364
   426
kpeter@364
   427
    private:
kpeter@365
   428
      const Graph& _graph;
kpeter@364
   429
      std::vector<T> _values;
kpeter@364
   430
      T _first_value;
kpeter@364
   431
      BF _bin_func;
kpeter@364
   432
    };
kpeter@364
   433
kpeter@364
   434
  };
kpeter@364
   435
kpeter@364
   436
}
kpeter@364
   437
kpeter@364
   438
#endif