# HG changeset patch
# User jacint
# Date 1091183045 0
# Node ID 8e933219691e3641593e27f1618993d0ed0c26a2
# Parent  a0e497db23ee95cdb960c0ae263f6f12fb17e382
bug fixing

diff -r a0e497db23ee -r 8e933219691e src/hugo/max_flow.h
--- a/src/hugo/max_flow.h	Thu Jul 29 17:23:55 2004 +0000
+++ b/src/hugo/max_flow.h	Fri Jul 30 10:24:05 2004 +0000
@@ -1,17 +1,15 @@
 // -*- C++ -*-
-#ifndef HUGO_MAX_FLOW_NO_STACK_H
-#define HUGO_MAX_FLOW_NO_STACK_H
+#ifndef HUGO_MAX_FLOW_H
+#define HUGO_MAX_FLOW_H
 
 #include <vector>
 #include <queue>
-//#include <stack>
 
 #include <hugo/graph_wrapper.h>
 #include <hugo/invalid.h>
 #include <hugo/maps.h>
 
 /// \file
-/// \brief The same as max_flow.h, but without using stl stack for the active nodes. Only for test.
 /// \ingroup galgs
 
 namespace hugo {
@@ -51,7 +49,6 @@
     typedef typename Graph::OutEdgeIt OutEdgeIt;
     typedef typename Graph::InEdgeIt InEdgeIt;
 
-    //    typedef typename std::vector<std::stack<Node> > VecStack;
     typedef typename std::vector<Node> VecFirst;
     typedef typename Graph::template NodeMap<Node> NNMap;
     typedef typename std::vector<Node> VecNode;
@@ -66,7 +63,6 @@
     //typedef ExpResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW;
     typedef typename ResGW::OutEdgeIt ResGWOutEdgeIt;
     typedef typename ResGW::Edge ResGWEdge;
-    //typedef typename ResGW::template NodeMap<bool> ReachedMap;
     typedef typename Graph::template NodeMap<int> ReachedMap;
 
 
@@ -110,7 +106,7 @@
       AFTER_PRE_FLOW_PHASE_2
     };
 
-    /// Don not needle this flag only if necessary.
+    /// Do not needle this flag only if necessary.
     StatusEnum status;
 
 //     int number_of_augmentations;
@@ -188,7 +184,7 @@
     ///The preflow algorithm consists of two phases, this method runs the
     ///first phase. After the first phase the maximum flow value and a
     ///minimum value cut can already be computed, though a maximum flow
-    ///is net yet obtained. So after calling this method \ref flowValue
+    ///is not yet obtained. So after calling this method \ref flowValue
     ///and \ref actMinCut gives proper results.
     ///\warning: \ref minCut, \ref minMinCut and \ref maxMinCut do not
     ///give minimum value cuts unless calling \ref preflowPhase2.
@@ -223,65 +219,9 @@
       VecNode level_list(n,INVALID);
       //List of the nodes in level i<n, set to n.
 
-      NodeIt v;
-      for(g->first(v); g->valid(v); g->next(v)) level.set(v,n);
-      //setting each node to level n
-
-      if ( fe == NO_FLOW ) {
-	EdgeIt e;
-	for(g->first(e); g->valid(e); g->next(e)) flow->set(e,0);
-      }
-
-      switch (fe) { //computing the excess
-      case PRE_FLOW:
-	{
-	  NodeIt v;
-	  for(g->first(v); g->valid(v); g->next(v)) {
-	    Num exc=0;
-
-	    InEdgeIt e;
-	    for(g->first(e,v); g->valid(e); g->next(e)) exc+=(*flow)[e];
-	    OutEdgeIt f;
-	    for(g->first(f,v); g->valid(f); g->next(f)) exc-=(*flow)[f];
-
-	    excess.set(v,exc);
-
-	    //putting the active nodes into the stack
-	    int lev=level[v];
-	    if ( exc > 0 && lev < n && v != t ) 
-	      {
-		next.set(v,first[lev]);
-		first[lev]=v;
-	      }
-	  }
-	  break;
-	}
-      case GEN_FLOW:
-	{
-	  NodeIt v;
-	  for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0);
-
-	  Num exc=0;
-	  InEdgeIt e;
-	  for(g->first(e,t); g->valid(e); g->next(e)) exc+=(*flow)[e];
-	  OutEdgeIt f;
-	  for(g->first(f,t); g->valid(f); g->next(f)) exc-=(*flow)[f];
-	  excess.set(t,exc);
-	  break;
-	}
-      case ZERO_FLOW:
-      case NO_FLOW:
-	{
-	  NodeIt v;
-	  for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0);
-	  break;
-	}
-      }
-
       preflowPreproc(fe, next, first, level_list, left, right);
       //End of preprocessing
 
