lemon/grid_graph.h
author Peter Kovacs <kpeter@inf.elte.hu>
Thu, 12 Nov 2009 23:26:13 +0100
changeset 806 fa6f37d7a25b
parent 735 853fcddcf282
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 GRID_GRAPH_H
    20 #define GRID_GRAPH_H
    21 
    22 #include <lemon/core.h>
    23 #include <lemon/bits/graph_extender.h>
    24 #include <lemon/dim2.h>
    25 #include <lemon/assert.h>
    26 
    27 ///\ingroup graphs
    28 ///\file
    29 ///\brief GridGraph class.
    30 
    31 namespace lemon {
    32 
    33   class GridGraphBase {
    34 
    35   public:
    36 
    37     typedef GridGraphBase Graph;
    38 
    39     class Node;
    40     class Edge;
    41     class Arc;
    42 
    43   public:
    44 
    45     GridGraphBase() {}
    46 
    47   protected:
    48 
    49     void construct(int width, int height) {
    50        _width = width; _height = height;
    51       _node_num = width * height;
    52       _edge_num = 2 * _node_num - width - height;
    53       _edge_limit = _node_num - _width;
    54     }
    55 
    56   public:
    57 
    58     Node operator()(int i, int j) const {
    59       LEMON_DEBUG(0 <= i && i < _width &&
    60                   0 <= j  && j < _height, "Index out of range");
    61       return Node(i + j * _width);
    62     }
    63 
    64     int col(Node n) const {
    65       return n._id % _width;
    66     }
    67 
    68     int row(Node n) const {
    69       return n._id / _width;
    70     }
    71 
    72     dim2::Point<int> pos(Node n) const {
    73       return dim2::Point<int>(col(n), row(n));
    74     }
    75 
    76     int width() const {
    77       return _width;
    78     }
    79 
    80     int height() const {
    81       return _height;
    82     }
    83 
    84     typedef True NodeNumTag;
    85     typedef True EdgeNumTag;
    86     typedef True ArcNumTag;
    87 
    88     int nodeNum() const { return _node_num; }
    89     int edgeNum() const { return _edge_num; }
    90     int arcNum() const { return 2 * _edge_num; }
    91 
    92     Node u(Edge edge) const {
    93       if (edge._id < _edge_limit) {
    94         return edge._id;
    95       } else {
    96         return (edge._id - _edge_limit) % (_width - 1) +
    97           (edge._id - _edge_limit) / (_width - 1) * _width;
    98       }
    99     }
   100 
   101     Node v(Edge edge) const {
   102       if (edge._id < _edge_limit) {
   103         return edge._id + _width;
   104       } else {
   105         return (edge._id - _edge_limit) % (_width - 1) +
   106           (edge._id - _edge_limit) / (_width - 1) * _width + 1;
   107       }
   108     }
   109 
   110     Node source(Arc arc) const {
   111       return (arc._id & 1) == 1 ? u(arc) : v(arc);
   112     }
   113 
   114     Node target(Arc arc) const {
   115       return (arc._id & 1) == 1 ? v(arc) : u(arc);
   116     }
   117 
   118     static int id(Node node) { return node._id; }
   119     static int id(Edge edge) { return edge._id; }
   120     static int id(Arc arc) { return arc._id; }
   121 
   122     int maxNodeId() const { return _node_num - 1; }
   123     int maxEdgeId() const { return _edge_num - 1; }
   124     int maxArcId() const { return 2 * _edge_num - 1; }
   125 
   126     static Node nodeFromId(int id) { return Node(id);}
   127     static Edge edgeFromId(int id) { return Edge(id);}
   128     static Arc arcFromId(int id) { return Arc(id);}
   129 
   130     typedef True FindEdgeTag;
   131     typedef True FindArcTag;
   132 
   133     Edge findEdge(Node u, Node v, Edge prev = INVALID) const {
   134       if (prev != INVALID) return INVALID;
   135       if (v._id > u._id) {
   136         if (v._id - u._id == _width)
   137           return Edge(u._id);
   138         if (v._id - u._id == 1 && u._id % _width < _width - 1) {
   139           return Edge(u._id / _width * (_width - 1) +
   140                       u._id % _width + _edge_limit);
   141         }
   142       } else {
   143         if (u._id - v._id == _width)
   144           return Edge(v._id);
   145         if (u._id - v._id == 1 && v._id % _width < _width - 1) {
   146           return Edge(v._id / _width * (_width - 1) +
   147                       v._id % _width + _edge_limit);
   148         }
   149       }
   150       return INVALID;
   151     }
   152 
   153     Arc findArc(Node u, Node v, Arc prev = INVALID) const {
   154       if (prev != INVALID) return INVALID;
   155       if (v._id > u._id) {
   156         if (v._id - u._id == _width)
   157           return Arc((u._id << 1) | 1);
   158         if (v._id - u._id == 1 && u._id % _width < _width - 1) {
   159           return Arc(((u._id / _width * (_width - 1) +
   160                        u._id % _width + _edge_limit) << 1) | 1);
   161         }
   162       } else {
   163         if (u._id - v._id == _width)
   164           return Arc(v._id << 1);
   165         if (u._id - v._id == 1 && v._id % _width < _width - 1) {
   166           return Arc((v._id / _width * (_width - 1) +
   167                        v._id % _width + _edge_limit) << 1);
   168         }
   169       }
   170       return INVALID;
   171     }
   172 
   173     class Node {
   174       friend class GridGraphBase;
   175 
   176     protected:
   177       int _id;
   178       Node(int id) : _id(id) {}
   179     public:
   180       Node() {}
   181       Node (Invalid) : _id(-1) {}
   182       bool operator==(const Node node) const {return _id == node._id;}
   183       bool operator!=(const Node node) const {return _id != node._id;}
   184       bool operator<(const Node node) const {return _id < node._id;}
   185     };
   186 
   187     class Edge {
   188       friend class GridGraphBase;
   189       friend class Arc;
   190 
   191     protected:
   192       int _id;
   193 
   194       Edge(int id) : _id(id) {}
   195 
   196     public:
   197       Edge() {}
   198       Edge (Invalid) : _id(-1) {}
   199       bool operator==(const Edge edge) const {return _id == edge._id;}
   200       bool operator!=(const Edge edge) const {return _id != edge._id;}
   201       bool operator<(const Edge edge) const {return _id < edge._id;}
   202     };
   203 
   204     class Arc {
   205       friend class GridGraphBase;
   206 
   207     protected:
   208       int _id;
   209 
   210       Arc(int id) : _id(id) {}
   211 
   212     public:
   213       Arc() {}
   214       Arc (Invalid) : _id(-1) {}
   215       operator Edge() const { return _id != -1 ? Edge(_id >> 1) : INVALID; }
   216       bool operator==(const Arc arc) const {return _id == arc._id;}
   217       bool operator!=(const Arc arc) const {return _id != arc._id;}
   218       bool operator<(const Arc arc) const {return _id < arc._id;}
   219     };
   220 
   221     static bool direction(Arc arc) {
   222       return (arc._id & 1) == 1;
   223     }
   224 
   225     static Arc direct(Edge edge, bool dir) {
   226       return Arc((edge._id << 1) | (dir ? 1 : 0));
   227     }
   228 
   229     void first(Node& node) const {
   230       node._id = _node_num - 1;
   231     }
   232 
   233     static void next(Node& node) {
   234       --node._id;
   235     }
   236 
   237     void first(Edge& edge) const {
   238       edge._id = _edge_num - 1;
   239     }
   240 
   241     static void next(Edge& edge) {
   242       --edge._id;
   243     }
   244 
   245     void first(Arc& arc) const {
   246       arc._id = 2 * _edge_num - 1;
   247     }
   248 
   249     static void next(Arc& arc) {
   250       --arc._id;
   251     }
   252 
   253     void firstOut(Arc& arc, const Node& node) const {
   254       if (node._id % _width < _width - 1) {
   255         arc._id = (_edge_limit + node._id % _width +
   256                    (node._id / _width) * (_width - 1)) << 1 | 1;
   257         return;
   258       }
   259       if (node._id < _node_num - _width) {
   260         arc._id = node._id << 1 | 1;
   261         return;
   262       }
   263       if (node._id % _width > 0) {
   264         arc._id = (_edge_limit + node._id % _width +
   265                    (node._id / _width) * (_width - 1) - 1) << 1;
   266         return;
   267       }
   268       if (node._id >= _width) {
   269         arc._id = (node._id - _width) << 1;
   270         return;
   271       }
   272       arc._id = -1;
   273     }
   274 
   275     void nextOut(Arc& arc) const {
   276       int nid = arc._id >> 1;
   277       if ((arc._id & 1) == 1) {
   278         if (nid >= _edge_limit) {
   279           nid = (nid - _edge_limit) % (_width - 1) +
   280             (nid - _edge_limit) / (_width - 1) * _width;
   281           if (nid < _node_num - _width) {
   282             arc._id = nid << 1 | 1;
   283             return;
   284           }
   285         }
   286         if (nid % _width > 0) {
   287           arc._id = (_edge_limit + nid % _width +
   288                      (nid / _width) * (_width - 1) - 1) << 1;
   289           return;
   290         }
   291         if (nid >= _width) {
   292           arc._id = (nid - _width) << 1;
   293           return;
   294         }
   295       } else {
   296         if (nid >= _edge_limit) {
   297           nid = (nid - _edge_limit) % (_width - 1) +
   298             (nid - _edge_limit) / (_width - 1) * _width + 1;
   299           if (nid >= _width) {
   300             arc._id = (nid - _width) << 1;
   301             return;
   302           }
   303         }
   304       }
   305       arc._id = -1;
   306     }
   307 
   308     void firstIn(Arc& arc, const Node& node) const {
   309       if (node._id % _width < _width - 1) {
   310         arc._id = (_edge_limit + node._id % _width +
   311                    (node._id / _width) * (_width - 1)) << 1;
   312         return;
   313       }
   314       if (node._id < _node_num - _width) {
   315         arc._id = node._id << 1;
   316         return;
   317       }
   318       if (node._id % _width > 0) {
   319         arc._id = (_edge_limit + node._id % _width +
   320                    (node._id / _width) * (_width - 1) - 1) << 1 | 1;
   321         return;
   322       }
   323       if (node._id >= _width) {
   324         arc._id = (node._id - _width) << 1 | 1;
   325         return;
   326       }
   327       arc._id = -1;
   328     }
   329 
   330     void nextIn(Arc& arc) const {
   331       int nid = arc._id >> 1;
   332       if ((arc._id & 1) == 0) {
   333         if (nid >= _edge_limit) {
   334           nid = (nid - _edge_limit) % (_width - 1) +
   335             (nid - _edge_limit) / (_width - 1) * _width;
   336           if (nid < _node_num - _width) {
   337             arc._id = nid << 1;
   338             return;
   339           }
   340         }
   341         if (nid % _width > 0) {
   342           arc._id = (_edge_limit + nid % _width +
   343                      (nid / _width) * (_width - 1) - 1) << 1 | 1;
   344           return;
   345         }
   346         if (nid >= _width) {
   347           arc._id = (nid - _width) << 1 | 1;
   348           return;
   349         }
   350       } else {
   351         if (nid >= _edge_limit) {
   352           nid = (nid - _edge_limit) % (_width - 1) +
   353             (nid - _edge_limit) / (_width - 1) * _width + 1;
   354           if (nid >= _width) {
   355             arc._id = (nid - _width) << 1 | 1;
   356             return;
   357           }
   358         }
   359       }
   360       arc._id = -1;
   361     }
   362 
   363     void firstInc(Edge& edge, bool& dir, const Node& node) const {
   364       if (node._id % _width < _width - 1) {
   365         edge._id = _edge_limit + node._id % _width +
   366           (node._id / _width) * (_width - 1);
   367         dir = true;
   368         return;
   369       }
   370       if (node._id < _node_num - _width) {
   371         edge._id = node._id;
   372         dir = true;
   373         return;
   374       }
   375       if (node._id % _width > 0) {
   376         edge._id = _edge_limit + node._id % _width +
   377           (node._id / _width) * (_width - 1) - 1;
   378         dir = false;
   379         return;
   380       }
   381       if (node._id >= _width) {
   382         edge._id = node._id - _width;
   383         dir = false;
   384         return;
   385       }
   386       edge._id = -1;
   387       dir = true;
   388     }
   389 
   390     void nextInc(Edge& edge, bool& dir) const {
   391       int nid = edge._id;
   392       if (dir) {
   393         if (nid >= _edge_limit) {
   394           nid = (nid - _edge_limit) % (_width - 1) +
   395             (nid - _edge_limit) / (_width - 1) * _width;
   396           if (nid < _node_num - _width) {
   397             edge._id = nid;
   398             return;
   399           }
   400         }
   401         if (nid % _width > 0) {
   402           edge._id = _edge_limit + nid % _width +
   403             (nid / _width) * (_width - 1) - 1;
   404           dir = false;
   405           return;
   406         }
   407         if (nid >= _width) {
   408           edge._id = nid - _width;
   409           dir = false;
   410           return;
   411         }
   412       } else {
   413         if (nid >= _edge_limit) {
   414           nid = (nid - _edge_limit) % (_width - 1) +
   415             (nid - _edge_limit) / (_width - 1) * _width + 1;
   416           if (nid >= _width) {
   417             edge._id = nid - _width;
   418             return;
   419           }
   420         }
   421       }
   422       edge._id = -1;
   423       dir = true;
   424     }
   425 
   426     Arc right(Node n) const {
   427       if (n._id % _width < _width - 1) {
   428         return Arc(((_edge_limit + n._id % _width +
   429                     (n._id / _width) * (_width - 1)) << 1) | 1);
   430       } else {
   431         return INVALID;
   432       }
   433     }
   434 
   435     Arc left(Node n) const {
   436       if (n._id % _width > 0) {
   437         return Arc((_edge_limit + n._id % _width +
   438                      (n._id / _width) * (_width - 1) - 1) << 1);
   439       } else {
   440         return INVALID;
   441       }
   442     }
   443 
   444     Arc up(Node n) const {
   445       if (n._id < _edge_limit) {
   446         return Arc((n._id << 1) | 1);
   447       } else {
   448         return INVALID;
   449       }
   450     }
   451 
   452     Arc down(Node n) const {
   453       if (n._id >= _width) {
   454         return Arc((n._id - _width) << 1);
   455       } else {
   456         return INVALID;
   457       }
   458     }
   459 
   460   private:
   461     int _width, _height;
   462     int _node_num, _edge_num;
   463     int _edge_limit;
   464   };
   465 
   466 
   467   typedef GraphExtender<GridGraphBase> ExtendedGridGraphBase;
   468 
   469   /// \ingroup graphs
   470   ///
   471   /// \brief Grid graph class
   472   ///
   473   /// GridGraph implements a special graph type. The nodes of the
   474   /// graph can be indexed by two integer values \c (i,j) where \c i is
   475   /// in the range <tt>[0..width()-1]</tt> and j is in the range
   476   /// <tt>[0..height()-1]</tt>. Two nodes are connected in the graph if
   477   /// the indices differ exactly on one position and the difference is
   478   /// also exactly one. The nodes of the graph can be obtained by position
   479   /// using the \c operator()() function and the indices of the nodes can
   480   /// be obtained using \c pos(), \c col() and \c row() members. The outgoing
   481   /// arcs can be retrieved with the \c right(), \c up(), \c left()
   482   /// and \c down() functions, where the bottom-left corner is the
   483   /// origin.
   484   ///
   485   /// This class is completely static and it needs constant memory space.
   486   /// Thus you can neither add nor delete nodes or edges, however
   487   /// the structure can be resized using resize().
   488   ///
   489   /// \image html grid_graph.png
   490   /// \image latex grid_graph.eps "Grid graph" width=\textwidth
   491   ///
   492   /// A short example about the basic usage:
   493   ///\code
   494   /// GridGraph graph(rows, cols);
   495   /// GridGraph::NodeMap<int> val(graph);
   496   /// for (int i = 0; i < graph.width(); ++i) {
   497   ///   for (int j = 0; j < graph.height(); ++j) {
   498   ///     val[graph(i, j)] = i + j;
   499   ///   }
   500   /// }
   501   ///\endcode
   502   ///
   503   /// This type fully conforms to the \ref concepts::Graph "Graph concept".
   504   /// Most of its member functions and nested classes are documented
   505   /// only in the concept class.
   506   ///
   507   /// This class provides constant time counting for nodes, edges and arcs.
   508   class GridGraph : public ExtendedGridGraphBase {
   509     typedef ExtendedGridGraphBase Parent;
   510 
   511   public:
   512 
   513     /// \brief Map to get the indices of the nodes as \ref dim2::Point
   514     /// "dim2::Point<int>".
   515     ///
   516     /// Map to get the indices of the nodes as \ref dim2::Point
   517     /// "dim2::Point<int>".
   518     class IndexMap {
   519     public:
   520       /// \brief The key type of the map
   521       typedef GridGraph::Node Key;
   522       /// \brief The value type of the map
   523       typedef dim2::Point<int> Value;
   524 
   525       /// \brief Constructor
   526       IndexMap(const GridGraph& graph) : _graph(graph) {}
   527 
   528       /// \brief The subscript operator
   529       Value operator[](Key key) const {
   530         return _graph.pos(key);
   531       }
   532 
   533     private:
   534       const GridGraph& _graph;
   535     };
   536 
   537     /// \brief Map to get the column of the nodes.
   538     ///
   539     /// Map to get the column of the nodes.
   540     class ColMap {
   541     public:
   542       /// \brief The key type of the map
   543       typedef GridGraph::Node Key;
   544       /// \brief The value type of the map
   545       typedef int Value;
   546 
   547       /// \brief Constructor
   548       ColMap(const GridGraph& graph) : _graph(graph) {}
   549 
   550       /// \brief The subscript operator
   551       Value operator[](Key key) const {
   552         return _graph.col(key);
   553       }
   554 
   555     private:
   556       const GridGraph& _graph;
   557     };
   558 
   559     /// \brief Map to get the row of the nodes.
   560     ///
   561     /// Map to get the row of the nodes.
   562     class RowMap {
   563     public:
   564       /// \brief The key type of the map
   565       typedef GridGraph::Node Key;
   566       /// \brief The value type of the map
   567       typedef int Value;
   568 
   569       /// \brief Constructor
   570       RowMap(const GridGraph& graph) : _graph(graph) {}
   571 
   572       /// \brief The subscript operator
   573       Value operator[](Key key) const {
   574         return _graph.row(key);
   575       }
   576 
   577     private:
   578       const GridGraph& _graph;
   579     };
   580 
   581     /// \brief Constructor
   582     ///
   583     /// Construct a grid graph with the given size.
   584     GridGraph(int width, int height) { construct(width, height); }
   585 
   586     /// \brief Resizes the graph
   587     ///
   588     /// This function resizes the graph. It fully destroys and
   589     /// rebuilds the structure, therefore the maps of the graph will be
   590     /// reallocated automatically and the previous values will be lost.
   591     void resize(int width, int height) {
   592       Parent::notifier(Arc()).clear();
   593       Parent::notifier(Edge()).clear();
   594       Parent::notifier(Node()).clear();
   595       construct(width, height);
   596       Parent::notifier(Node()).build();
   597       Parent::notifier(Edge()).build();
   598       Parent::notifier(Arc()).build();
   599     }
   600 
   601     /// \brief The node on the given position.
   602     ///
   603     /// Gives back the node on the given position.
   604     Node operator()(int i, int j) const {
   605       return Parent::operator()(i, j);
   606     }
   607 
   608     /// \brief The column index of the node.
   609     ///
   610     /// Gives back the column index of the node.
   611     int col(Node n) const {
   612       return Parent::col(n);
   613     }
   614 
   615     /// \brief The row index of the node.
   616     ///
   617     /// Gives back the row index of the node.
   618     int row(Node n) const {
   619       return Parent::row(n);
   620     }
   621 
   622     /// \brief The position of the node.
   623     ///
   624     /// Gives back the position of the node, ie. the <tt>(col,row)</tt> pair.
   625     dim2::Point<int> pos(Node n) const {
   626       return Parent::pos(n);
   627     }
   628 
   629     /// \brief The number of the columns.
   630     ///
   631     /// Gives back the number of the columns.
   632     int width() const {
   633       return Parent::width();
   634     }
   635 
   636     /// \brief The number of the rows.
   637     ///
   638     /// Gives back the number of the rows.
   639     int height() const {
   640       return Parent::height();
   641     }
   642 
   643     /// \brief The arc goes right from the node.
   644     ///
   645     /// Gives back the arc goes right from the node. If there is not
   646     /// outgoing arc then it gives back INVALID.
   647     Arc right(Node n) const {
   648       return Parent::right(n);
   649     }
   650 
   651     /// \brief The arc goes left from the node.
   652     ///
   653     /// Gives back the arc goes left from the node. If there is not
   654     /// outgoing arc then it gives back INVALID.
   655     Arc left(Node n) const {
   656       return Parent::left(n);
   657     }
   658 
   659     /// \brief The arc goes up from the node.
   660     ///
   661     /// Gives back the arc goes up from the node. If there is not
   662     /// outgoing arc then it gives back INVALID.
   663     Arc up(Node n) const {
   664       return Parent::up(n);
   665     }
   666 
   667     /// \brief The arc goes down from the node.
   668     ///
   669     /// Gives back the arc goes down from the node. If there is not
   670     /// outgoing arc then it gives back INVALID.
   671     Arc down(Node n) const {
   672       return Parent::down(n);
   673     }
   674 
   675     /// \brief Index map of the grid graph
   676     ///
   677     /// Just returns an IndexMap for the grid graph.
   678     IndexMap indexMap() const {
   679       return IndexMap(*this);
   680     }
   681 
   682     /// \brief Row map of the grid graph
   683     ///
   684     /// Just returns a RowMap for the grid graph.
   685     RowMap rowMap() const {
   686       return RowMap(*this);
   687     }
   688 
   689     /// \brief Column map of the grid graph
   690     ///
   691     /// Just returns a ColMap for the grid graph.
   692     ColMap colMap() const {
   693       return ColMap(*this);
   694     }
   695 
   696   };
   697 
   698 }
   699 #endif