1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/src/work/jacint/preflow_excess.h Tue Apr 27 10:27:34 2004 +0000
1.3 @@ -0,0 +1,664 @@
1.4 +// -*- C++ -*-
1.5 +
1.6 +//run gyorsan tudna adni a minmincutot a 2 fazis elejen , ne vegyuk be konstruktorba egy cutmapet?
1.7 +//constzero jo igy?
1.8 +
1.9 +//majd marci megmondja betegyem-e bfs-t meg resgraphot
1.10 +
1.11 +//constzero helyett az kell hogy flow-e vagy csak preflow, ha flow akor csak
1.12 +//excess[t]-t kell szmaolni
1.13 +
1.14 +/*
1.15 +Heuristics:
1.16 + 2 phase
1.17 + gap
1.18 + list 'level_list' on the nodes on level i implemented by hand
1.19 + stack 'active' on the active nodes on level i implemented by hand
1.20 + runs heuristic 'highest label' for H1*n relabels
1.21 + runs heuristic 'bound decrease' for H0*n relabels, starts with 'highest label'
1.22 +
1.23 +Parameters H0 and H1 are initialized to 20 and 10.
1.24 +
1.25 +Constructors:
1.26 +
1.27 +Preflow(Graph, Node, Node, CapMap, FlowMap, bool) : bool must be false if
1.28 + FlowMap is not constant zero, and should be true if it is
1.29 +
1.30 +Members:
1.31 +
1.32 +void run()
1.33 +
1.34 +T flowValue() : returns the value of a maximum flow
1.35 +
1.36 +void minMinCut(CutMap& M) : sets M to the characteristic vector of the
1.37 + minimum min cut. M should be a map of bools initialized to false.
1.38 +
1.39 +void maxMinCut(CutMap& M) : sets M to the characteristic vector of the
1.40 + maximum min cut. M should be a map of bools initialized to false.
1.41 +
1.42 +void minCut(CutMap& M) : sets M to the characteristic vector of
1.43 + a min cut. M should be a map of bools initialized to false.
1.44 +
1.45 +FIXME reset
1.46 +
1.47 +*/
1.48 +
1.49 +#ifndef HUGO_PREFLOW_H
1.50 +#define HUGO_PREFLOW_H
1.51 +
1.52 +#define H0 20
1.53 +#define H1 1
1.54 +
1.55 +#include <vector>
1.56 +#include <queue>
1.57 +#include <stack>
1.58 +
1.59 +namespace hugo {
1.60 +
1.61 + template <typename Graph, typename T,
1.62 + typename CapMap=typename Graph::template EdgeMap<T>,
1.63 + typename FlowMap=typename Graph::template EdgeMap<T> >
1.64 + class Preflow {
1.65 +
1.66 + typedef typename Graph::Node Node;
1.67 + typedef typename Graph::Edge Edge;
1.68 + typedef typename Graph::NodeIt NodeIt;
1.69 + typedef typename Graph::OutEdgeIt OutEdgeIt;
1.70 + typedef typename Graph::InEdgeIt InEdgeIt;
1.71 +
1.72 + const Graph& G;
1.73 + Node s;
1.74 + Node t;
1.75 + const CapMap& capacity;
1.76 + FlowMap& flow;
1.77 + T value;
1.78 + bool constzero;
1.79 + bool isflow;
1.80 +
1.81 + public:
1.82 + Preflow(Graph& _G, Node _s, Node _t, CapMap& _capacity,
1.83 + FlowMap& _flow, bool _constzero, bool _isflow ) :
1.84 + G(_G), s(_s), t(_t), capacity(_capacity), flow(_flow), constzero(_constzero), isflow(_isflow) {}
1.85 +
1.86 +
1.87 + void run() {
1.88 +
1.89 + value=0; //for the subsequent runs
1.90 +
1.91 + bool phase=0; //phase 0 is the 1st phase, phase 1 is the 2nd
1.92 + int n=G.nodeNum();
1.93 + int heur0=(int)(H0*n); //time while running 'bound decrease'
1.94 + int heur1=(int)(H1*n); //time while running 'highest label'
1.95 + int heur=heur1; //starting time interval (#of relabels)
1.96 + bool what_heur=1;
1.97 + /*
1.98 + what_heur is 0 in case 'bound decrease'
1.99 + and 1 in case 'highest label'
1.100 + */
1.101 + bool end=false;
1.102 + /*
1.103 + Needed for 'bound decrease', 'true'
1.104 + means no active nodes are above bound b.
1.105 + */
1.106 + int relabel=0;
1.107 + int k=n-2; //bound on the highest level under n containing a node
1.108 + int b=k; //bound on the highest level under n of an active node
1.109 +
1.110 + typename Graph::template NodeMap<int> level(G,n);
1.111 + typename Graph::template NodeMap<T> excess(G);
1.112 +
1.113 + std::vector<std::stack<Node> > active(n);
1.114 + /* std::vector<Node> active(n-1,INVALID);
1.115 + typename Graph::template NodeMap<Node> next(G,INVALID);
1.116 + //Stack of the active nodes in level i < n.
1.117 + //We use it in both phases.*/
1.118 +
1.119 + typename Graph::template NodeMap<Node> left(G,INVALID);
1.120 + typename Graph::template NodeMap<Node> right(G,INVALID);
1.121 + std::vector<Node> level_list(n,INVALID);
1.122 + /*
1.123 + List of the nodes in level i<n.
1.124 + */
1.125 +
1.126 +
1.127 + if ( constzero ) {
1.128 +
1.129 + /*Reverse_bfs from t, to find the starting level.*/
1.130 + level.set(t,0);
1.131 + std::queue<Node> bfs_queue;
1.132 + bfs_queue.push(t);
1.133 +
1.134 + while (!bfs_queue.empty()) {
1.135 +
1.136 + Node v=bfs_queue.front();
1.137 + bfs_queue.pop();
1.138 + int l=level[v]+1;
1.139 +
1.140 + InEdgeIt e;
1.141 + for(G.first(e,v); G.valid(e); G.next(e)) {
1.142 + Node w=G.tail(e);
1.143 + if ( level[w] == n && w != s ) {
1.144 + bfs_queue.push(w);
1.145 + Node first=level_list[l];
1.146 + if ( G.valid(first) ) left.set(first,w);
1.147 + right.set(w,first);
1.148 + level_list[l]=w;
1.149 + level.set(w, l);
1.150 + }
1.151 + }
1.152 + }
1.153 +
1.154 + //the starting flow
1.155 + OutEdgeIt e;
1.156 + for(G.first(e,s); G.valid(e); G.next(e))
1.157 + {
1.158 + T c=capacity[e];
1.159 + if ( c == 0 ) continue;
1.160 + Node w=G.head(e);
1.161 + if ( level[w] < n ) {
1.162 + if ( excess[w] == 0 && w!=t ) active[level[w]].push(w);
1.163 + flow.set(e, c);
1.164 + excess.set(w, excess[w]+c);
1.165 + }
1.166 + }
1.167 + }
1.168 + else
1.169 + {
1.170 +
1.171 + /*
1.172 + Reverse_bfs from t in the residual graph,
1.173 + to find the starting level.
1.174 + */
1.175 + level.set(t,0);
1.176 + std::queue<Node> bfs_queue;
1.177 + bfs_queue.push(t);
1.178 +
1.179 + while (!bfs_queue.empty()) {
1.180 +
1.181 + Node v=bfs_queue.front();
1.182 + bfs_queue.pop();
1.183 + int l=level[v]+1;
1.184 +
1.185 + InEdgeIt e;
1.186 + for(G.first(e,v); G.valid(e); G.next(e)) {
1.187 + if ( capacity[e] == flow[e] ) continue;
1.188 + Node w=G.tail(e);
1.189 + if ( level[w] == n && w != s ) {
1.190 + bfs_queue.push(w);
1.191 + Node first=level_list[l];
1.192 + if ( G.valid(first) ) left.set(first,w);
1.193 + right.set(w,first);
1.194 + level_list[l]=w;
1.195 + level.set(w, l);
1.196 + }
1.197 + }
1.198 +
1.199 + OutEdgeIt f;
1.200 + for(G.first(f,v); G.valid(f); G.next(f)) {
1.201 + if ( 0 == flow[f] ) continue;
1.202 + Node w=G.head(f);
1.203 + if ( level[w] == n && w != s ) {
1.204 + bfs_queue.push(w);
1.205 + Node first=level_list[l];
1.206 + if ( G.valid(first) ) left.set(first,w);
1.207 + right.set(w,first);
1.208 + level_list[l]=w;
1.209 + level.set(w, l);
1.210 + }
1.211 + }
1.212 + }
1.213 +
1.214 +
1.215 + /*
1.216 + Counting the excess
1.217 + */
1.218 +
1.219 + if ( !isflow ) {
1.220 + NodeIt v;
1.221 + for(G.first(v); G.valid(v); G.next(v)) {
1.222 + T exc=0;
1.223 +
1.224 + InEdgeIt e;
1.225 + for(G.first(e,v); G.valid(e); G.next(e)) exc+=flow[e];
1.226 + OutEdgeIt f;
1.227 + for(G.first(f,v); G.valid(f); G.next(f)) exc-=flow[f];
1.228 +
1.229 + excess.set(v,exc);
1.230 +
1.231 + //putting the active nodes into the stack
1.232 + int lev=level[v];
1.233 + if ( exc > 0 && lev < n && v != t ) active[lev].push(v);
1.234 + }
1.235 + } else {
1.236 + T exc=0;
1.237 +
1.238 + InEdgeIt e;
1.239 + for(G.first(e,t); G.valid(e); G.next(e)) exc+=flow[e];
1.240 + OutEdgeIt f;
1.241 + for(G.first(f,t); G.valid(f); G.next(f)) exc-=flow[f];
1.242 +
1.243 + excess.set(t,exc);
1.244 + }
1.245 +
1.246 +
1.247 + //the starting flow
1.248 + OutEdgeIt e;
1.249 + for(G.first(e,s); G.valid(e); G.next(e))
1.250 + {
1.251 + T rem=capacity[e]-flow[e];
1.252 + if ( rem == 0 ) continue;
1.253 + Node w=G.head(e);
1.254 + if ( level[w] < n ) {
1.255 + if ( excess[w] == 0 && w!=t ) active[level[w]].push(w);
1.256 + flow.set(e, capacity[e]);
1.257 + excess.set(w, excess[w]+rem);
1.258 + }
1.259 + }
1.260 +
1.261 + InEdgeIt f;
1.262 + for(G.first(f,s); G.valid(f); G.next(f))
1.263 + {
1.264 + if ( flow[f] == 0 ) continue;
1.265 + Node w=G.tail(f);
1.266 + if ( level[w] < n ) {
1.267 + if ( excess[w] == 0 && w!=t ) active[level[w]].push(w);
1.268 + excess.set(w, excess[w]+flow[f]);
1.269 + flow.set(f, 0);
1.270 + }
1.271 + }
1.272 + }
1.273 +
1.274 +
1.275 +
1.276 +
1.277 + /*
1.278 + End of preprocessing
1.279 + */
1.280 +
1.281 +
1.282 +
1.283 + /*
1.284 + Push/relabel on the highest level active nodes.
1.285 + */
1.286 + while ( true ) {
1.287 +
1.288 + if ( b == 0 ) {
1.289 + if ( phase ) break;
1.290 +
1.291 + if ( !what_heur && !end && k > 0 ) {
1.292 + b=k;
1.293 + end=true;
1.294 + } else {
1.295 + phase=1;
1.296 + level.set(s,0);
1.297 + std::queue<Node> bfs_queue;
1.298 + bfs_queue.push(s);
1.299 +
1.300 + while (!bfs_queue.empty()) {
1.301 +
1.302 + Node v=bfs_queue.front();
1.303 + bfs_queue.pop();
1.304 + int l=level[v]+1;
1.305 +
1.306 + InEdgeIt e;
1.307 + for(G.first(e,v); G.valid(e); G.next(e)) {
1.308 + if ( capacity[e] == flow[e] ) continue;
1.309 + Node u=G.tail(e);
1.310 + if ( level[u] >= n ) {
1.311 + bfs_queue.push(u);
1.312 + level.set(u, l);
1.313 + if ( excess[u] > 0 ) active[l].push(u);
1.314 + }
1.315 + }
1.316 +
1.317 + OutEdgeIt f;
1.318 + for(G.first(f,v); G.valid(f); G.next(f)) {
1.319 + if ( 0 == flow[f] ) continue;
1.320 + Node u=G.head(f);
1.321 + if ( level[u] >= n ) {
1.322 + bfs_queue.push(u);
1.323 + level.set(u, l);
1.324 + if ( excess[u] > 0 ) active[l].push(u);
1.325 + }
1.326 + }
1.327 + }
1.328 + b=n-2;
1.329 + }
1.330 +
1.331 + }
1.332 +
1.333 +
1.334 + ///
1.335 + if ( active[b].empty() ) --b;
1.336 + else {
1.337 + end=false;
1.338 +
1.339 + Node w=active[b].top();
1.340 + active[b].pop();
1.341 + int lev=level[w];
1.342 + T exc=excess[w];
1.343 + int newlevel=n; //bound on the next level of w
1.344 +
1.345 + OutEdgeIt e;
1.346 + for(G.first(e,w); G.valid(e); G.next(e)) {
1.347 +
1.348 + if ( flow[e] == capacity[e] ) continue;
1.349 + Node v=G.head(e);
1.350 + //e=wv
1.351 +
1.352 + if( lev > level[v] ) {
1.353 + /*Push is allowed now*/
1.354 +
1.355 + if ( excess[v]==0 && v!=t && v!=s ) {
1.356 + int lev_v=level[v];
1.357 + active[lev_v].push(v);
1.358 + }
1.359 +
1.360 + T cap=capacity[e];
1.361 + T flo=flow[e];
1.362 + T remcap=cap-flo;
1.363 +
1.364 + if ( remcap >= exc ) {
1.365 + /*A nonsaturating push.*/
1.366 +
1.367 + flow.set(e, flo+exc);
1.368 + excess.set(v, excess[v]+exc);
1.369 + exc=0;
1.370 + break;
1.371 +
1.372 + } else {
1.373 + /*A saturating push.*/
1.374 +
1.375 + flow.set(e, cap);
1.376 + excess.set(v, excess[v]+remcap);
1.377 + exc-=remcap;
1.378 + }
1.379 + } else if ( newlevel > level[v] ){
1.380 + newlevel = level[v];
1.381 + }
1.382 +
1.383 + } //for out edges wv
1.384 +
1.385 +
1.386 + if ( exc > 0 ) {
1.387 + InEdgeIt e;
1.388 + for(G.first(e,w); G.valid(e); G.next(e)) {
1.389 +
1.390 + if( flow[e] == 0 ) continue;
1.391 + Node v=G.tail(e);
1.392 + //e=vw
1.393 +
1.394 + if( lev > level[v] ) {
1.395 + /*Push is allowed now*/
1.396 +
1.397 + if ( excess[v]==0 && v!=t && v!=s ) {
1.398 + int lev_v=level[v];
1.399 + active[lev_v].push(v);
1.400 + }
1.401 +
1.402 + T flo=flow[e];
1.403 +
1.404 + if ( flo >= exc ) {
1.405 + /*A nonsaturating push.*/
1.406 +
1.407 + flow.set(e, flo-exc);
1.408 + excess.set(v, excess[v]+exc);
1.409 + exc=0;
1.410 + break;
1.411 + } else {
1.412 + /*A saturating push.*/
1.413 +
1.414 + excess.set(v, excess[v]+flo);
1.415 + exc-=flo;
1.416 + flow.set(e,0);
1.417 + }
1.418 + } else if ( newlevel > level[v] ) {
1.419 + newlevel = level[v];
1.420 + }
1.421 + } //for in edges vw
1.422 +
1.423 + } // if w still has excess after the out edge for cycle
1.424 +
1.425 + excess.set(w, exc);
1.426 + /// push
1.427 +
1.428 +
1.429 + /*
1.430 + Relabel
1.431 + */
1.432 +
1.433 +
1.434 + if ( exc > 0 ) {
1.435 + //now 'lev' is the old level of w
1.436 +
1.437 + if ( phase ) {
1.438 + level.set(w,++newlevel);
1.439 + active[newlevel].push(w);
1.440 + b=newlevel;
1.441 + } else {
1.442 + //unlacing starts
1.443 + Node right_n=right[w];
1.444 + Node left_n=left[w];
1.445 +
1.446 + if ( G.valid(right_n) ) {
1.447 + if ( G.valid(left_n) ) {
1.448 + right.set(left_n, right_n);
1.449 + left.set(right_n, left_n);
1.450 + } else {
1.451 + level_list[lev]=right_n;
1.452 + left.set(right_n, INVALID);
1.453 + }
1.454 + } else {
1.455 + if ( G.valid(left_n) ) {
1.456 + right.set(left_n, INVALID);
1.457 + } else {
1.458 + level_list[lev]=INVALID;
1.459 + }
1.460 + }
1.461 + //unlacing ends
1.462 +
1.463 + if ( !G.valid(level_list[lev]) ) {
1.464 +
1.465 + //gapping starts
1.466 + for (int i=lev; i!=k ; ) {
1.467 + Node v=level_list[++i];
1.468 + while ( G.valid(v) ) {
1.469 + level.set(v,n);
1.470 + v=right[v];
1.471 + }
1.472 + level_list[i]=INVALID;
1.473 + if ( !what_heur ) {
1.474 + while ( !active[i].empty() ) {
1.475 + active[i].pop(); //FIXME: ezt szebben kene
1.476 + }
1.477 + }
1.478 + }
1.479 +
1.480 + level.set(w,n);
1.481 + b=lev-1;
1.482 + k=b;
1.483 + //gapping ends
1.484 +
1.485 + } else {
1.486 +
1.487 + if ( newlevel == n ) level.set(w,n);
1.488 + else {
1.489 + level.set(w,++newlevel);
1.490 + active[newlevel].push(w);
1.491 + if ( what_heur ) b=newlevel;
1.492 + if ( k < newlevel ) ++k; //now k=newlevel
1.493 + Node first=level_list[newlevel];
1.494 + if ( G.valid(first) ) left.set(first,w);
1.495 + right.set(w,first);
1.496 + left.set(w,INVALID);
1.497 + level_list[newlevel]=w;
1.498 + }
1.499 + }
1.500 +
1.501 +
1.502 + ++relabel;
1.503 + if ( relabel >= heur ) {
1.504 + relabel=0;
1.505 + if ( what_heur ) {
1.506 + what_heur=0;
1.507 + heur=heur0;
1.508 + end=false;
1.509 + } else {
1.510 + what_heur=1;
1.511 + heur=heur1;
1.512 + b=k;
1.513 + }
1.514 + }
1.515 + } //phase 0
1.516 +
1.517 +
1.518 + } // if ( exc > 0 )
1.519 +
1.520 +
1.521 + } // if stack[b] is nonempty
1.522 +
1.523 + } // while(true)
1.524 +
1.525 +
1.526 + value = excess[t];
1.527 + /*Max flow value.*/
1.528 +
1.529 + } //void run()
1.530 +
1.531 +
1.532 +
1.533 +
1.534 +
1.535 + /*
1.536 + Returns the maximum value of a flow.
1.537 + */
1.538 +
1.539 + T flowValue() {
1.540 + return value;
1.541 + }
1.542 +
1.543 +
1.544 + FlowMap Flow() {
1.545 + return flow;
1.546 + }
1.547 +
1.548 +
1.549 + void Flow(FlowMap& _flow ) {
1.550 + NodeIt v;
1.551 + for(G.first(v) ; G.valid(v); G.next(v))
1.552 + _flow.set(v,flow[v]);
1.553 + }
1.554 +
1.555 +
1.556 +
1.557 + /*
1.558 + Returns the minimum min cut, by a bfs from s in the residual graph.
1.559 + */
1.560 +
1.561 + template<typename _CutMap>
1.562 + void minMinCut(_CutMap& M) {
1.563 +
1.564 + std::queue<Node> queue;
1.565 +
1.566 + M.set(s,true);
1.567 + queue.push(s);
1.568 +
1.569 + while (!queue.empty()) {
1.570 + Node w=queue.front();
1.571 + queue.pop();
1.572 +
1.573 + OutEdgeIt e;
1.574 + for(G.first(e,w) ; G.valid(e); G.next(e)) {
1.575 + Node v=G.head(e);
1.576 + if (!M[v] && flow[e] < capacity[e] ) {
1.577 + queue.push(v);
1.578 + M.set(v, true);
1.579 + }
1.580 + }
1.581 +
1.582 + InEdgeIt f;
1.583 + for(G.first(f,w) ; G.valid(f); G.next(f)) {
1.584 + Node v=G.tail(f);
1.585 + if (!M[v] && flow[f] > 0 ) {
1.586 + queue.push(v);
1.587 + M.set(v, true);
1.588 + }
1.589 + }
1.590 + }
1.591 + }
1.592 +
1.593 +
1.594 +
1.595 + /*
1.596 + Returns the maximum min cut, by a reverse bfs
1.597 + from t in the residual graph.
1.598 + */
1.599 +
1.600 + template<typename _CutMap>
1.601 + void maxMinCut(_CutMap& M) {
1.602 +
1.603 + std::queue<Node> queue;
1.604 +
1.605 + M.set(t,true);
1.606 + queue.push(t);
1.607 +
1.608 + while (!queue.empty()) {
1.609 + Node w=queue.front();
1.610 + queue.pop();
1.611 +
1.612 +
1.613 + InEdgeIt e;
1.614 + for(G.first(e,w) ; G.valid(e); G.next(e)) {
1.615 + Node v=G.tail(e);
1.616 + if (!M[v] && flow[e] < capacity[e] ) {
1.617 + queue.push(v);
1.618 + M.set(v, true);
1.619 + }
1.620 + }
1.621 +
1.622 + OutEdgeIt f;
1.623 + for(G.first(f,w) ; G.valid(f); G.next(f)) {
1.624 + Node v=G.head(f);
1.625 + if (!M[v] && flow[f] > 0 ) {
1.626 + queue.push(v);
1.627 + M.set(v, true);
1.628 + }
1.629 + }
1.630 + }
1.631 +
1.632 + NodeIt v;
1.633 + for(G.first(v) ; G.valid(v); G.next(v)) {
1.634 + M.set(v, !M[v]);
1.635 + }
1.636 +
1.637 + }
1.638 +
1.639 +
1.640 +
1.641 + template<typename CutMap>
1.642 + void minCut(CutMap& M) {
1.643 + minMinCut(M);
1.644 + }
1.645 +
1.646 +
1.647 + void resetTarget (Node _t) {t=_t;}
1.648 + void resetSource (Node _s) {s=_s;}
1.649 +
1.650 + void resetCap (CapMap _cap) {capacity=_cap;}
1.651 +
1.652 + void resetFlow (FlowMap _flow, bool _constzero) {
1.653 + flow=_flow;
1.654 + constzero=_constzero;
1.655 + }
1.656 +
1.657 +
1.658 +
1.659 + };
1.660 +
1.661 +} //namespace hugo
1.662 +
1.663 +#endif //PREFLOW_H
1.664 +
1.665 +
1.666 +
1.667 +
2.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
2.2 +++ b/src/work/jacint/preflow_excess_test.cc Tue Apr 27 10:27:34 2004 +0000
2.3 @@ -0,0 +1,179 @@
2.4 +/*
2.5 +The only difference between preflow.h and preflow_res.h is that the latter
2.6 +uses the ResGraphWrapper, while the first does not. (Bfs is implemented by
2.7 +hand in both.) This test program runs Preflow and PreflowRes on the same
2.8 +graph, tests the result of these implementations and writes the running time
2.9 +of them. */
2.10 +#include <iostream>
2.11 +
2.12 +#include <smart_graph.h>
2.13 +#include <dimacs.h>
2.14 +#include <preflow_excess.h>
2.15 +#include <time_measure.h>
2.16 +
2.17 +using namespace hugo;
2.18 +
2.19 +int main(int, char **) {
2.20 +
2.21 + typedef SmartGraph Graph;
2.22 +
2.23 + typedef Graph::Node Node;
2.24 + typedef Graph::EdgeIt EdgeIt;
2.25 +
2.26 + Graph G;
2.27 + Node s, t;
2.28 + Graph::EdgeMap<int> cap(G);
2.29 + readDimacsMaxFlow(std::cin, G, s, t, cap);
2.30 + Timer ts;
2.31 +
2.32 + std::cout <<
2.33 + "\n Are we slower?"
2.34 + <<std::endl;
2.35 + std::cout <<
2.36 + "\n Running preflow.h on a graph with " <<
2.37 + G.nodeNum() << " nodes and " << G.edgeNum() << " edges..."
2.38 + << std::endl<<std::endl;
2.39 +
2.40 +
2.41 +
2.42 +
2.43 +
2.44 +
2.45 +
2.46 + Graph::EdgeMap<int> flow(G,0);
2.47 + Preflow<Graph, int> max_flow_test(G, s, t, cap, flow, 0 , 0);
2.48 + ts.reset();
2.49 + max_flow_test.run();
2.50 + std::cout << "Elapsed time from a preflow: " << std::endl
2.51 + <<ts << std::endl;
2.52 +
2.53 + Graph::NodeMap<bool> mincut(G);
2.54 + max_flow_test.minMinCut(mincut);
2.55 + int min_min_cut_value=0;
2.56 + EdgeIt e;
2.57 + for(G.first(e); G.valid(e); G.next(e)) {
2.58 + if (mincut[G.tail(e)] && !mincut[G.head(e)]) min_min_cut_value+=cap[e];
2.59 + }
2.60 +
2.61 + Graph::NodeMap<bool> cut(G);
2.62 + max_flow_test.minCut(cut);
2.63 + int min_cut_value=0;
2.64 + for(G.first(e); G.valid(e); G.next(e)) {
2.65 + if (cut[G.tail(e)] && !cut[G.head(e)])
2.66 + min_cut_value+=cap[e];
2.67 + }
2.68 +
2.69 + Graph::NodeMap<bool> maxcut(G);
2.70 + max_flow_test.maxMinCut(maxcut);
2.71 + int max_min_cut_value=0;
2.72 + for(G.first(e); G.valid(e); G.next(e)) {
2.73 + if (maxcut[G.tail(e)] && !maxcut[G.head(e)])
2.74 + max_min_cut_value+=cap[e];
2.75 + }
2.76 +
2.77 + std::cout << "\n Checking the result: " <<std::endl;
2.78 + std::cout << "Flow value: "<< max_flow_test.flowValue() << std::endl;
2.79 + std::cout << "Min cut value: "<< min_cut_value << std::endl;
2.80 + std::cout << "Min min cut value: "<< min_min_cut_value << std::endl;
2.81 + std::cout << "Max min cut value: "<< max_min_cut_value <<
2.82 + std::endl;
2.83 +
2.84 + if ( max_flow_test.flowValue() == min_cut_value &&
2.85 + min_cut_value == min_min_cut_value &&
2.86 + min_min_cut_value == max_min_cut_value )
2.87 + std::cout << "They are equal! " <<std::endl<< std::endl<<"\n";
2.88 +
2.89 +
2.90 +
2.91 +
2.92 +
2.93 +
2.94 + Graph::EdgeMap<int> flow2(G,0);
2.95 + Preflow<Graph, int> max_flow_test2(G, s, t, cap, flow2, 0 , 1);
2.96 + ts.reset();
2.97 + max_flow_test2.run();
2.98 + std::cout << "Elapsed time from a flow: " << std::endl
2.99 + << ts << std::endl;
2.100 +
2.101 + Graph::NodeMap<bool> mincut2(G);
2.102 + max_flow_test2.minMinCut(mincut2);
2.103 + int min_min_cut_value2=0;
2.104 + for(G.first(e); G.valid(e); G.next(e)) {
2.105 + if (mincut2[G.tail(e)] && !mincut2[G.head(e)]) min_min_cut_value2+=cap[e];
2.106 + }
2.107 +
2.108 + Graph::NodeMap<bool> cut2(G);
2.109 + max_flow_test2.minCut(cut2);
2.110 + int min_cut_value2=0;
2.111 + for(G.first(e); G.valid(e); G.next(e)) {
2.112 + if (cut2[G.tail(e)] && !cut2[G.head(e)])
2.113 + min_cut_value2+=cap[e];
2.114 + }
2.115 +
2.116 + Graph::NodeMap<bool> maxcut2(G);
2.117 + max_flow_test2.maxMinCut(maxcut2);
2.118 + int max_min_cut_value2=0;
2.119 + for(G.first(e); G.valid(e); G.next(e)) {
2.120 + if (maxcut2[G.tail(e)] && !maxcut2[G.head(e)])
2.121 + max_min_cut_value2+=cap[e];
2.122 + }
2.123 +
2.124 + std::cout << "\n Checking the result: " <<std::endl;
2.125 + std::cout << "Flow value: "<< max_flow_test2.flowValue() << std::endl;
2.126 + std::cout << "Min cut value: "<< min_cut_value2 << std::endl;
2.127 + std::cout << "Min min cut value: "<< min_min_cut_value2 << std::endl;
2.128 + std::cout << "Max min cut value: "<< max_min_cut_value2 <<
2.129 + std::endl;
2.130 + if ( max_flow_test.flowValue() == min_cut_value &&
2.131 + min_cut_value == min_min_cut_value &&
2.132 + min_min_cut_value == max_min_cut_value )
2.133 + std::cout << "They are equal! " <<std::endl;
2.134 +
2.135 +
2.136 +
2.137 +
2.138 +
2.139 + Graph::EdgeMap<int> flow3(G,0);
2.140 + Preflow<Graph, int> max_flow_test3(G, s, t, cap, flow3, 1 , 1);
2.141 + ts.reset();
2.142 + max_flow_test3.run();
2.143 + std::cout << "Elapsed time from a const zero flow: " << std::endl
2.144 + <<ts << std::endl;
2.145 +
2.146 + Graph::NodeMap<bool> mincut3(G);
2.147 + max_flow_test3.minMinCut(mincut3);
2.148 + int min_min_cut_value3=0;
2.149 + for(G.first(e); G.valid(e); G.next(e)) {
2.150 + if (mincut3[G.tail(e)] && !mincut3[G.head(e)]) min_min_cut_value3+=cap[e];
2.151 + }
2.152 +
2.153 + Graph::NodeMap<bool> cut3(G);
2.154 + max_flow_test3.minCut(cut3);
2.155 + int min_cut_value3=0;
2.156 + for(G.first(e); G.valid(e); G.next(e)) {
2.157 + if (cut3[G.tail(e)] && !cut3[G.head(e)])
2.158 + min_cut_value3+=cap[e];
2.159 + }
2.160 +
2.161 + Graph::NodeMap<bool> maxcut3(G);
2.162 + max_flow_test3.maxMinCut(maxcut3);
2.163 + int max_min_cut_value3=0;
2.164 + for(G.first(e); G.valid(e); G.next(e)) {
2.165 + if (maxcut3[G.tail(e)] && !maxcut3[G.head(e)])
2.166 + max_min_cut_value3+=cap[e];
2.167 + }
2.168 +
2.169 + std::cout << "\n Checking the result: " <<std::endl;
2.170 + std::cout << "Flow value: "<< max_flow_test3.flowValue() << std::endl;
2.171 + std::cout << "Min cut value: "<< min_cut_value3 << std::endl;
2.172 + std::cout << "Min min cut value: "<< min_min_cut_value3 << std::endl;
2.173 + std::cout << "Max min cut value: "<< max_min_cut_value3 <<
2.174 + std::endl;
2.175 +
2.176 + if ( max_flow_test3.flowValue() == min_cut_value3 &&
2.177 + min_cut_value3 == min_min_cut_value3 &&
2.178 + min_min_cut_value3 == max_min_cut_value3 )
2.179 + std::cout << "They are equal! " <<std::endl<< std::endl<<"\n";
2.180 +
2.181 + return 0;
2.182 +}