-
       //Push/relabel on the highest level active nodes.
       while ( true ) {
 	if ( b == 0 ) {
@@ -394,7 +334,7 @@
 	    first[newlevel]=w;
 	    b=newlevel;
 	  }
-	}  // if stack[b] is nonempty
+	} 
       } // while(true)
 
       status=AFTER_PRE_FLOW_PHASE_2;
@@ -413,15 +353,7 @@
       return a;
       //marci figyu: excess[t] epp ezt adja preflow 1. fazisa utan   
     }
-    Num flowValue2() const {
-      return excess[t];
-//       Num a=0;
-//       for(InEdgeIt e(*g,t);g->valid(e);g->next(e)) a+=(*flow)[e];
-//       for(OutEdgeIt e(*g,t);g->valid(e);g->next(e)) a-=(*flow)[e];
-//       return a;
-//       //marci figyu: excess[t] epp ezt adja preflow 1. fazisa utan  
-      
-    }
+
 
     ///Returns a minimum value cut after calling \ref preflowPhase1.
 
@@ -652,13 +584,51 @@
     }
 
 
+
     void preflowPreproc(FlowEnum fe, NNMap& next, VecFirst& first,
 			VecNode& level_list, NNMap& left, NNMap& right)
     {
+      switch (fe) { //setting excess
+	case NO_FLOW: 
+	{
+	  EdgeIt e;
+	  for(g->first(e); g->valid(e); g->next(e)) flow->set(e,0);
+	  
+	  NodeIt v;
+	  for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0);
+	  break;
+	}
+	case ZERO_FLOW: 
+	{
+	  NodeIt v;
+	  for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0);
+	  break;
+	}
+	case GEN_FLOW:
+	{
+	  NodeIt v;
+	  for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0);
+
+	  Num exc=0;
+	  InEdgeIt e;
+	  for(g->first(e,t); g->valid(e); g->next(e)) exc+=(*flow)[e];
+	  OutEdgeIt f;
+	  for(g->first(f,t); g->valid(f); g->next(f)) exc-=(*flow)[f];
+	  excess.set(t,exc);
+	  break;
+	}
+	default: break;
+      }
+
+      NodeIt v;
+      for(g->first(v); g->valid(v); g->next(v)) level.set(v,n);
+      //setting each node to level n
+      
       std::queue<Node> bfs_queue;
 
