Changeset 472:052af4060f3e in lemon-0.x for src

Timestamp:

04/29/04 17:58:34 (20 years ago)

Author:

marci

Branch:

default

Phase:

public

Convert:

svn:c9d7d8f5-90d6-0310-b91f-818b3a526b0e/lemon/trunk@629

Message:

preflow, maxflow

Location:

src/work

Files:

• jacint/preflow.h (modified) (5 diffs)
• marci/edmonds_karp_demo.cc (modified) (6 diffs)

src/work/jacint/preflow.h

@@ -45 +45 @@
 
 #include <graph_wrapper.h>
+#include <bfs_iterator.h>
+#include <invalid.h>
+#include <maps.h>
+#include <for_each_macros.h>
+
 
 namespace hugo {
@@ -112 +117 @@
     bool augmentOnBlockingFlow2();
 
-    //Returns the maximum value of a flow.
-    Num flowValue() {
-      return excess[t];
+    /// Returns the actual flow value.
+    /// More precisely, it returns the negative excess of s, thus 
+    /// this works also for preflows.
+    Num flowValue() { 
+      Num a=0;
+      FOR_EACH_INC_LOC(OutEdgeIt, e, *g, s) a+=(*flow)[e];
+      FOR_EACH_INC_LOC(InEdgeIt, e, *g, s) a-=(*flow)[e];
+      return a;
     }
 
@@ -434 +444 @@
     void relabel(Node w, int newlevel, VecStack& active,  
                  VecNode& level_list, NNMap& left, 
-                 NNMap& right, int& b, int& k, bool what_heur ) {
+                 NNMap& right, int& b, int& k, bool what_heur ) 
+    {
 
       Num lev=level[w]; 
@@ -498 +509 @@
       
     } //relabel
+
+
+    template<typename MapGraphWrapper> 
+    class DistanceMap {
+    protected:
+      const MapGraphWrapper* g;
+      typename MapGraphWrapper::template NodeMap<int> dist; 
+    public:
+      DistanceMap(MapGraphWrapper& _g) : g(&_g), dist(*g, g->nodeNum()) { }
+      void set(const typename MapGraphWrapper::Node& n, int a) { 
+        dist.set(n, a); 
+      }
+      int operator[](const typename MapGraphWrapper::Node& n) 
+        { return dist[n]; }
+//       int get(const typename MapGraphWrapper::Node& n) const { 
+//      return dist[n]; }
+//       bool get(const typename MapGraphWrapper::Edge& e) const { 
+//      return (dist.get(g->tail(e))<dist.get(g->head(e))); }
+      bool operator[](const typename MapGraphWrapper::Edge& e) const { 
+        return (dist[g->tail(e)]<dist[g->head(e)]); 
+      }
+    };
     
   };
@@ -675 +708 @@
 
 
-
-
-
-
-
+  template <typename Graph, typename Num, typename CapMap, typename FlowMap>
+  bool Preflow<Graph, Num, CapMap, FlowMap>::augmentOnShortestPath() 
+  {
+      ResGW res_graph(*g, *capacity, *flow);
+      bool _augment=false;
+      
+      //ReachedMap level(res_graph);
+      FOR_EACH_LOC(typename Graph::NodeIt, e, *g) level.set(e, 0);
+      BfsIterator<ResGW, ReachedMap> bfs(res_graph, level);
+      bfs.pushAndSetReached(s);
+        
+      typename ResGW::template NodeMap<ResGWEdge> pred(res_graph); 
+      pred.set(s, INVALID);
+      
+      typename ResGW::template NodeMap<Num> free(res_graph);
+        
+      //searching for augmenting path
+      while ( !bfs.finished() ) { 
+        ResGWOutEdgeIt e=bfs;
+        if (res_graph.valid(e) && bfs.isBNodeNewlyReached()) {
+          Node v=res_graph.tail(e);
+          Node w=res_graph.head(e);
+          pred.set(w, e);
+          if (res_graph.valid(pred[v])) {
+            free.set(w, std::min(free[v], res_graph.resCap(e)));
+          } else {
+            free.set(w, res_graph.resCap(e)); 
+          }
+          if (res_graph.head(e)==t) { _augment=true; break; }
+        }
+        
+        ++bfs;
+      } //end of searching augmenting path
+
+      if (_augment) {
+        Node n=t;
+        Num augment_value=free[t];
+        while (res_graph.valid(pred[n])) { 
+          ResGWEdge e=pred[n];
+          res_graph.augment(e, augment_value); 
+          n=res_graph.tail(e);
+        }
+      }
+
+      return _augment;
+    }
+
+
+
+
+
+
+
+
+
+  template <typename Graph, typename Num, typename CapMap, typename FlowMap>
+  template<typename MutableGraph> 
+  bool Preflow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow() 
+  {      
+      typedef MutableGraph MG;
+      bool _augment=false;
+
+      ResGW res_graph(*g, *capacity, *flow);
+
+      //bfs for distances on the residual graph
+      //ReachedMap level(res_graph);
+      FOR_EACH_LOC(typename Graph::NodeIt, e, *g) level.set(e, 0);
+      BfsIterator<ResGW, ReachedMap> bfs(res_graph, level);
+      bfs.pushAndSetReached(s);
+      typename ResGW::template NodeMap<int> 
+        dist(res_graph); //filled up with 0's
+
+      //F will contain the physical copy of the residual graph
+      //with the set of edges which are on shortest paths
+      MG F;
+      typename ResGW::template NodeMap<typename MG::Node> 
+        res_graph_to_F(res_graph);
+      {
+        typename ResGW::NodeIt n;
+        for(res_graph.first(n); res_graph.valid(n); res_graph.next(n)) {
+          res_graph_to_F.set(n, F.addNode());
+        }
+      }
+
+      typename MG::Node sF=res_graph_to_F[s];
+      typename MG::Node tF=res_graph_to_F[t];
+      typename MG::template EdgeMap<ResGWEdge> original_edge(F);
+      typename MG::template EdgeMap<Num> residual_capacity(F);
+
+      while ( !bfs.finished() ) { 
+        ResGWOutEdgeIt e=bfs;
+        if (res_graph.valid(e)) {
+          if (bfs.isBNodeNewlyReached()) {
+            dist.set(res_graph.head(e), dist[res_graph.tail(e)]+1);
+            typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.tail(e)], res_graph_to_F[res_graph.head(e)]);
+            original_edge.update();
+            original_edge.set(f, e);
+            residual_capacity.update();
+            residual_capacity.set(f, res_graph.resCap(e));
+          } else {
+            if (dist[res_graph.head(e)]==(dist[res_graph.tail(e)]+1)) {
+              typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.tail(e)], res_graph_to_F[res_graph.head(e)]);
+              original_edge.update();
+              original_edge.set(f, e);
+              residual_capacity.update();
+              residual_capacity.set(f, res_graph.resCap(e));
+            }
+          }
+        }
+        ++bfs;
+      } //computing distances from s in the residual graph
+
+      bool __augment=true;
+
+      while (__augment) {
+        __augment=false;
+        //computing blocking flow with dfs
+        DfsIterator< MG, typename MG::template NodeMap<bool> > dfs(F);
+        typename MG::template NodeMap<typename MG::Edge> pred(F);
+        pred.set(sF, INVALID);
+        //invalid iterators for sources
+
+        typename MG::template NodeMap<Num> free(F);
+
+        dfs.pushAndSetReached(sF);      
+        while (!dfs.finished()) {
+          ++dfs;
+          if (F.valid(/*typename MG::OutEdgeIt*/(dfs))) {
+            if (dfs.isBNodeNewlyReached()) {
+              typename MG::Node v=F.aNode(dfs);
+              typename MG::Node w=F.bNode(dfs);
+              pred.set(w, dfs);
+              if (F.valid(pred[v])) {
+                free.set(w, std::min(free[v], residual_capacity[dfs]));
+              } else {
+                free.set(w, residual_capacity[dfs]); 
+              }
+              if (w==tF) { 
+                __augment=true; 
+                _augment=true;
+                break; 
+              }
+              
+            } else {
+              F.erase(/*typename MG::OutEdgeIt*/(dfs));
+            }
+          } 
+        }
+
+        if (__augment) {
+          typename MG::Node n=tF;
+          Num augment_value=free[tF];
+          while (F.valid(pred[n])) { 
+            typename MG::Edge e=pred[n];
+            res_graph.augment(original_edge[e], augment_value); 
+            n=F.tail(e);
+            if (residual_capacity[e]==augment_value) 
+              F.erase(e); 
+            else 
+              residual_capacity.set(e, residual_capacity[e]-augment_value);
+          }
+        }
+        
+      }
+            
+      return _augment;
+    }
+
+
+
+
+
+
+  template <typename Graph, typename Num, typename CapMap, typename FlowMap>
+  bool Preflow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow2() 
+  {
+      bool _augment=false;
+
+      ResGW res_graph(*g, *capacity, *flow);
+      
+      //ReachedMap level(res_graph);
+      FOR_EACH_LOC(typename Graph::NodeIt, e, *g) level.set(e, 0);
+      BfsIterator<ResGW, ReachedMap> bfs(res_graph, level);
+
+      bfs.pushAndSetReached(s);
+      DistanceMap<ResGW> dist(res_graph);
+      while ( !bfs.finished() ) { 
+        ResGWOutEdgeIt e=bfs;
+        if (res_graph.valid(e) && bfs.isBNodeNewlyReached()) {
+          dist.set(res_graph.head(e), dist[res_graph.tail(e)]+1);
+        }
+        ++bfs;
+      } //computing distances from s in the residual graph
+
+      //Subgraph containing the edges on some shortest paths
+      ConstMap<typename ResGW::Node, bool> true_map(true);
+      typedef SubGraphWrapper<ResGW, ConstMap<typename ResGW::Node, bool>, 
+        DistanceMap<ResGW> > FilterResGW;
+      FilterResGW filter_res_graph(res_graph, true_map, dist);
+
+      //Subgraph, which is able to delete edges which are already 
+      //met by the dfs
+      typename FilterResGW::template NodeMap<typename FilterResGW::OutEdgeIt> 
+        first_out_edges(filter_res_graph);
+      typename FilterResGW::NodeIt v;
+      for(filter_res_graph.first(v); filter_res_graph.valid(v); 
+          filter_res_graph.next(v)) 
+      {
+        typename FilterResGW::OutEdgeIt e;
+        filter_res_graph.first(e, v);
+        first_out_edges.set(v, e);
+      }
+      typedef ErasingFirstGraphWrapper<FilterResGW, typename FilterResGW::
+        template NodeMap<typename FilterResGW::OutEdgeIt> > ErasingResGW;
+      ErasingResGW erasing_res_graph(filter_res_graph, first_out_edges);
+
+      bool __augment=true;
+
+      while (__augment) {
+
+        __augment=false;
+        //computing blocking flow with dfs
+        DfsIterator< ErasingResGW, 
+          typename ErasingResGW::template NodeMap<bool> > 
+          dfs(erasing_res_graph);
+        typename ErasingResGW::
+          template NodeMap<typename ErasingResGW::OutEdgeIt> 
+          pred(erasing_res_graph); 
+        pred.set(s, INVALID);
+        //invalid iterators for sources
+
+        typename ErasingResGW::template NodeMap<Num> 
+          free1(erasing_res_graph);
+
+        dfs.pushAndSetReached(
+          typename ErasingResGW::Node(
+            typename FilterResGW::Node(
+              typename ResGW::Node(s)
+              )
+            )
+          );
+        while (!dfs.finished()) {
+          ++dfs;
+          if (erasing_res_graph.valid(
+                typename ErasingResGW::OutEdgeIt(dfs))) 
+          { 
+            if (dfs.isBNodeNewlyReached()) {
+          
+              typename ErasingResGW::Node v=erasing_res_graph.aNode(dfs);
+              typename ErasingResGW::Node w=erasing_res_graph.bNode(dfs);
+
+              pred.set(w, /*typename ErasingResGW::OutEdgeIt*/(dfs));
+              if (erasing_res_graph.valid(pred[v])) {
+                free1.set(w, std::min(free1[v], res_graph.resCap(
+                                       typename ErasingResGW::OutEdgeIt(dfs))));
+              } else {
+                free1.set(w, res_graph.resCap(
+                           typename ErasingResGW::OutEdgeIt(dfs))); 
+              }
+              
+              if (w==t) { 
+                __augment=true; 
+                _augment=true;
+                break; 
+              }
+            } else {
+              erasing_res_graph.erase(dfs);
+            }
+          }
+        }       
+
+        if (__augment) {
+          typename ErasingResGW::Node n=typename FilterResGW::Node(typename ResGW::Node(t));
+//        typename ResGW::NodeMap<Num> a(res_graph);
+//        typename ResGW::Node b;
+//        Num j=a[b];
+//        typename FilterResGW::NodeMap<Num> a1(filter_res_graph);
+//        typename FilterResGW::Node b1;
+//        Num j1=a1[b1];
+//        typename ErasingResGW::NodeMap<Num> a2(erasing_res_graph);
+//        typename ErasingResGW::Node b2;
+//        Num j2=a2[b2];
+          Num augment_value=free1[n];
+          while (erasing_res_graph.valid(pred[n])) { 
+            typename ErasingResGW::OutEdgeIt e=pred[n];
+            res_graph.augment(e, augment_value);
+            n=erasing_res_graph.tail(e);
+            if (res_graph.resCap(e)==0)
+              erasing_res_graph.erase(e);
+        }
+      }
+      
+      } //while (__augment) 
+            
+      return _augment;
+    }
 
 

