src/work/jacint/preflowproba.h
author marci
Fri, 23 Apr 2004 08:08:41 +0000
changeset 380 6399494e30b1
parent 374 0fc9cd9b854a
permissions -rw-r--r--
.
jacint@370
     1
// -*- C++ -*-
jacint@370
     2
jacint@370
     3
//run gyorsan tudna adni a minmincutot a 2 fazis elejen , ne vegyuk be konstruktorba egy cutmapet?
jacint@370
     4
//constzero jo igy?
jacint@370
     5
jacint@370
     6
//majd marci megmondja betegyem-e bfs-t meg resgraphot
jacint@370
     7
jacint@370
     8
/*
jacint@370
     9
Heuristics: 
jacint@370
    10
 2 phase
jacint@370
    11
 gap
jacint@370
    12
 list 'level_list' on the nodes on level i implemented by hand
jacint@370
    13
 stack 'active' on the active nodes on level i implemented by hand
jacint@370
    14
 runs heuristic 'highest label' for H1*n relabels
jacint@370
    15
 runs heuristic 'bound decrease' for H0*n relabels, starts with 'highest label'
jacint@370
    16
 
jacint@370
    17
Parameters H0 and H1 are initialized to 20 and 10.
jacint@370
    18
jacint@370
    19
Constructors:
jacint@370
    20
jacint@370
    21
Preflow(Graph, Node, Node, CapMap, FlowMap, bool) : bool must be false if 
jacint@370
    22
     FlowMap is not constant zero, and should be true if it is
jacint@370
    23
jacint@370
    24
Members:
jacint@370
    25
jacint@370
    26
void run()
jacint@370
    27
jacint@370
    28
T flowValue() : returns the value of a maximum flow
jacint@370
    29
jacint@370
    30
void minMinCut(CutMap& M) : sets M to the characteristic vector of the 
jacint@370
    31
     minimum min cut. M should be a map of bools initialized to false.
jacint@370
    32
jacint@370
    33
void maxMinCut(CutMap& M) : sets M to the characteristic vector of the 
jacint@370
    34
     maximum min cut. M should be a map of bools initialized to false.
jacint@370
    35
jacint@370
    36
void minCut(CutMap& M) : sets M to the characteristic vector of 
jacint@370
    37
     a min cut. M should be a map of bools initialized to false.
jacint@370
    38
jacint@370
    39
FIXME reset
jacint@370
    40
jacint@370
    41
*/
jacint@370
    42
marci@376
    43
#ifndef HUGO_PREFLOW_PROBA_H
marci@376
    44
#define HUGO_PREFLOW_PROBA_H
jacint@370
    45
jacint@370
    46
#define H0 20
jacint@370
    47
#define H1 1
jacint@370
    48
jacint@370
    49
#include <vector>
jacint@370
    50
#include <queue>
jacint@370
    51
#include<graph_wrapper.h>
jacint@370
    52
jacint@370
    53
