[Lemon-commits] Peter Kovacs: Entirely rework CapacityScaling (#...

Lemon HG hg at lemon.cs.elte.hu
Mon Dec 14 06:17:44 CET 2009


details:   http://lemon.cs.elte.hu/hg/lemon/rev/fa6f37d7a25b
changeset: 872:fa6f37d7a25b
user:      Peter Kovacs <kpeter [at] inf.elte.hu>
date:      Thu Nov 12 23:26:13 2009 +0100
description:
	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.

diffstat:

 lemon/capacity_scaling.h |  1103 +++++++++++++++++++++++++++--------------------
 1 files changed, 635 insertions(+), 468 deletions(-)

diffs (truncated from 1323 to 300 lines):

diff --git a/lemon/capacity_scaling.h b/lemon/capacity_scaling.h
--- a/lemon/capacity_scaling.h
+++ b/lemon/capacity_scaling.h
@@ -19,396 +19,417 @@
 #ifndef LEMON_CAPACITY_SCALING_H
 #define LEMON_CAPACITY_SCALING_H
 
-/// \ingroup min_cost_flow
+/// \ingroup min_cost_flow_algs
 ///
 /// \file
-/// \brief Capacity scaling algorithm for finding a minimum cost flow.
+/// \brief Capacity Scaling algorithm for finding a minimum cost flow.
 
 #include <vector>
+#include <limits>
+#include <lemon/core.h>
 #include <lemon/bin_heap.h>
 
 namespace lemon {
 
-  /// \addtogroup min_cost_flow
+  /// \addtogroup min_cost_flow_algs
   /// @{
 
-  /// \brief Implementation of the capacity scaling algorithm for
-  /// finding a minimum cost flow.
+  /// \brief Implementation of the Capacity Scaling algorithm for
+  /// finding a \ref min_cost_flow "minimum cost flow".
   ///
   /// \ref CapacityScaling implements the capacity scaling version
-  /// of the successive shortest path algorithm for finding a minimum
-  /// cost flow.
+  /// of the successive shortest path algorithm for finding a
+  /// \ref min_cost_flow "minimum cost flow". It is an efficient dual
+  /// solution method.
   ///
-  /// \tparam Digraph The digraph type the algorithm runs on.
-  /// \tparam LowerMap The type of the lower bound map.
-  /// \tparam CapacityMap The type of the capacity (upper bound) map.
-  /// \tparam CostMap The type of the cost (length) map.
-  /// \tparam SupplyMap The type of the supply map.
+  /// Most of the parameters of the problem (except for the digraph)
+  /// can be given using separate functions, and the algorithm can be
+  /// executed using the \ref run() function. If some parameters are not
+  /// specified, then default values will be used.
   ///
-  /// \warning
-  /// - Arc capacities and costs should be \e non-negative \e integers.
-  /// - Supply values should be \e signed \e integers.
-  /// - The value types of the maps should be convertible to each other.
-  /// - \c CostMap::Value must be signed type.
+  /// \tparam GR The digraph type the algorithm runs on.
+  /// \tparam V The value type used for flow amounts, capacity bounds
+  /// and supply values in the algorithm. By default it is \c int.
+  /// \tparam C The value type used for costs and potentials in the
+  /// algorithm. By default it is the same as \c V.
   ///
-  /// \author Peter Kovacs
-  template < typename Digraph,
-             typename LowerMap = typename Digraph::template ArcMap<int>,
-             typename CapacityMap = typename Digraph::template ArcMap<int>,
-             typename CostMap = typename Digraph::template ArcMap<int>,
-             typename SupplyMap = typename Digraph::template NodeMap<int> >
+  /// \warning Both value types must be signed and all input data must
+  /// be integer.
+  /// \warning This algorithm does not support negative costs for such
+  /// arcs that have infinite upper bound.
+  template <typename GR, typename V = int, typename C = V>
   class CapacityScaling
   {
-    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
+  public:
 
-    typedef typename CapacityMap::Value Capacity;
-    typedef typename CostMap::Value Cost;
-    typedef typename SupplyMap::Value Supply;
-    typedef typename Digraph::template ArcMap<Capacity> CapacityArcMap;
-    typedef typename Digraph::template NodeMap<Supply> SupplyNodeMap;
-    typedef typename Digraph::template NodeMap<Arc> PredMap;
+    /// The type of the flow amounts, capacity bounds and supply values
+    typedef V Value;
+    /// The type of the arc costs
+    typedef C Cost;
 
   public:
 
-    /// The type of the flow map.
-    typedef typename Digraph::template ArcMap<Capacity> FlowMap;
-    /// The type of the potential map.
-    typedef typename Digraph::template NodeMap<Cost> PotentialMap;
+    /// \brief Problem type constants for the \c run() function.
+    ///
+    /// Enum type containing the problem type constants that can be
+    /// returned by the \ref run() function of the algorithm.
+    enum ProblemType {
+      /// The problem has no feasible solution (flow).
+      INFEASIBLE,
+      /// The problem has optimal solution (i.e. it is feasible and
+      /// bounded), and the algorithm has found optimal flow and node
+      /// potentials (primal and dual solutions).
+      OPTIMAL,
+      /// The digraph contains an arc of negative cost and infinite
+      /// upper bound. It means that the objective function is unbounded
+      /// on that arc, however note that it could actually be bounded
+      /// over the feasible flows, but this algroithm cannot handle
+      /// these cases.
+      UNBOUNDED
+    };
+  
+  private:
+
+    TEMPLATE_DIGRAPH_TYPEDEFS(GR);
+
+    typedef std::vector<Arc> ArcVector;
+    typedef std::vector<Node> NodeVector;
+    typedef std::vector<int> IntVector;
+    typedef std::vector<bool> BoolVector;
+    typedef std::vector<Value> ValueVector;
+    typedef std::vector<Cost> CostVector;
 
   private:
 
-    /// \brief Special implementation of the \ref Dijkstra algorithm
-    /// for finding shortest paths in the residual network.
+    // Data related to the underlying digraph
+    const GR &_graph;
+    int _node_num;
+    int _arc_num;
+    int _res_arc_num;
+    int _root;
+
+    // Parameters of the problem
+    bool _have_lower;
+    Value _sum_supply;
+
+    // Data structures for storing the digraph
+    IntNodeMap _node_id;
+    IntArcMap _arc_idf;
+    IntArcMap _arc_idb;
+    IntVector _first_out;
+    BoolVector _forward;
+    IntVector _source;
+    IntVector _target;
+    IntVector _reverse;
+
+    // Node and arc data
+    ValueVector _lower;
+    ValueVector _upper;
+    CostVector _cost;
+    ValueVector _supply;
+
+    ValueVector _res_cap;
+    CostVector _pi;
+    ValueVector _excess;
+    IntVector _excess_nodes;
+    IntVector _deficit_nodes;
+
+    Value _delta;
+    int _phase_num;
+    IntVector _pred;
+
+  public:
+  
+    /// \brief Constant for infinite upper bounds (capacities).
     ///
-    /// \ref ResidualDijkstra is a special implementation of the
-    /// \ref Dijkstra algorithm for finding shortest paths in the
-    /// residual network of the digraph with respect to the reduced arc
-    /// costs and modifying the node potentials according to the
-    /// distance of the nodes.
+    /// Constant for infinite upper bounds (capacities).
+    /// It is \c std::numeric_limits<Value>::infinity() if available,
+    /// \c std::numeric_limits<Value>::max() otherwise.
+    const Value INF;
+
+  private:
+
+    // Special implementation of the Dijkstra algorithm for finding
+    // shortest paths in the residual network of the digraph with
+    // respect to the reduced arc costs and modifying the node
+    // potentials according to the found distance labels.
     class ResidualDijkstra
     {
-      typedef typename Digraph::template NodeMap<int> HeapCrossRef;
+      typedef RangeMap<int> HeapCrossRef;
       typedef BinHeap<Cost, HeapCrossRef> Heap;
 
     private:
 
-      // The digraph the algorithm runs on
-      const Digraph &_graph;
-
-      // The main maps
-      const FlowMap &_flow;
-      const CapacityArcMap &_res_cap;
-      const CostMap &_cost;
-      const SupplyNodeMap &_excess;
-      PotentialMap &_potential;
-
-      // The distance map
-      PotentialMap _dist;
-      // The pred arc map
-      PredMap &_pred;
-      // The processed (i.e. permanently labeled) nodes
-      std::vector<Node> _proc_nodes;
-
+      int _node_num;
+      const IntVector &_first_out;
+      const IntVector &_target;
+      const CostVector &_cost;
+      const ValueVector &_res_cap;
+      const ValueVector &_excess;
+      CostVector &_pi;
+      IntVector &_pred;
+      
+      IntVector _proc_nodes;
+      CostVector _dist;
+      
     public:
 
-      /// Constructor.
-      ResidualDijkstra( const Digraph &digraph,
-                        const FlowMap &flow,
-                        const CapacityArcMap &res_cap,
-                        const CostMap &cost,
-                        const SupplyMap &excess,
-                        PotentialMap &potential,
-                        PredMap &pred ) :
-        _graph(digraph), _flow(flow), _res_cap(res_cap), _cost(cost),
-        _excess(excess), _potential(potential), _dist(digraph),
-        _pred(pred)
+      ResidualDijkstra(CapacityScaling& cs) :
+        _node_num(cs._node_num), _first_out(cs._first_out), 
+        _target(cs._target), _cost(cs._cost), _res_cap(cs._res_cap),
+        _excess(cs._excess), _pi(cs._pi), _pred(cs._pred),
+        _dist(cs._node_num)
       {}
 
-      /// Run the algorithm from the given source node.
-      Node run(Node s, Capacity delta = 1) {
-        HeapCrossRef heap_cross_ref(_graph, Heap::PRE_HEAP);
+      int run(int s, Value delta = 1) {
+        HeapCrossRef heap_cross_ref(_node_num, Heap::PRE_HEAP);
         Heap heap(heap_cross_ref);
         heap.push(s, 0);
-        _pred[s] = INVALID;
+        _pred[s] = -1;
         _proc_nodes.clear();
 
-        // Processing nodes
+        // Process nodes
         while (!heap.empty() && _excess[heap.top()] > -delta) {
-          Node u = heap.top(), v;
-          Cost d = heap.prio() + _potential[u], nd;
+          int u = heap.top(), v;
+          Cost d = heap.prio() + _pi[u], dn;
           _dist[u] = heap.prio();
+          _proc_nodes.push_back(u);
           heap.pop();
-          _proc_nodes.push_back(u);
 
-          // Traversing outgoing arcs
-          for (OutArcIt e(_graph, u); e != INVALID; ++e) {
-            if (_res_cap[e] >= delta) {
-              v = _graph.target(e);
-              switch(heap.state(v)) {
+          // Traverse outgoing residual arcs
+          for (int a = _first_out[u]; a != _first_out[u+1]; ++a) {
+            if (_res_cap[a] < delta) continue;
+            v = _target[a];
+            switch (heap.state(v)) {
               case Heap::PRE_HEAP:
-                heap.push(v, d + _cost[e] - _potential[v]);
-                _pred[v] = e;
+                heap.push(v, d + _cost[a] - _pi[v]);
+                _pred[v] = a;
                 break;
               case Heap::IN_HEAP:
-                nd = d + _cost[e] - _potential[v];
-                if (nd < heap[v]) {
-                  heap.decrease(v, nd);
-                  _pred[v] = e;
+                dn = d + _cost[a] - _pi[v];
+                if (dn < heap[v]) {
+                  heap.decrease(v, dn);
+                  _pred[v] = a;
                 }
                 break;
               case Heap::POST_HEAP:
                 break;
-              }
-            }
-          }
-
-          // Traversing incoming arcs
-          for (InArcIt e(_graph, u); e != INVALID; ++e) {
-            if (_flow[e] >= delta) {
-              v = _graph.source(e);
-              switch(heap.state(v)) {



More information about the Lemon-commits mailing list