src/work/marci/edmonds_karp_demo.cc

@@ -6 +6 @@
 #include <smart_graph.h>
 #include <dimacs.h>
-#include <edmonds_karp.h>
+//#include <edmonds_karp.h>
 #include <time_measure.h>
 //#include <graph_wrapper.h>
 #include <preflow.h>
-#include <preflow_res.h>
+//#include <preflow_res.h>
 #include <for_each_macros.h>
 
@@ -73 +73 @@
   Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > 
     pre_flow_test(G, s, t, cap, flow/*, true*/);
-  Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > 
-    pre_flow_ize(G, s, t, cap, flow/*, false*/);
+  //  Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > 
+  //  pre_flow_ize(G, s, t, cap, flow/*, false*/);
 //   PreflowRes<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > 
 //     pre_flow_res(G, s, t, cap, flow/*, true*/);
-  MaxFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > 
-    max_flow_test(G, s, t, cap, flow);
+  //MaxFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > 
+  //  max_flow_test(G, s, t, cap, flow);
 
   {
@@ -92 +92 @@
     FOR_EACH_LOC(Graph::EdgeIt, e, G) flow.set(e, 0);
     ts.reset();
-    pre_flow_ize.preflow(Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >::GEN_FLOW);
+    pre_flow_test.preflow(Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >::GEN_FLOW);
     std::cout << "elapsed time: " << ts << std::endl;
-    std::cout << "flow value: "<< pre_flow_ize.flowValue() << std::endl;
+    std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl;
   }
 
