lemon-0.x: comparison src/work/jacint/preflow

equal deleted inserted replaced

-:723ab01137f3
+:435e37fdebf6
 // -*- C++ -*-
 /*
-preflow_hl4.h
+preflow_h5.h
 by jacint.
-Runs the two phase highest label preflow push algorithm. In phase 0
+Heuristics:
-we maintain in a list the nodes in level i < n, and we maintain a
+2 phase
-bound k on the max level i < n containing a node, so we can do
+gap
-the gap heuristic fast. Phase 1 is the same. (The algorithm is the
+list 'level_list' on the nodes on level i implemented by hand
-same as preflow.hl3, the only diff is that here we use the gap
+highest label
-heuristic with the list of the nodes on level i, and not a bfs form the
+relevel: in phase 0, after BFS*n relabels, it runs a reverse
-upgraded node.)
+bfs from t in the res graph to relevel the nodes reachable from t.
+BFS is initialized to 20
-In phase 1 we shift everything downwards by n.
+Due to the last heuristic, this algorithm is quite fast on very
-Member functions:
+sparse graphs, but relatively bad on even the dense graphs.
-void run() : runs the algorithm
+'NodeMap<bool> cut' is a member, in this way we can count it fast, after
+the algorithm was run.
-The following functions should be used after run() was already run.
+The constructor runs the algorithm.
-T maxflow() : returns the value of a maximum flow
+Members:
-T flowonedge(EdgeIt e) : for a fixed maximum flow x it returns x(e)
+T maxFlow() : returns the value of a maximum flow
-FlowMap allflow() : returns a maximum flow
+T flowOnEdge(EdgeIt e) : for a fixed maximum flow x it returns x(e)
-void allflow(FlowMap& _flow ) : returns a maximum flow
+FlowMap Flow() : returns the fixed maximum flow x
-void mincut(CutMap& M) : sets M to the characteristic vector of a
-minimum cut. M should be a map of bools initialized to false.
+void Flow(FlowMap& _flow ) : returns the fixed maximum flow x
-void min_mincut(CutMap& M) : sets M to the characteristic vector of the
+void minMinCut(CutMap& M) : sets M to the characteristic vector of the
 minimum min cut. M should be a map of bools initialized to false.
-void max_mincut(CutMap& M) : sets M to the characteristic vector of the
+void maxMinCut(CutMap& M) : sets M to the characteristic vector of the
 maximum min cut. M should be a map of bools initialized to false.
+void minCut(CutMap& M) : fast function, sets M to the characteristic
+vector of a minimum cut.
+Different member from the other preflow_hl-s (here we have a member
+'NodeMap<bool> cut').
+CutMap minCut() : fast function, giving the characteristic
+vector of a minimum cut.
 */
 #ifndef PREFLOW_HL4_H
 #define PREFLOW_HL4_H
+#define BFS 20
 #include <vector>
-#include <stack>
 #include <queue>