namespace hugo {
jacint@370
    54
jacint@370
    55
  template <typename Graph, typename T, 
jacint@370
    56
	    typename CapMap=typename Graph::EdgeMap<T>, 
jacint@370
    57
            typename FlowMap=typename Graph::EdgeMap<T> >
marci@376
    58
  class PreflowProba {
jacint@370
    59
    
jacint@370
    60
    typedef typename Graph::Node Node;
jacint@370
    61
    typedef typename Graph::Edge Edge;
jacint@370
    62
    typedef typename Graph::NodeIt NodeIt;
jacint@370
    63
    typedef typename Graph::OutEdgeIt OutEdgeIt;
jacint@370
    64
    typedef typename Graph::InEdgeIt InEdgeIt;
jacint@370
    65
    
jacint@370
    66
    const Graph& G;
jacint@370
    67
    Node s;
jacint@370
    68
    Node t;
jacint@370
    69
    const CapMap& capacity;  
jacint@370
    70
    FlowMap& flow;
jacint@370
    71
    T value;
jacint@370
    72
    bool constzero;
jacint@374
    73
    bool res;
jacint@370
    74
jacint@370
    75
    typedef ResGraphWrapper<const Graph, T, CapMap, FlowMap> ResGW;
jacint@370
    76
    typedef typename ResGW::OutEdgeIt ResOutEdgeIt;
jacint@370
    77
    typedef typename ResGW::InEdgeIt ResInEdgeIt;
jacint@370
    78
    typedef typename ResGW::Edge ResEdge;
jacint@370
    79
 
jacint@370
    80
  public:
marci@376
    81
    PreflowProba(Graph& _G, Node _s, Node _t, CapMap& _capacity, 
jacint@374
    82
	    FlowMap& _flow, bool _constzero, bool _res ) :
jacint@374
    83
      G(_G), s(_s), t(_t), capacity(_capacity), flow(_flow), constzero(_constzero), res(_res) {}
jacint@370
    84
    
jacint@370
    85
    
jacint@370
    86
    void run() {
jacint@370
    87
jacint@370
    88
      ResGW res_graph(G, capacity, flow);
jacint@370
    89
jacint@370
    90
      value=0;                //for the subsequent runs
jacint@370
    91
jacint@370
    92
      bool phase=0;        //phase 0 is the 1st phase, phase 1 is the 2nd
jacint@370
    93
      int n=G.nodeNum(); 
jacint@370
    94
      int heur0=(int)(H0*n);  //time while running 'bound decrease' 
jacint@370
    95
      int heur1=(int)(H1*n);  //time while running 'highest label'
jacint@370
    96
      int heur=heur1;         //starting time interval (#of relabels)
jacint@370
    97
      bool what_heur=1;       
jacint@370
    98
      /*
jacint@370
    99
	what_heur is 0 in case 'bound decrease' 
jacint@370
   100
	and 1 in case 'highest label'
jacint@370
   101
      */
jacint@370
   102
      bool end=false;     
jacint@370
   103
      /*
jacint@370
   104
	Needed for 'bound decrease', 'true'
jacint@370
   105
	means no active nodes are above bound b.
jacint@370
   106
      */
jacint@370
   107
      int relabel=0;
jacint@370
   108
      int k=n-2;  //bound on the highest level under n containing a node
jacint@370
   109
      int b=k;    //bound on the highest level under n of an active node
jacint@370
   110
      
jacint@370
   111
      typename Graph::NodeMap<int> level(G,n);      
jacint@370
   112
      typename Graph::NodeMap<T> excess(G); 
jacint@370
   113
jacint@370
   114
      std::vector<Node> active(n-1,INVALID);
jacint@370
   115
      typename Graph::NodeMap<Node> next(G,INVALID);
jacint@370
   116
      //Stack of the active nodes in level i < n.
jacint@370
   117
      //We use it in both phases.
jacint@370
   118
jacint@370
   119
      typename Graph::NodeMap<Node> left(G,INVALID);
jacint@370
   120
      typename Graph::NodeMap<Node> right(G,INVALID);
jacint@370
   121
      std::vector<Node> level_list(n,INVALID);
jacint@370
   122
      /*
jacint@370
   123
	List of the nodes in level i<n.
jacint@370
   124
      */
jacint@370
   125
jacint@370
   126
jacint@370
   127
      if ( constzero ) {
jacint@370
   128
     
jacint@370
   129
	/*Reverse_bfs from t, to find the starting level.*/
jacint@370
   130
	level.set(t,0);
jacint@370
   131
	std::queue<Node> bfs_queue;
jacint@370
   132
	bfs_queue.push(t);
jacint@370
   133
	
jacint@370
   134
	while (!bfs_queue.empty()) {
jacint@370
   135
	  
jacint@370
   136
	  Node v=bfs_queue.front();	
jacint@370
   137
	  bfs_queue.pop();
jacint@370
   138
	  int l=level[v]+1;
jacint@370
   139
	  
jacint@370
   140
	  InEdgeIt e;
jacint@370
   141
	  for(G.first(e,v); G.valid(e); G.next(e)) {
jacint@370
   142
	    Node w=G.tail(e);
jacint@370
   143
	    if ( level[w] == n && w != s ) {
jacint@370
   144
	      bfs_queue.push(w);
jacint@370
   145
	      Node first=level_list[l];
jacint@370
   146
	      if ( G.valid(first) ) left.set(first,w);
jacint@370
   147
	      right.set(w,first);
jacint@370
   148
	      level_list[l]=w;
jacint@370
   149
	      level.set(w, l);
jacint@370
   150
	    }
jacint@370
   151
	  }
jacint@370
   152
	}
jacint@370
   153
jacint@370
   154
	//the starting flow
jacint@370
   155
	OutEdgeIt e;
jacint@370
   156
	for(G.first(e,s); G.valid(e); G.next(e)) 
jacint@370
   157
	{
jacint@370
   158
	  T c=capacity[e];
jacint@370
   159
	  if ( c == 0 ) continue;
jacint@370
   160
	  Node w=G.head(e);
jacint@370
   161
	  if ( level[w] < n ) {	  
jacint@370
   162
	    if ( excess[w] == 0 && w!=t ) {
jacint@370
   163
	      next.set(w,active[level[w]]);
jacint@370
   164
	      active[level[w]]=w;
jacint@370
   165
	    }
jacint@370
   166
	    flow.set(e, c); 
jacint@370
   167
	    excess.set(w, excess[w]+c);
jacint@370
   168
	  }
jacint@370
   169
	}
jacint@370
   170
      }
jacint@370
   171
      else 
jacint@370
   172
      {
jacint@370
   173
	
jacint@370
   174
	/*
jacint@370
   175
	  Reverse_bfs from t in the residual graph, 
jacint@370
   176
	  to find the starting level.
jacint@370
   177
	*/
jacint@370
   178
	level.set(t,0);
jacint@370
   179
	std::queue<Node> bfs_queue;
jacint@370
   180
	bfs_queue.push(t);
jacint@370
   181
	
jacint@370
   182
	while (!bfs_queue.empty()) {
jacint@370
   183
	  
jacint@370
   184
	  Node v=bfs_queue.front();	
jacint@370
   185
	  bfs_queue.pop();
jacint@370
   186
	  int l=level[v]+1;
jacint@370
   187
	  
jacint@370
   188
	  InEdgeIt e;
jacint@370
   189
	  for(G.first(e,v); G.valid(e); G.next(e)) {
jacint@370
   190
	    if ( capacity[e] == flow[e] ) continue;
jacint@370
   191
	    Node w=G.tail(e);
jacint@370
   192
	    if ( level[w] == n && w != s ) {
jacint@370
   193
	      bfs_queue.push(w);
jacint@370
   194
	      Node first=level_list[l];
jacint@370
   195
	      if ( G.valid(first) ) left.set(first,w);
jacint@370
   196
	      right.set(w,first);
jacint@370
   197
	      level_list[l]=w;
jacint@370
   198
	      level.set(w, l);
jacint@370
   199
	    }
jacint@370
   200
	  }
jacint@370
   201
	    
jacint@370
   202
	  OutEdgeIt f;
jacint@370
   203
	  for(G.first(f,v); G.valid(f); G.next(f)) {
jacint@370
   204
	    if ( 0 == flow[f] ) continue;
jacint@370
   205
	    Node w=G.head(f);
jacint@370
   206
	    if ( level[w] == n && w != s ) {
jacint@370
   207
	      bfs_queue.push(w);
jacint@370
   208
	      Node first=level_list[l];
jacint@370
   209
	      if ( G.valid(first) ) left.set(first,w);
jacint@370
   210
	      right.set(w,first);
jacint@370
   211
	      level_list[l]=w;
jacint@370
   212
	      level.set(w, l);
jacint@370
   213
	    }
jacint@370
   214
	    }
jacint@370
   215
	}
jacint@370
   216
      
jacint@370
   217
	
jacint@370
   218
	/*
jacint@370
   219
	  Counting the excess
jacint@370
   220
	*/
jacint@370
   221
	NodeIt v;
jacint@370
   222
	for(G.first(v); G.valid(v); G.next(v)) {
jacint@370
   223
	  T exc=0;
jacint@370
   224
jacint@370
   225
	  InEdgeIt e;
jacint@370
   226
	  for(G.first(e,v); G.valid(e); G.next(e)) exc+=flow[e];
jacint@370
   227
	  OutEdgeIt f;
jacint@370
   228
	  for(G.first(f,v); G.valid(f); G.next(f)) exc-=flow[e];
jacint@370
   229
jacint@370
   230
	  excess.set(v,exc);	  
jacint@370
   231
jacint@370
   232
	  //putting the active nodes into the stack
jacint@370
   233
	  int lev=level[v];
jacint@370
   234
	  if ( exc > 0 && lev < n ) {
jacint@370
   235
	    next.set(v,active[lev]);
jacint@370
   236
	    active[lev]=v;
jacint@370
   237
	  }
jacint@370
   238
	}
jacint@370
   239
jacint@370
   240
jacint@370
   241
	//the starting flow
jacint@370
   242
	OutEdgeIt e;
jacint@370
   243
	for(G.first(e,s); G.valid(e); G.next(e)) 
jacint@370
   244
	{
jacint@370
   245
	  T rem=capacity[e]-flow[e];
jacint@370
   246
	  if ( rem == 0 ) continue;
jacint@370
   247
	  Node w=G.head(e);
jacint@370
   248
	  if ( level[w] < n ) {	  
jacint@370
   249
	    if ( excess[w] == 0 && w!=t ) {
jacint@370
   250
	      next.set(w,active[level[w]]);
jacint@370
   251
	      active[level[w]]=w;
jacint@370
   252
	    }
jacint@370
   253
	    flow.set(e, capacity[e]); 
jacint@370
   254
	    excess.set(w, excess[w]+rem);
jacint@370
   255
	  }
jacint@370
   256
	}
jacint@370
   257
	
jacint@370
   258
	InEdgeIt f;
jacint@370
   259
	for(G.first(f,s); G.valid(f); G.next(f)) 
jacint@370
   260
	{
jacint@370
   261
	  if ( flow[f] == 0 ) continue;
jacint@370
   262
	  Node w=G.head(f);
jacint@370
   263
	  if ( level[w] < n ) {	  
jacint@370
   264
	    if ( excess[w] == 0 && w!=t ) {
jacint@370
   265
	      next.set(w,active[level[w]]);
jacint@370
   266
	      active[level[w]]=w;
jacint@370
   267
	    }
jacint@370
   268
	    excess.set(w, excess[w]+flow[f]);
jacint@370
   269
	    flow.set(f, 0); 
jacint@370
   270
	  }
jacint@370
   271
	}
jacint@370
   272
      }
jacint@370
   273
jacint@370
   274
jacint@370
   275
jacint@370
   276
jacint@370
   277
      /* 
jacint@370
   278
	 End of preprocessing 
jacint@370
   279
      */
jacint@370
   280
jacint@370
   281
jacint@370
   282
jacint@370
   283
      /*
jacint@370
   284
	Push/relabel on the highest level active nodes.
jacint@370
   285
      */	
jacint@370
   286
      while ( true ) {
jacint@370
   287
	
jacint@370
   288
	if ( b == 0 ) {
jacint@370
   289
	  if ( phase ) break;
jacint@370
   290
	  
jacint@370
   291
	  if ( !what_heur && !end && k > 0 ) {
jacint@370
   292
	    b=k;
jacint@370
   293
	    end=true;
jacint@370
   294
	  } else {
jacint@370
   295
	    phase=1;
jacint@370
   296
	    level.set(s,0);
jacint@370
   297
	    std::queue<Node> bfs_queue;
jacint@370
   298
	    bfs_queue.push(s);
jacint@370
   299
	    
jacint@370
   300
	    while (!bfs_queue.empty()) {
jacint@370
   301
	      
jacint@370
   302
	      Node v=bfs_queue.front();	
jacint@370
   303
	      bfs_queue.pop();
jacint@370
   304
	      int l=level[v]+1;
jacint@370
   305
	      
jacint@374
   306
	      if (res){	
jacint@374
   307
		ResInEdgeIt e;
jacint@374
   308
		for(res_graph.first(e,v); res_graph.valid(e); 
jacint@374
   309
		    res_graph.next(e)) {
jacint@374
   310
		  Node u=res_graph.tail(e);
jacint@374
   311
		  if ( level[u] >= n ) { 
jacint@374
   312
		    bfs_queue.push(u);
jacint@374
   313
		    level.set(u, l);
jacint@374
   314
		    if ( excess[u] > 0 ) {
jacint@374
   315
		      next.set(u,active[l]);
jacint@374
   316
		      active[l]=u;
jacint@374
   317
		    }
jacint@370
   318
		  }
jacint@370
   319
		}
jacint@374
   320
	      } else {
jacint@374
   321
	        InEdgeIt e;
jacint@374
   322
		for(G.first(e,v); G.valid(e); G.next(e)) {
jacint@374
   323
		  if ( capacity[e] == flow[e] ) continue;
jacint@374
   324
		  Node u=G.tail(e);
jacint@374
   325
		  if ( level[u] >= n ) { 
jacint@374
   326
		    bfs_queue.push(u);
jacint@374
   327
		    level.set(u, l);
jacint@374
   328
		    if ( excess[u] > 0 ) {
jacint@374
   329
		      next.set(u,active[l]);
jacint@374
   330
		      active[l]=u;
jacint@374
   331
		    }
jacint@374
   332
		  }
jacint@374
   333
		}
jacint@374
   334
		
jacint@374
   335
		OutEdgeIt f;
jacint@374
   336
		for(G.first(f,v); G.valid(f); G.next(f)) {
jacint@374
   337
		  if ( 0 == flow[f] ) continue;
jacint@374
   338
		  Node u=G.head(f);
jacint@374
   339
		  if ( level[u] >= n ) { 
jacint@374
   340
		    bfs_queue.push(u);
jacint@374
   341
		    level.set(u, l);
jacint@374
   342
		    if ( excess[u] > 0 ) {
jacint@374
   343
		      next.set(u,active[l]);
jacint@374
   344
		      active[l]=u;
jacint@374
   345
		    }
jacint@370
   346
		  }
jacint@370
   347
		}
jacint@370
   348
	      }
jacint@370
   349
	    }
jacint@370
   350
	    b=n-2;
jacint@370
   351
	    }
jacint@370
   352
	    
jacint@370
   353
	}
jacint@370
   354
	  
jacint@370
   355
	  
jacint@370
   356
	if ( !G.valid(active[b]) ) --b; 
jacint@370
   357
	else {
jacint@370
   358
	  end=false;  
jacint@370
   359
jacint@370
   360
	  Node w=active[b];
jacint@370
   361
	  active[b]=next[w];
jacint@370
   362
	  int lev=level[w];
jacint@370
   363
	  T exc=excess[w];
jacint@370
   364
	  int newlevel=n;       //bound on the next level of w
jacint@370
   365
	  
jacint@370
   366
	  OutEdgeIt e;
jacint@370
   367
	  for(G.first(e,w); G.valid(e); G.next(e)) {
jacint@370
   368
	    
jacint@370
   369
	    if ( flow[e] == capacity[e] ) continue; 
jacint@370
   370
	    Node v=G.head(e);            
jacint@370
   371
	    //e=wv	    
jacint@370
   372
	    
jacint@370
   373
	    if( lev > level[v] ) {      
jacint@370
   374
	      /*Push is allowed now*/
jacint@370
   375
	      
jacint@370
   376
	      if ( excess[v]==0 && v!=t && v!=s ) {
jacint@370
   377
		int lev_v=level[v];
jacint@370
   378
		next.set(v,active[lev_v]);
jacint@370
   379
		active[lev_v]=v;
jacint@370
   380
	      }
jacint@370
   381
	      
jacint@370
   382
	      T cap=capacity[e];
jacint@370
   383
	      T flo=flow[e];
jacint@370
   384
	      T remcap=cap-flo;
jacint@370
   385
	      
jacint@370
   386
	      if ( remcap >= exc ) {       
jacint@370
   387
		/*A nonsaturating push.*/
jacint@370
   388
		
jacint@370
   389
		flow.set(e, flo+exc);
jacint@370
   390
		excess.set(v, excess[v]+exc);
jacint@370
   391
		exc=0;
jacint@370
   392
		break; 
jacint@370
   393
		
jacint@370
   394
	      } else { 
jacint@370
   395
		/*A saturating push.*/
jacint@370
   396
		
jacint@370
   397
		flow.set(e, cap);
jacint@370
   398
		excess.set(v, excess[v]+remcap);
jacint@370
   399
		exc-=remcap;
jacint@370
   400
	      }
jacint@370
   401
	    } else if ( newlevel > level[v] ){
jacint@370
   402
	      newlevel = level[v];
jacint@370
   403
	    }	    
jacint@370
   404
	    
jacint@370
   405
	  } //for out edges wv 
jacint@370
   406
	
jacint@370
   407
	
jacint@370
   408
	if ( exc > 0 ) {	
jacint@370
   409
	  InEdgeIt e;
jacint@370
   410
	  for(G.first(e,w); G.valid(e); G.next(e)) {
jacint@370
   411
	    
jacint@370
   412
	    if( flow[e] == 0 ) continue; 
jacint@370
   413
	    Node v=G.tail(e);  
jacint@370
   414
	    //e=vw
jacint@370
   415
	    
jacint@370
   416
	    if( lev > level[v] ) {  
jacint@370
   417
	      /*Push is allowed now*/
jacint@370
   418
	      
jacint@370
   419
	      if ( excess[v]==0 && v!=t && v!=s ) {
jacint@370
   420
		int lev_v=level[v];
jacint@370
   421
		next.set(v,active[lev_v]);
jacint@370
   422
		active[lev_v]=v;
jacint@370
   423
	      }
jacint@370
   424
	      
jacint@370
   425
	      T flo=flow[e];
jacint@370
   426
	      
jacint@370
   427
	      if ( flo >= exc ) { 
jacint@370
   428
		/*A nonsaturating push.*/
jacint@370
   429
		
jacint@370
   430
		flow.set(e, flo-exc);
jacint@370
   431
		excess.set(v, excess[v]+exc);
jacint@370
   432
		exc=0;
jacint@370
   433
		break; 
jacint@370
   434
	      } else {                                               
jacint@370
   435
		/*A saturating push.*/
jacint@370
   436
		
jacint@370
   437
		excess.set(v, excess[v]+flo);
jacint@370
   438
		exc-=flo;
jacint@370
   439
		flow.set(e,0);
jacint@370
   440
	      }  
jacint@370
   441
	    } else if ( newlevel > level[v] ) {
jacint@370
   442
	      newlevel = level[v];
jacint@370
   443
	    }	    
jacint@370
   444
	  } //for in edges vw
jacint@370
   445
	  
jacint@370
   446
	} // if w still has excess after the out edge for cycle
jacint@370
   447
	
jacint@370
   448
	excess.set(w, exc);
jacint@370
   449
	 
jacint@370
   450
	/*
jacint@370
   451
	  Relabel
jacint@370
   452
	*/
jacint@370
   453
	
jacint@370
   454
jacint@370
   455
	if ( exc > 0 ) {
jacint@370
   456
	  //now 'lev' is the old level of w
jacint@370
   457
	
jacint@370
   458
	  if ( phase ) {
jacint@370
   459
	    level.set(w,++newlevel);
jacint@370
   460
	    next.set(w,active[newlevel]);
jacint@370
   461
	    active[newlevel]=w;
jacint@370
   462
	    b=newlevel;
jacint@370
   463
	  } else {
jacint@370
   464
	    //unlacing starts
jacint@370
   465
	    Node right_n=right[w];
jacint@370
   466
	    Node left_n=left[w];
jacint@370
   467
jacint@370
   468
	    if ( G.valid(right_n) ) {
jacint@370
   469
	      if ( G.valid(left_n) ) {
jacint@370
   470
		right.set(left_n, right_n);
jacint@370
   471
		left.set(right_n, left_n);
jacint@370
   472
	      } else {
jacint@370
   473
		level_list[lev]=right_n;   
jacint@370
   474
		left.set(right_n, INVALID);
jacint@370
   475
	      } 
jacint@370
   476
	    } else {
jacint@370
   477
	      if ( G.valid(left_n) ) {
jacint@370
   478
		right.set(left_n, INVALID);
jacint@370
   479
	      } else { 
jacint@370
   480
		level_list[lev]=INVALID;   
jacint@370
   481
	      } 
jacint@370
   482
	    } 
jacint@370
   483
	    //unlacing ends
jacint@370
   484
		
jacint@370
   485
	    if ( !G.valid(level_list[lev]) ) {
jacint@370
   486
	      
jacint@370
   487
	       //gapping starts
jacint@370
   488
	      for (int i=lev; i!=k ; ) {
jacint@370
   489
		Node v=level_list[++i];
jacint@370
   490
		while ( G.valid(v) ) {
jacint@370
   491
		  level.set(v,n);
jacint@370
   492
		  v=right[v];
jacint@370
   493
		}
jacint@370
   494
		level_list[i]=INVALID;
jacint@370
   495
		if ( !what_heur ) active[i]=INVALID;
jacint@370
   496
	      }	     
jacint@370
   497
jacint@370
   498
	      level.set(w,n);
jacint@370
   499
	      b=lev-1;
jacint@370
   500
	      k=b;
jacint@370
   501
	      //gapping ends
jacint@370
   502
	    
jacint@370
   503
	    } else {
jacint@370
   504
	      
jacint@370
   505
	      if ( newlevel == n ) level.set(w,n); 
jacint@370
   506
	      else {
jacint@370
   507
		level.set(w,++newlevel);
jacint@370
   508
		next.set(w,active[newlevel]);
jacint@370
   509
		active[newlevel]=w;
jacint@370
   510
		if ( what_heur ) b=newlevel;
jacint@370
   511
		if ( k < newlevel ) ++k;      //now k=newlevel
jacint@370
   512
		Node first=level_list[newlevel];
jacint@370
   513
		if ( G.valid(first) ) left.set(first,w);
jacint@370
   514
		right.set(w,first);
jacint@370
   515
		left.set(w,INVALID);
jacint@370
   516
		level_list[newlevel]=w;
jacint@370
   517
	      }
jacint@370
   518
	    }
jacint@370
   519
jacint@370
   520
jacint@370
   521
	    ++relabel; 
jacint@370
   522
	    if ( relabel >= heur ) {
jacint@370
   523
	      relabel=0;
jacint@370
   524
	      if ( what_heur ) {
jacint@370
   525
		what_heur=0;
jacint@370
   526
		heur=heur0;
jacint@370
   527
		end=false;
jacint@370
   528
	      } else {
jacint@370
   529
		what_heur=1;
jacint@370
   530
		heur=heur1;
jacint@370
   531
		b=k; 
jacint@370
   532
	      }
jacint@370
   533
	    }
jacint@370
   534
	  } //phase 0
jacint@370
   535
	  
jacint@370
   536
	  
jacint@370
   537
	} // if ( exc > 0 )
jacint@370
   538
	  
jacint@370
   539
	
jacint@370
   540
	}  // if stack[b] is nonempty
jacint@370
   541
	
jacint@370
   542
      } // while(true)
jacint@370
   543
jacint@370
   544
jacint@370
   545
      value = excess[t];
jacint@370
   546
      /*Max flow value.*/
jacint@370
   547
     
jacint@370
   548
    } //void run()
jacint@370
   549
jacint@370
   550
jacint@370
   551
jacint@370
   552
jacint@370
   553
jacint@370
   554
    /*
jacint@370
   555
      Returns the maximum value of a flow.
jacint@370
   556
     */
jacint@370
   557
jacint@370
   558
    T flowValue() {
jacint@370
   559
      return value;
jacint@370
   560
    }
jacint@370
   561
jacint@370
   562
jacint@370
   563
    FlowMap Flow() {
jacint@370
   564
      return flow;
jacint@370
   565
      }
jacint@370
   566
jacint@370
   567
jacint@370
   568
    
jacint@370
   569
    void Flow(FlowMap& _flow ) {
jacint@370
   570
      NodeIt v;
jacint@370
   571
      for(G.first(v) ; G.valid(v); G.next(v))
jacint@370
   572
	_flow.set(v,flow[v]);
jacint@370
   573
    }
jacint@370
   574
jacint@370
   575
jacint@370
   576
jacint@370
   577
    /*
jacint@370
   578
      Returns the minimum min cut, by a bfs from s in the residual graph.
jacint@370
   579
    */
jacint@370
   580
   
jacint@370
   581
    template<typename _CutMap>
jacint@370
   582
    void minMinCut(_CutMap& M) {
jacint@370
   583
    
jacint@370
   584
      std::queue<Node> queue;
jacint@370
   585
      
jacint@370
   586
      M.set(s,true);      
jacint@370
   587
      queue.push(s);
jacint@370
   588
jacint@370
   589
      while (!queue.empty()) {
jacint@370
   590
        Node w=queue.front();
jacint@370
   591
	queue.pop();
jacint@370
   592
jacint@370
   593
	OutEdgeIt e;
jacint@370
   594
	for(G.first(e,w) ; G.valid(e); G.next(e)) {
jacint@370
   595
	  Node v=G.head(e);
jacint@370
   596
	  if (!M[v] && flow[e] < capacity[e] ) {
jacint@370
   597
	    queue.push(v);
jacint@370
   598
	    M.set(v, true);
jacint@370
   599
	  }
jacint@370
   600
	} 
jacint@370
   601
jacint@370
   602
	InEdgeIt f;
jacint@370
   603
	for(G.first(f,w) ; G.valid(f); G.next(f)) {
jacint@370
   604
	  Node v=G.tail(f);
jacint@370
   605
	  if (!M[v] && flow[f] > 0 ) {
jacint@370
   606
	    queue.push(v);
jacint@370
   607
	    M.set(v, true);
jacint@370
   608
	  }
jacint@370
   609
	} 
jacint@370
   610
      }
jacint@370
   611
    }
jacint@370
   612
jacint@370
   613
jacint@370
   614
  
jacint@370
   615
    /*
jacint@370
   616
      Returns the maximum min cut, by a reverse bfs 
jacint@370
   617
      from t in the residual graph.
jacint@370
   618
    */
jacint@370
   619
    
jacint@370
   620
    template<typename _CutMap>
jacint@370
   621
    void maxMinCut(_CutMap& M) {
jacint@370
   622
    
jacint@370
   623
      std::queue<Node> queue;
jacint@370
   624
      
jacint@370
   625
      M.set(t,true);        
jacint@370
   626
      queue.push(t);
jacint@370
   627
jacint@370
   628
      while (!queue.empty()) {
jacint@370
   629
        Node w=queue.front();
jacint@370
   630
	queue.pop();
jacint@370
   631
jacint@370
   632
jacint@370
   633
	InEdgeIt e;
jacint@370
   634
	for(G.first(e,w) ; G.valid(e); G.next(e)) {
jacint@370
   635
	  Node v=G.tail(e);
jacint@370
   636
	  if (!M[v] && flow[e] < capacity[e] ) {
jacint@370
   637
	    queue.push(v);
jacint@370
   638
	    M.set(v, true);
jacint@370
   639
	  }
jacint@370
   640
	}
jacint@370
   641
	
jacint@370
   642
	OutEdgeIt f;
jacint@370
   643
	for(G.first(f,w) ; G.valid(f); G.next(f)) {
jacint@370
   644
	  Node v=G.head(f);
jacint@370
   645
	  if (!M[v] && flow[f] > 0 ) {
jacint@370
   646
	    queue.push(v);
jacint@370
   647
	    M.set(v, true);
jacint@370
   648
	  }
jacint@370
   649
	}
jacint@370
   650
      }
jacint@370
   651
jacint@370
   652
      NodeIt v;
jacint@370
   653
      for(G.first(v) ; G.valid(v); G.next(v)) {
jacint@370
   654
	M.set(v, !M[v]);
jacint@370
   655
      }
jacint@370
   656
jacint@370
   657
    }
jacint@370
   658
jacint@370
   659
jacint@370
   660
jacint@370
   661
    template<typename CutMap>
jacint@370
   662
    void minCut(CutMap& M) {
jacint@370
   663
      minMinCut(M);
jacint@370
   664
    }
jacint@370
   665
jacint@370
   666
    
jacint@370
   667
    void reset_target (Node _t) {t=_t;}
jacint@370
   668
    void reset_source (Node _s) {s=_s;}
jacint@370
   669
   
jacint@370
   670
    template<typename _CapMap>   
jacint@370
   671
    void reset_cap (_CapMap _cap) {capacity=_cap;}
jacint@370
   672
jacint@370
   673
    template<typename _FlowMap>   
jacint@370
   674
    void reset_cap (_FlowMap _flow, bool _constzero) {
jacint@370
   675
      flow=_flow;
jacint@370
   676
      constzero=_constzero;
jacint@370
   677
    }
jacint@370
   678
jacint@370
   679
jacint@370
   680
jacint@370
   681
  };
jacint@370
   682
jacint@370
   683
} //namespace hugo
jacint@370
   684
marci@376
   685
#endif //PREFLOW_PROBA_H
jacint@370
   686
jacint@370
   687
jacint@370
   688
jacint@370
   689