@@ -111 +111 @@
     ts.reset();
     int i=0;
-    while (max_flow_test.augmentOnBlockingFlow<MutableGraph>()) { ++i; }
+    while (pre_flow_test.augmentOnBlockingFlow<MutableGraph>()) { ++i; }
     std::cout << "elapsed time: " << ts << std::endl;
     std::cout << "number of augmentation phases: " << i << std::endl; 
-    std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
+    std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl;
   }
 
-  {
-    std::cout << "faster physical blocking flow augmentation ..." << std::endl;
-    FOR_EACH_LOC(Graph::EdgeIt, e, G) flow.set(e, 0);
-    ts.reset();
-    int i=0;
-    while (max_flow_test.augmentOnBlockingFlow1<MutableGraph>()) { ++i; }
-    std::cout << "elapsed time: " << ts << std::endl;
-    std::cout << "number of augmentation phases: " << i << std::endl; 
-    std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
-  }
+//   {
+//     std::cout << "faster physical blocking flow augmentation ..." << std::endl;
+//     FOR_EACH_LOC(Graph::EdgeIt, e, G) flow.set(e, 0);
+//     ts.reset();
+//     int i=0;
+//     while (max_flow_test.augmentOnBlockingFlow1<MutableGraph>()) { ++i; }
+//     std::cout << "elapsed time: " << ts << std::endl;
+//     std::cout << "number of augmentation phases: " << i << std::endl; 
+//     std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
+//   }
 
   {
@@ -133 +133 @@
     ts.reset();
     int i=0;
-    while (max_flow_test.augmentOnBlockingFlow2()) { ++i; }
+    while (pre_flow_test.augmentOnBlockingFlow2()) { ++i; }
     std::cout << "elapsed time: " << ts << std::endl;
     std::cout << "number of augmentation phases: " << i << std::endl; 
-    std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
+    std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl;
   }
 
@@ -144 +144 @@
     ts.reset();
     int i=0;
-    while (max_flow_test.augmentOnShortestPath()) { ++i; }
+    while (pre_flow_test.augmentOnShortestPath()) { ++i; }
     std::cout << "elapsed time: " << ts << std::endl;
     std::cout << "number of augmentation phases: " << i << std::endl; 
-    std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
+    std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl;
   }