+#include <time_measure.h>  //for test
 namespace hugo {
 template <typename Graph, typename T,
 typename FlowMap=typename Graph::EdgeMap<T>,
+typename CutMap=typename Graph::NodeMap<bool>,
 typename CapMap=typename Graph::EdgeMap<T> >
 class preflow_hl4 {
 typedef typename Graph::NodeIt NodeIt;
 typedef typename Graph::EdgeIt EdgeIt;
 Graph& G;
 NodeIt s;
 NodeIt t;
 FlowMap flow;
 CapMap& capacity;
+CutMap cut;
 T value;
 public:
+double time;
 preflow_hl4(Graph& _G, NodeIt _s, NodeIt _t, CapMap& _capacity) :
-G(_G), s(_s), t(_t), flow(_G, 0), capacity(_capacity) { }
+G(_G), s(_s), t(_t), flow(_G, 0), capacity(_capacity),
+cut(G, false) {
-void run() {
+bool phase=0;   //phase 0 is the 1st phase, phase 1 is the 2nd
-bool phase=0;
 int n=G.nodeNum();
+int relabel=0;
+int heur=(int)BFS*n;
 int k=n-2;
 int b=k;
 /*
 	b is a bound on the highest level of the stack.
 	k is a bound on the highest nonempty level i < n.
 */
 typename Graph::NodeMap<int> level(G,n);
 typename Graph::NodeMap<T> excess(G);
-std::vector<std::stack<NodeIt> > stack(n);
+std::vector<NodeIt> active(n);
+typename Graph::NodeMap<NodeIt> next(G);
 //Stack of the active nodes in level i < n.
 //We use it in both phases.
 typename Graph::NodeMap<NodeIt> left(G);
 typename Graph::NodeMap<NodeIt> right(G);
 	bfs_queue.pop();
 	int l=level.get(v)+1;
 	for(InEdgeIt e=G.template first<InEdgeIt>(v); e.valid(); ++e) {
 	  NodeIt w=G.tail(e);
-	  if ( level.get(w) == n ) {
+	  if ( level.get(w) == n && w !=s ) {
 	    bfs_queue.push(w);
 	    NodeIt first=level_list[l];
 	    if ( first != 0 ) left.set(first,w);
 	    right.set(w,first);
 	    level_list[l]=w;
 	{
 	  T c=capacity.get(e);
 	  if ( c == 0 ) continue;
 	  NodeIt w=G.head(e);
 	  if ( level.get(w) < n ) {
-	    if ( excess.get(w) == 0 && w!=t ) stack[level.get(w)].push(w);
+	    if ( excess.get(w) == 0 && w!=t ) {
+	      next.set(w,active[level.get(w)]);
+	      active[level.get(w)]=w;
+	    }
 	    flow.set(e, c);
 	    excess.set(w, excess.get(w)+c);
 	  }
 	}
 /*
 	  /*
 	    In the end of phase 0 we apply a bfs from s in
 	    the residual graph.
 	  */
 	  phase=1;
+	  //Now have a min cut.
+	  for( EachNodeIt v=G.template first<EachNodeIt>();
+	       v.valid(); ++v)
+	    if (level.get(v) >= n ) cut.set(v,true);
+	  time=currTime();
 	  level.set(s,0);
 	  std::queue<NodeIt> bfs_queue;
 	  bfs_queue.push(s);
 	  while (!bfs_queue.empty()) {
 	      if ( capacity.get(e) == flow.get(e) ) continue;
 	      NodeIt u=G.tail(e);
 	      if ( level.get(u) >= n ) {
 		bfs_queue.push(u);
 		level.set(u, l);
-		if ( excess.get(u) > 0 ) stack[l].push(u);
+		if ( excess.get(u) > 0 ) {
+		    next.set(u,active[l]);
+		    active[l]=u;
+		}
 	      }
 	    }
 	    for(OutEdgeIt e=G.template first<OutEdgeIt>(v); e.valid(); ++e) {
 	      if ( 0 == flow.get(e) ) continue;
 	      NodeIt u=G.head(e);
 	      if ( level.get(u) >= n ) {
 		bfs_queue.push(u);
 		level.set(u, l);
-		if ( excess.get(u) > 0 ) stack[l].push(u);
+		if ( excess.get(u) > 0 ) {
+		    next.set(u,active[l]);
+		    active[l]=u;
+		}
 	      }
 	    }
 	  }
 	  b=n-2;
 	}
-	if ( stack[b].empty() ) --b;
+	if ( active[b] == 0 ) --b;
 	else {
-	  NodeIt w=stack[b].top();        //w is a highest label active node.
+	  NodeIt w=active[b];
-	  stack[b].pop();
+	  active[b]=next.get(w);
 	  int lev=level.get(w);
 	  T exc=excess.get(w);
-	  int newlevel=n;          //In newlevel we bound the next level of w.
+	  int newlevel=n;          //bound on the next level of w.
 	  for(OutEdgeIt e=G.template first<OutEdgeIt>(w); e.valid(); ++e) {
 	    if ( flow.get(e) == capacity.get(e) ) continue;
 	    NodeIt v=G.head(e);
 	    //e=wv
 	    if( lev > level.get(v) ) {
 	      /*Push is allowed now*/
-	      if ( excess.get(v)==0 && v!=t && v!=s )
+	      if ( excess.get(v)==0 && v!=t && v!=s )  {
-		stack[level.get(v)].push(v);
+		int lev_v=level.get(v);
-	      /*v becomes active.*/
+		next.set(v,active[lev_v]);
+		active[lev_v]=v;
+	      }
 	      T cap=capacity.get(e);
 	      T flo=flow.get(e);
 	      T remcap=cap-flo;
 	    //e=vw
 	    if( lev > level.get(v) ) {
 	      /*Push is allowed now*/
-	      if ( excess.get(v)==0 && v!=t && v!=s )
+	      if ( excess.get(v)==0 && v!=t && v!=s ) {
-		stack[level.get(v)].push(v);
+		int lev_v=level.get(v);
-	      /*v becomes active.*/
+		next.set(v,active[lev_v]);
+		active[lev_v]=v;
+	      }
 	      T flo=flow.get(e);
 	      if ( flo >= exc ) {
 		/*A nonsaturating push.*/
 	if ( exc > 0 ) {
 	  //now 'lev' is the old level of w
 	  if ( phase ) {
 	    level.set(w,++newlevel);
-	    stack[newlevel].push(w);
+	    next.set(w,active[newlevel]);
+	    active[newlevel]=w;
 	    b=newlevel;
 	  } else {
 	    //unlacing
 	    NodeIt right_n=right.get(w);
 	    NodeIt left_n=left.get(w);
 		right.set(left_n, 0);
 	      } else {
 		level_list[lev]=0;
 	      }
 	    }
+	    //unlacing ends
 	    if ( level_list[lev]==0 ) {
 	      for (int i=lev; i!=k ; ) {
 	      if ( newlevel == n ) {
 		level.set(w,n);
 	      } else {
 		level.set(w,++newlevel);
-		stack[newlevel].push(w);
+		next.set(w,active[newlevel]);
+		active[newlevel]=w;
 		b=newlevel;
 		if ( k < newlevel ) ++k;
 		NodeIt first=level_list[newlevel];
 		if ( first != 0 ) left.set(first,w);
 		right.set(w,first);
 		left.set(w,0);
 		level_list[newlevel]=w;
 	      }
 	    }
+	    ++relabel;
+	    if ( relabel >= heur ) {
+	      relabel=0;
+	      b=n-2;
+	      k=b;
+	      for ( int i=1; i!=n; ++i ) {
+		active[i]=0;
+		level_list[i]=0;
+	      }
+	      //bfs from t in the res graph to relevel the nodes
+	      for( EachNodeIt v=G.template first<EachNodeIt>();
+		   v.valid(); ++v) level.set(v,n);
+	      level.set(t,0);
+	      std::queue<NodeIt> bfs_queue;
+	      bfs_queue.push(t);
+	      while (!bfs_queue.empty()) {
+		NodeIt v=bfs_queue.front();
+		bfs_queue.pop();
+		int l=level.get(v)+1;
+		for(InEdgeIt e=G.template first<InEdgeIt>(v);
+		    e.valid(); ++e) {
+		  if ( capacity.get(e) == flow.get(e) ) continue;
+		  NodeIt u=G.tail(e);
+		  if ( level.get(u) == n ) {
+		    bfs_queue.push(u);
+		    level.set(u, l);
+		    if ( excess.get(u) > 0 ) {
+		      next.set(u,active[l]);
+		      active[l]=u;
+		    }
+		    NodeIt first=level_list[l];
+		    if ( first != 0 ) left.set(first,w);
+		    right.set(w,first);
+		    left.set(w,0);
+		    level_list[l]=w;
+		  }
+		}
+		for(OutEdgeIt e=G.template first<OutEdgeIt>(v);
+		    e.valid(); ++e) {
+		  if ( 0 == flow.get(e) ) continue;
+		  NodeIt u=G.head(e);
+		  if ( level.get(u) == n ) {
+		    bfs_queue.push(u);
+		    level.set(u, l);
+		    if ( excess.get(u) > 0 ) {
+		      next.set(u,active[l]);
+		      active[l]=u;
+		    }
+		    NodeIt first=level_list[l];
+		    if ( first != 0 ) left.set(first,w);
+		    right.set(w,first);
+		    left.set(w,0);
+		    level_list[l]=w;
+		  }
+		}
+	      }
+	      level.set(s,n);
+	    }
 	  } //phase 0
 	} // if ( exc > 0 )
 	} // if stack[b] is nonempty
 } // while(true)
 /*
 Returns the maximum value of a flow.
 */
-T maxflow() {
+T maxFlow() {
 return value;
 }
 /*
 For the maximum flow x found by the algorithm, it returns the flow value on Edge e, i.e. x(e).
 */
-T flowonedge(EdgeIt e) {
+T flowOnEdge(EdgeIt e) {
 return flow.get(e);
 }
-FlowMap allflow() {
+FlowMap Flow() {
 return flow;
 }
-void allflow(FlowMap& _flow ) {
+void Flow(FlowMap& _flow ) {
 for(EachNodeIt v=G.template first<EachNodeIt>() ; v.valid(); ++v)
 	_flow.set(v,flow.get(v));
 }
 /*
 Returns the minimum min cut, by a bfs from s in the residual graph.
 */
-template<typename CutMap>
+template<typename _CutMap>
-void mincut(CutMap& M) {
+void minMinCut(_CutMap& M) {
 std::queue<NodeIt> queue;
 M.set(s,true);
 queue.push(s);
 /*
 Returns the maximum min cut, by a reverse bfs
 from t in the residual graph.
 */
-template<typename CutMap>
+template<typename _CutMap>
-void max_mincut(CutMap& M) {
+void maxMinCut(_CutMap& M) {
 std::queue<NodeIt> queue;
 M.set(t,true);
 queue.push(t);
 }
 }
+template<typename _CutMap>
-template<typename CutMap>
+void minCut(_CutMap& M) {
-void min_mincut(CutMap& M) {
+for( EachNodeIt v=G.template first<EachNodeIt>();
-mincut(M);
+	   v.valid(); ++v)
-}
+	M.set(v, cut.get(v));
+}
+CutMap minCut() {
+return cut;
+}
 };
-}//namespace hugo
+}//namespace marci
 #endif

changeset 113	cf7b01232d86
parent 105	a3c73e9b9b2e