+
       switch (fe) {
-      case NO_FLOW:   //flow is already set to const zero in this case
+      case NO_FLOW:   //flow is already set to const zero
       case ZERO_FLOW:
 	{
 	  //Reverse_bfs from t, to find the starting level.
@@ -693,8 +663,8 @@
 	      if ( c <= 0 ) continue;
 	      Node w=g->head(e);
 	      if ( level[w] < n ) {
-		if ( excess[w] <= 0 && w!=t ) 
-		  {
+		if ( excess[w] <= 0 && w!=t ) //putting into the stack
+		  { 
 		    next.set(w,first[level[w]]);
 		    first[level[w]]=w;
 		  }
@@ -706,6 +676,83 @@
 	}
 
       case GEN_FLOW:
+	{
+	  //Reverse_bfs from t in the residual graph,
+	  //to find the starting level.
+	  level.set(t,0);
+	  bfs_queue.push(t);
+
+	  while (!bfs_queue.empty()) {
+
+	    Node v=bfs_queue.front();
+	    bfs_queue.pop();
+	    int l=level[v]+1;
+
+	    InEdgeIt e;
+	    for(g->first(e,v); g->valid(e); g->next(e)) {
+	      if ( (*capacity)[e] <= (*flow)[e] ) continue;
+	      Node w=g->tail(e);
+	      if ( level[w] == n && w != s ) {
+		bfs_queue.push(w);
+		Node z=level_list[l];
+		if ( g->valid(z) ) left.set(z,w);
+		right.set(w,z);
+		level_list[l]=w;
+		level.set(w, l);
+	      }
+	    }
+
+	    OutEdgeIt f;
+	    for(g->first(f,v); g->valid(f); g->next(f)) {
+	      if ( 0 >= (*flow)[f] ) continue;
+	      Node w=g->head(f);
+	      if ( level[w] == n && w != s ) {
+		bfs_queue.push(w);
+		Node z=level_list[l];
+		if ( g->valid(z) ) left.set(z,w);
+		right.set(w,z);
+		level_list[l]=w;
+		level.set(w, l);
+	      }
+	    }
+	  }
+
+	  //the starting flow
+	  OutEdgeIt e;
+	  for(g->first(e,s); g->valid(e); g->next(e))
+	    {
+	      Num rem=(*capacity)[e]-(*flow)[e];
+	      if ( rem <= 0 ) continue;
+	      Node w=g->head(e);
+	      if ( level[w] < n ) {
+		if ( excess[w] <= 0 && w!=t ) //putting into the stack
+		  {
+		    next.set(w,first[level[w]]);
+		    first[level[w]]=w;
+		  }   
+		flow->set(e, (*capacity)[e]);
+		excess.set(w, excess[w]+rem);
+	      }
+	    }
+
+	  InEdgeIt f;
+	  for(g->first(f,s); g->valid(f); g->next(f))
+	    {
+	      if ( (*flow)[f] <= 0 ) continue;
+	      Node w=g->tail(f);
+	      if ( level[w] < n ) {
+		if ( excess[w] <= 0 && w!=t )
+		  {
+		    next.set(w,first[level[w]]);
+		    first[level[w]]=w;
+		  }  
+		excess.set(w, excess[w]+(*flow)[f]);
+		flow->set(f, 0);
+	      }
+	    }
+	  break;
+	} //case GEN_FLOW
+    
       case PRE_FLOW:
 	{
 	  //Reverse_bfs from t in the residual graph,
@@ -757,11 +804,6 @@
 	      if ( rem <= 0 ) continue;
 	      Node w=g->head(e);
 	      if ( level[w] < n ) {
-		if ( excess[w] <= 0 && w!=t )
-		  {
-		    next.set(w,first[level[w]]);
-		    first[level[w]]=w;
-		  }   
 		flow->set(e, (*capacity)[e]);
 		excess.set(w, excess[w]+rem);
 	      }
@@ -773,22 +815,36 @@
 	      if ( (*flow)[f] <= 0 ) continue;
 	      Node w=g->tail(f);
 	      if ( level[w] < n ) {
-		if ( excess[w] <= 0 && w!=t )
-		  {
-		    next.set(w,first[level[w]]);
-		    first[level[w]]=w;
-		  }   
 		excess.set(w, excess[w]+(*flow)[f]);
 		flow->set(f, 0);
 	      }
 	    }
+	  
+	  NodeIt w; //computing the excess
+	  for(g->first(w); g->valid(w); g->next(w)) {
+	    Num exc=0;
+
+	    InEdgeIt e;
+	    for(g->first(e,w); g->valid(e); g->next(e)) exc+=(*flow)[e];
+	    OutEdgeIt f;
+	    for(g->first(f,w); g->valid(f); g->next(f)) exc-=(*flow)[f];
+
+	    excess.set(w,exc);
+
+	    //putting the active nodes into the stack
+	    int lev=level[w];
+	    if ( exc > 0 && lev < n && w != t ) 
+	      {
+		next.set(w,first[lev]);
+		first[lev]=w;
+	      }
+	  }
 	  break;
 	} //case PRE_FLOW
       }
     } //preflowPreproc
 
 
-
     void relabel(Node w, int newlevel, NNMap& next, VecFirst& first,
 		 VecNode& level_list, NNMap& left,
 		 NNMap& right, int& b, int& k, bool what_heur )
@@ -852,6 +908,35 @@
 	}
       }
     } //relabel
+
+    void printexcess() {////
+      std::cout << "Excesses:" <<std::endl;
+
+      NodeIt v;
+      for(g->first(v); g->valid(v); g->next(v)) {
+	std::cout << 1+(g->id(v)) << ":" << excess[v]<<std::endl; 
+      }
+    }
+
+ void printlevel() {////
+      std::cout << "Levels:" <<std::endl;
+
+      NodeIt v;
+      for(g->first(v); g->valid(v); g->next(v)) {
+	std::cout << 1+(g->id(v)) << ":" << level[v]<<std::endl; 
+      }
+    }
+
+void printactive() {////
+      std::cout << "Levels:" <<std::endl;
+
+      NodeIt v;
+      for(g->first(v); g->valid(v); g->next(v)) {
+	std::cout << 1+(g->id(v)) << ":" << level[v]<<std::endl; 
+      }
+    }
+
+
   };  //class MaxFlow
 } //namespace hugo