src/work/jacint/preflow_push_hl.h
author alpar
Fri, 20 Feb 2004 21:57:39 +0000
changeset 106 0508d63fcc96
parent 101 d2ac583ed195
permissions -rw-r--r--
.
     1 // -*- C++ -*-
     2 
     3 //kerdesek: nem tudom lehet-e a 
     4 //kieleket csak a legf n szintu pontokra nezni.
     5 
     6 /*
     7 preflow_push_hl.h
     8 by jacint. 
     9 Runs the highest label variant of the preflow push algorithm with 
    10 running time O(n^2\sqrt(m)), and with the 'empty level' heuristic. 
    11 
    12 'A' is a parameter for the empty_level heuristic
    13 
    14 Member functions:
    15 
    16 void run() : runs the algorithm
    17 
    18  The following functions should be used after run() was already run.
    19 
    20 T maxflow() : returns the value of a maximum flow
    21 
    22 T flowonedge(EdgeIt e) : for a fixed maximum flow x it returns x(e) 
    23 
    24 FlowMap allflow() : returns the fixed maximum flow x
    25 
    26 void mincut(CutMap& M) : sets M to the characteristic vector of a 
    27      minimum cut. M should be a map of bools initialized to false.
    28 
    29 void min_mincut(CutMap& M) : sets M to the characteristic vector of the 
    30      minimum min cut. M should be a map of bools initialized to false.
    31 
    32 void max_mincut(CutMap& M) : sets M to the characteristic vector of the 
    33      maximum min cut. M should be a map of bools initialized to false.
    34 
    35 */
    36 
    37 #ifndef PREFLOW_PUSH_HL_H
    38 #define PREFLOW_PUSH_HL_H
    39 
    40 #define A 1
    41 
    42 #include <vector>
    43 #include <stack>
    44 #include <queue>
    45 
    46 namespace hugo {
    47 
    48   template <typename Graph, typename T, 
    49     typename FlowMap=typename Graph::EdgeMap<T>, typename CapMap=typename Graph::EdgeMap<T>, 
    50     typename IntMap=typename Graph::NodeMap<int>, typename TMap=typename Graph::NodeMap<T> >
    51   class preflow_push_hl {
    52     
    53     typedef typename Graph::NodeIt NodeIt;
    54     typedef typename Graph::EdgeIt EdgeIt;
    55     typedef typename Graph::EachNodeIt EachNodeIt;
    56     typedef typename Graph::OutEdgeIt OutEdgeIt;
    57     typedef typename Graph::InEdgeIt InEdgeIt;
    58     
    59     Graph& G;
    60     NodeIt s;
    61     NodeIt t;
    62     FlowMap flow;
    63     CapMap& capacity;  
    64     T value;
    65     
    66   public:
    67 
    68     preflow_push_hl(Graph& _G, NodeIt _s, NodeIt _t, CapMap& _capacity) :
    69       G(_G), s(_s), t(_t), flow(_G, 0), capacity(_capacity) { }
    70 
    71 
    72     
    73 
    74     void run() {
    75  
    76       int n=G.nodeNum(); 
    77       int b=n-2; 
    78       /*
    79 	b is a bound on the highest level of an active node. 
    80       */
    81 
    82       IntMap level(G,n);      
    83       TMap excess(G); 
    84 
    85       std::vector<int> numb(n);    
    86       /*
    87 	The number of nodes on level i < n. It is
    88 	initialized to n+1, because of the reverse_bfs-part.
    89       */
    90 
    91       std::vector<std::stack<NodeIt> > stack(2*n-1);    
    92       //Stack of the active nodes in level i.
    93 
    94 
    95       /*Reverse_bfs from t, to find the starting level.*/
    96       level.set(t,0);
    97       std::queue<NodeIt> bfs_queue;
    98       bfs_queue.push(t);
    99 
   100       while (!bfs_queue.empty()) {
   101 
   102 	NodeIt v=bfs_queue.front();	
   103 	bfs_queue.pop();
   104 	int l=level.get(v)+1;
   105 
   106 	for(InEdgeIt e=G.template first<InEdgeIt>(v); e.valid(); ++e) {
   107 	  NodeIt w=G.tail(e);
   108 	  if ( level.get(w) == n ) {
   109 	    bfs_queue.push(w);
   110 	    ++numb[l];
   111 	    level.set(w, l);
   112 	  }
   113 	}
   114       }
   115 	
   116       level.set(s,n);
   117 
   118 
   119 
   120       /* Starting flow. It is everywhere 0 at the moment. */     
   121       for(OutEdgeIt e=G.template first<OutEdgeIt>(s); e.valid(); ++e) 
   122 	{
   123 	  T c=capacity.get(e);
   124 	  if ( c == 0 ) continue;
   125 	  NodeIt w=G.head(e);
   126 	  if ( w!=s ) {	  
   127 	    if ( excess.get(w) == 0 && w!=t ) stack[level.get(w)].push(w); 
   128 	    flow.set(e, c); 
   129 	    excess.set(w, excess.get(w)+c);
   130 	  }
   131 	}
   132       
   133       /* 
   134 	 End of preprocessing 
   135       */
   136 
   137 
   138       /*
   139 	Push/relabel on the highest level active nodes.
   140       */
   141       /*While there exists an active node.*/
   142       while (b) { 
   143 	if ( stack[b].empty() ) { 
   144 	  --b;
   145 	  continue;
   146 	} 
   147 	
   148 	NodeIt w=stack[b].top();        //w is a highest label active node.
   149 	stack[b].pop();           
   150 	int lev=level.get(w);
   151 	int exc=excess.get(w);
   152 	int newlevel=2*n;      //In newlevel we bound the next level of w.
   153 	//vagy MAXINT	
   154 
   155 	//  if ( level.get(w) < n ) { //Nem tudom ez mukodik-e
   156 	  for(OutEdgeIt e=G.template first<OutEdgeIt>(w); e.valid(); ++e) {
   157 	    
   158 	    if ( flow.get(e) == capacity.get(e) ) continue; 
   159 	    NodeIt v=G.head(e);            
   160 	    //e=wv	    
   161 	    
   162 	    if( lev > level.get(v) ) {      
   163 	      /*Push is allowed now*/
   164 	      
   165 	      if ( excess.get(v)==0 && v != s && v !=t ) 
   166 		stack[level.get(v)].push(v); 
   167 	      /*v becomes active.*/
   168 	      
   169 	      int cap=capacity.get(e);
   170 	      int flo=flow.get(e);
   171 	      int remcap=cap-flo;
   172 	      
   173 	      if ( remcap >= exc ) {       
   174 		/*A nonsaturating push.*/
   175 		
   176 		flow.set(e, flo+exc);
   177 		excess.set(v, excess.get(v)+exc);
   178 		exc=0;
   179 		break; 
   180 		
   181 	      } else { 
   182 		/*A saturating push.*/
   183 		
   184 		flow.set(e, cap );
   185 		excess.set(v, excess.get(v)+remcap);
   186 		exc-=remcap;
   187 	      }
   188 	    } else if ( newlevel > level.get(v) ) newlevel = level.get(v);
   189 	    
   190 	  } //for out edges wv 
   191 	
   192 	
   193 	if ( exc > 0 ) {	
   194 	  for( InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) {
   195 	    
   196 	    if( flow.get(e) == 0 ) continue; 
   197 	    NodeIt v=G.tail(e);  
   198 	    //e=vw
   199 	    
   200 	    if( lev > level.get(v) ) {  
   201 	      /*Push is allowed now*/
   202 	      
   203 	      if ( excess.get(v)==0 && v != s && v !=t) 
   204 		stack[level.get(v)].push(v); 
   205 	      /*v becomes active.*/
   206 	      
   207 	      int flo=flow.get(e);
   208 	      
   209 	      if ( flo >= exc ) { 
   210 		/*A nonsaturating push.*/
   211 		
   212 		flow.set(e, flo-exc);
   213 		excess.set(v, excess.get(v)+exc);
   214 		exc=0;
   215 		break; 
   216 	      } else {                                               
   217 		/*A saturating push.*/
   218 		
   219 		excess.set(v, excess.get(v)+flo);
   220 		exc-=flo;
   221 		flow.set(e,0);
   222 	      }  
   223 	    } else if ( newlevel > level.get(v) ) newlevel = level.get(v);
   224 	    
   225 	  } //for in edges vw
   226 	  
   227 	} // if w still has excess after the out edge for cycle
   228 	
   229 	excess.set(w, exc);
   230 	
   231 
   232 	
   233 
   234 	/*
   235 	  Relabel
   236 	*/
   237 	
   238 	if ( exc > 0 ) {
   239 	  //now 'lev' is the old level of w
   240 	  level.set(w,++newlevel);
   241 	  
   242 	  if ( lev < n ) {
   243 	    --numb[lev];
   244 	    
   245 	    if ( !numb[lev] && lev < A*n ) {  //If the level of w gets empty. 
   246 	      
   247 	      for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid() ; ++v) {
   248 		if (level.get(v) > lev && level.get(v) < n ) level.set(v,n);  
   249 	      }
   250 	      for (int i=lev+1 ; i!=n ; ++i) numb[i]=0; 
   251 	      if ( newlevel < n ) newlevel=n; 
   252 	    } else { 
   253 	      if ( newlevel < n ) ++numb[newlevel]; 
   254 	    }
   255 	  } 
   256 	  
   257 	  stack[newlevel].push(w);
   258 	  b=newlevel;
   259 	  
   260 	}
   261 	
   262       } // while(b)
   263       
   264       
   265       value = excess.get(t);
   266       /*Max flow value.*/
   267 
   268 
   269     } //void run()
   270 
   271 
   272 
   273 
   274 
   275     /*
   276       Returns the maximum value of a flow.
   277      */
   278 
   279     T maxflow() {
   280       return value;
   281     }
   282 
   283 
   284 
   285     /*
   286       For the maximum flow x found by the algorithm, it returns the flow value on Edge e, i.e. x(e). 
   287     */
   288 
   289     T flowonedge(const EdgeIt e) {
   290       return flow.get(e);
   291     }
   292 
   293 
   294 
   295     /*
   296       Returns the maximum flow x found by the algorithm.
   297     */
   298 
   299     FlowMap allflow() {
   300       return flow;
   301     }
   302 
   303 
   304 
   305 
   306     /*
   307       Returns the minimum min cut, by a bfs from s in the residual graph.
   308     */
   309     
   310     template<typename CutMap>
   311     void mincut(CutMap& M) {
   312     
   313       std::queue<NodeIt> queue;
   314       
   315       M.set(s,true);      
   316       queue.push(s);
   317 
   318       while (!queue.empty()) {
   319         NodeIt w=queue.front();
   320 	queue.pop();
   321 	
   322 	for(OutEdgeIt e=G.template first<OutEdgeIt>(w) ; e.valid(); ++e) {
   323 	  NodeIt v=G.head(e);
   324 	  if (!M.get(v) && flow.get(e) < capacity.get(e) ) {
   325 	    queue.push(v);
   326 	    M.set(v, true);
   327 	  }
   328 	} 
   329 
   330 	for(InEdgeIt e=G.template first<InEdgeIt>(w) ; e.valid(); ++e) {
   331 	  NodeIt v=G.tail(e);
   332 	  if (!M.get(v) && flow.get(e) > 0 ) {
   333 	    queue.push(v);
   334 	    M.set(v, true);
   335 	  }
   336 	}
   337 
   338       }
   339     }
   340 
   341 
   342 
   343     /*
   344       Returns the maximum min cut, by a reverse bfs 
   345       from t in the residual graph.
   346     */
   347     
   348     template<typename CutMap>
   349     void max_mincut(CutMap& M) {
   350     
   351       std::queue<NodeIt> queue;
   352       
   353       M.set(t,true);        
   354       queue.push(t);
   355 
   356       while (!queue.empty()) {
   357         NodeIt w=queue.front();
   358 	queue.pop();
   359 
   360 	for(InEdgeIt e=G.template first<InEdgeIt>(w) ; e.valid(); ++e) {
   361 	  NodeIt v=G.tail(e);
   362 	  if (!M.get(v) && flow.get(e) < capacity.get(e) ) {
   363 	    queue.push(v);
   364 	    M.set(v, true);
   365 	  }
   366 	}
   367 
   368 	for(OutEdgeIt e=G.template first<OutEdgeIt>(w) ; e.valid(); ++e) {
   369 	  NodeIt v=G.head(e);
   370 	  if (!M.get(v) && flow.get(e) > 0 ) {
   371 	    queue.push(v);
   372 	    M.set(v, true);
   373 	  }
   374 	}
   375       }
   376 
   377       for(EachNodeIt v=G.template first<EachNodeIt>() ; v.valid(); ++v) {
   378 	M.set(v, !M.get(v));
   379       }
   380 
   381     }
   382 
   383 
   384 
   385     template<typename CutMap>
   386     void min_mincut(CutMap& M) {
   387       mincut(M);
   388     }
   389 
   390 
   391 
   392   };
   393 }//namespace hugo
   394 #endif 
   395 
   396 
   397 
   398