changeset 835 | eb9587f09b42 |
parent 773 | ce9438c5a82d |
7:9f1a4a2a8c63 | 8:9a220d72378d |
---|---|
3 #define HUGO_MAX_FLOW_H |
3 #define HUGO_MAX_FLOW_H |
4 |
4 |
5 #include <vector> |
5 #include <vector> |
6 #include <queue> |
6 #include <queue> |
7 |
7 |
8 #include <hugo/graph_wrapper.h> |
8 //#include <hugo/graph_wrapper.h> |
9 #include <hugo/invalid.h> |
9 #include <hugo/invalid.h> |
10 #include <hugo/maps.h> |
10 #include <hugo/maps.h> |
11 |
11 |
12 /// \file |
12 /// \file |
13 /// \ingroup flowalgs |
13 /// \ingroup flowalgs |
60 Node s; |
60 Node s; |
61 Node t; |
61 Node t; |
62 const CapMap* capacity; |
62 const CapMap* capacity; |
63 FlowMap* flow; |
63 FlowMap* flow; |
64 int n; //the number of nodes of G |
64 int n; //the number of nodes of G |
65 typedef ResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW; |
65 // typedef ResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW; |
66 //typedef ExpResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW; |
66 //typedef ExpResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW; |
67 typedef typename ResGW::OutEdgeIt ResGWOutEdgeIt; |
67 // typedef typename ResGW::OutEdgeIt ResGWOutEdgeIt; |
68 typedef typename ResGW::Edge ResGWEdge; |
68 // typedef typename ResGW::Edge ResGWEdge; |
69 typedef typename Graph::template NodeMap<int> ReachedMap; |
69 typedef typename Graph::template NodeMap<int> ReachedMap; |
70 |
70 |
71 |
71 |
72 //level works as a bool map in augmenting path algorithms and is |
72 //level works as a bool map in augmenting path algorithms and is |
73 //used by bfs for storing reached information. In preflow, it |
73 //used by bfs for storing reached information. In preflow, it |
110 }; |
110 }; |
111 |
111 |
112 /// Do not needle this flag only if necessary. |
112 /// Do not needle this flag only if necessary. |
113 StatusEnum status; |
113 StatusEnum status; |
114 |
114 |
115 // int number_of_augmentations; |
115 // int number_of_augmentations; |
116 |
116 |
117 |
117 |
118 // template<typename IntMap> |
118 // template<typename IntMap> |
119 // class TrickyReachedMap { |
119 // class TrickyReachedMap { |
120 // protected: |
120 // protected: |
121 // IntMap* map; |
121 // IntMap* map; |
122 // int* number_of_augmentations; |
122 // int* number_of_augmentations; |
123 // public: |
123 // public: |
124 // TrickyReachedMap(IntMap& _map, int& _number_of_augmentations) : |
124 // TrickyReachedMap(IntMap& _map, int& _number_of_augmentations) : |
125 // map(&_map), number_of_augmentations(&_number_of_augmentations) { } |
125 // map(&_map), number_of_augmentations(&_number_of_augmentations) { } |
126 // void set(const Node& n, bool b) { |
126 // void set(const Node& n, bool b) { |
127 // if (b) |
127 // if (b) |
128 // map->set(n, *number_of_augmentations); |
128 // map->set(n, *number_of_augmentations); |
129 // else |
129 // else |
130 // map->set(n, *number_of_augmentations-1); |
130 // map->set(n, *number_of_augmentations-1); |
131 // } |
131 // } |
132 // bool operator[](const Node& n) const { |
132 // bool operator[](const Node& n) const { |
133 // return (*map)[n]==*number_of_augmentations; |
133 // return (*map)[n]==*number_of_augmentations; |
134 // } |
134 // } |
135 // }; |
135 // }; |
136 |
136 |
137 ///Constructor |
137 ///Constructor |
138 |
138 |
139 ///\todo Document, please. |
139 ///\todo Document, please. |
140 /// |
140 /// |
232 b=k; |
232 b=k; |
233 end=true; |
233 end=true; |
234 } else break; |
234 } else break; |
235 } |
235 } |
236 |
236 |
237 if ( !g->valid(first[b]) ) --b; |
237 if ( first[b]==INVALID ) --b; |
238 else { |
238 else { |
239 end=false; |
239 end=false; |
240 Node w=first[b]; |
240 Node w=first[b]; |
241 first[b]=next[w]; |
241 first[b]=next[w]; |
242 int newlevel=push(w, next, first); |
242 int newlevel=push(w, next, first); |
287 |
287 |
288 Node v=bfs_queue.front(); |
288 Node v=bfs_queue.front(); |
289 bfs_queue.pop(); |
289 bfs_queue.pop(); |
290 int l=level[v]+1; |
290 int l=level[v]+1; |
291 |
291 |
292 InEdgeIt e; |
292 for(InEdgeIt e(*g,v); e!=INVALID; ++e) { |
293 for(g->first(e,v); g->valid(e); g->next(e)) { |
|
294 if ( (*capacity)[e] <= (*flow)[e] ) continue; |
293 if ( (*capacity)[e] <= (*flow)[e] ) continue; |
295 Node u=g->tail(e); |
294 Node u=g->tail(e); |
296 if ( level[u] >= n ) { |
295 if ( level[u] >= n ) { |
297 bfs_queue.push(u); |
296 bfs_queue.push(u); |
298 level.set(u, l); |
297 level.set(u, l); |
301 first[l]=u; |
300 first[l]=u; |
302 } |
301 } |
303 } |
302 } |
304 } |
303 } |
305 |
304 |
306 OutEdgeIt f; |
305 for(OutEdgeIt e(*g,v); e!=INVALID; ++e) { |
307 for(g->first(f,v); g->valid(f); g->next(f)) { |
306 if ( 0 >= (*flow)[e] ) continue; |
308 if ( 0 >= (*flow)[f] ) continue; |
307 Node u=g->head(e); |
309 Node u=g->head(f); |
|
310 if ( level[u] >= n ) { |
308 if ( level[u] >= n ) { |
311 bfs_queue.push(u); |
309 bfs_queue.push(u); |
312 level.set(u, l); |
310 level.set(u, l); |
313 if ( excess[u] > 0 ) { |
311 if ( excess[u] > 0 ) { |
314 next.set(u,first[l]); |
312 next.set(u,first[l]); |
321 |
319 |
322 while ( true ) { |
320 while ( true ) { |
323 |
321 |
324 if ( b == 0 ) break; |
322 if ( b == 0 ) break; |
325 |
323 |
326 if ( !g->valid(first[b]) ) --b; |
324 if ( first[b]==INVALID ) --b; |
327 else { |
325 else { |
328 |
326 |
329 Node w=first[b]; |
327 Node w=first[b]; |
330 first[b]=next[w]; |
328 first[b]=next[w]; |
331 int newlevel=push(w,next, first/*active*/); |
329 int newlevel=push(w,next, first/*active*/); |
349 /// Returns the excess of the target node \ref t. |
347 /// Returns the excess of the target node \ref t. |
350 /// After running \ref preflowPhase1, this is the value of |
348 /// After running \ref preflowPhase1, this is the value of |
351 /// the maximum flow. |
349 /// the maximum flow. |
352 /// It can be called already after running \ref preflowPhase1. |
350 /// It can be called already after running \ref preflowPhase1. |
353 Num flowValue() const { |
351 Num flowValue() const { |
354 // Num a=0; |
352 // Num a=0; |
355 // for(InEdgeIt e(*g,t);g->valid(e);g->next(e)) a+=(*flow)[e]; |
353 // for(InEdgeIt e(*g,t);g->valid(e);g->next(e)) a+=(*flow)[e]; |
356 // for(OutEdgeIt e(*g,t);g->valid(e);g->next(e)) a-=(*flow)[e]; |
354 // for(OutEdgeIt e(*g,t);g->valid(e);g->next(e)) a-=(*flow)[e]; |
357 // return a; |
355 // return a; |
358 return excess[t]; |
356 return excess[t]; |
359 //marci figyu: excess[t] epp ezt adja preflow 1. fazisa utan |
357 //marci figyu: excess[t] epp ezt adja preflow 1. fazisa utan |
360 } |
358 } |
361 |
359 |
362 |
360 |
372 /// actual state |
370 /// actual state |
373 /// of the class. This enables us to determine which methods are valid |
371 /// of the class. This enables us to determine which methods are valid |
374 /// for MinCut computation |
372 /// for MinCut computation |
375 template<typename _CutMap> |
373 template<typename _CutMap> |
376 void actMinCut(_CutMap& M) const { |
374 void actMinCut(_CutMap& M) const { |
377 NodeIt v; |
|
378 switch (status) { |
375 switch (status) { |
379 case AFTER_PRE_FLOW_PHASE_1: |
376 case AFTER_PRE_FLOW_PHASE_1: |
380 for(g->first(v); g->valid(v); g->next(v)) { |
377 for(NodeIt v(*g); v!=INVALID; ++v) { |
381 if (level[v] < n) { |
378 if (level[v] < n) { |
382 M.set(v, false); |
379 M.set(v, false); |
383 } else { |
380 } else { |
384 M.set(v, true); |
381 M.set(v, true); |
385 } |
382 } |
386 } |
383 } |
387 break; |
384 break; |
388 case AFTER_PRE_FLOW_PHASE_2: |
385 case AFTER_PRE_FLOW_PHASE_2: |
389 case AFTER_NOTHING: |
386 case AFTER_NOTHING: |
390 case AFTER_AUGMENTING: |
387 case AFTER_AUGMENTING: |
391 case AFTER_FAST_AUGMENTING: |
388 case AFTER_FAST_AUGMENTING: |
392 minMinCut(M); |
389 minMinCut(M); |
393 break; |
390 break; |
394 } |
391 } |
395 } |
392 } |
396 |
393 |
410 |
407 |
411 while (!queue.empty()) { |
408 while (!queue.empty()) { |
412 Node w=queue.front(); |
409 Node w=queue.front(); |
413 queue.pop(); |
410 queue.pop(); |
414 |
411 |
415 OutEdgeIt e; |
412 for(OutEdgeIt e(*g,w) ; e!=INVALID; ++e) { |
416 for(g->first(e,w) ; g->valid(e); g->next(e)) { |
|
417 Node v=g->head(e); |
413 Node v=g->head(e); |
418 if (!M[v] && (*flow)[e] < (*capacity)[e] ) { |
414 if (!M[v] && (*flow)[e] < (*capacity)[e] ) { |
419 queue.push(v); |
415 queue.push(v); |
420 M.set(v, true); |
416 M.set(v, true); |
421 } |
417 } |
422 } |
418 } |
423 |
419 |
424 InEdgeIt f; |
420 for(InEdgeIt e(*g,w) ; e!=INVALID; ++e) { |
425 for(g->first(f,w) ; g->valid(f); g->next(f)) { |
421 Node v=g->tail(e); |
426 Node v=g->tail(f); |
422 if (!M[v] && (*flow)[e] > 0 ) { |
427 if (!M[v] && (*flow)[f] > 0 ) { |
|
428 queue.push(v); |
423 queue.push(v); |
429 M.set(v, true); |
424 M.set(v, true); |
430 } |
425 } |
431 } |
426 } |
432 } |
427 } |
440 ///\pre M should be a node map of bools initialized to false. |
435 ///\pre M should be a node map of bools initialized to false. |
441 ///\pre \c flow must be a maximum flow. |
436 ///\pre \c flow must be a maximum flow. |
442 template<typename _CutMap> |
437 template<typename _CutMap> |
443 void maxMinCut(_CutMap& M) const { |
438 void maxMinCut(_CutMap& M) const { |
444 |
439 |
445 NodeIt v; |
440 for(NodeIt v(*g) ; v!=INVALID; ++v) M.set(v, true); |
446 for(g->first(v) ; g->valid(v); g->next(v)) { |
|
447 M.set(v, true); |
|
448 } |
|
449 |
441 |
450 std::queue<Node> queue; |
442 std::queue<Node> queue; |
451 |
443 |
452 M.set(t,false); |
444 M.set(t,false); |
453 queue.push(t); |
445 queue.push(t); |
454 |
446 |
455 while (!queue.empty()) { |
447 while (!queue.empty()) { |
456 Node w=queue.front(); |
448 Node w=queue.front(); |
457 queue.pop(); |
449 queue.pop(); |
458 |
450 |
459 InEdgeIt e; |
451 for(InEdgeIt e(*g,w) ; e!=INVALID; ++e) { |
460 for(g->first(e,w) ; g->valid(e); g->next(e)) { |
|
461 Node v=g->tail(e); |
452 Node v=g->tail(e); |
462 if (M[v] && (*flow)[e] < (*capacity)[e] ) { |
453 if (M[v] && (*flow)[e] < (*capacity)[e] ) { |
463 queue.push(v); |
454 queue.push(v); |
464 M.set(v, false); |
455 M.set(v, false); |
465 } |
456 } |
466 } |
457 } |
467 |
458 |
468 OutEdgeIt f; |
459 for(OutEdgeIt e(*g,w) ; e!=INVALID; ++e) { |
469 for(g->first(f,w) ; g->valid(f); g->next(f)) { |
460 Node v=g->head(e); |
470 Node v=g->head(f); |
461 if (M[v] && (*flow)[e] > 0 ) { |
471 if (M[v] && (*flow)[f] > 0 ) { |
|
472 queue.push(v); |
462 queue.push(v); |
473 M.set(v, false); |
463 M.set(v, false); |
474 } |
464 } |
475 } |
465 } |
476 } |
466 } |
516 |
506 |
517 int lev=level[w]; |
507 int lev=level[w]; |
518 Num exc=excess[w]; |
508 Num exc=excess[w]; |
519 int newlevel=n; //bound on the next level of w |
509 int newlevel=n; //bound on the next level of w |
520 |
510 |
521 OutEdgeIt e; |
511 for(OutEdgeIt e(*g,w) ; e!=INVALID; ++e) { |
522 for(g->first(e,w); g->valid(e); g->next(e)) { |
|
523 |
|
524 if ( (*flow)[e] >= (*capacity)[e] ) continue; |
512 if ( (*flow)[e] >= (*capacity)[e] ) continue; |
525 Node v=g->head(e); |
513 Node v=g->head(e); |
526 |
514 |
527 if( lev > level[v] ) { //Push is allowed now |
515 if( lev > level[v] ) { //Push is allowed now |
528 |
516 |
529 if ( excess[v]<=0 && v!=t && v!=s ) { |
517 if ( excess[v]<=0 && v!=t && v!=s ) { |
530 next.set(v,first[level[v]]); |
518 next.set(v,first[level[v]]); |
531 first[level[v]]=v; |
519 first[level[v]]=v; |
532 } |
520 } |
533 |
521 |
534 Num cap=(*capacity)[e]; |
522 Num cap=(*capacity)[e]; |
535 Num flo=(*flow)[e]; |
523 Num flo=(*flow)[e]; |
536 Num remcap=cap-flo; |
524 Num remcap=cap-flo; |
537 |
525 |
538 if ( remcap >= exc ) { //A nonsaturating push. |
526 if ( remcap >= exc ) { //A nonsaturating push. |
539 |
527 |
540 flow->set(e, flo+exc); |
528 flow->set(e, flo+exc); |
541 excess.set(v, excess[v]+exc); |
529 excess.set(v, excess[v]+exc); |
542 exc=0; |
530 exc=0; |
543 break; |
531 break; |
544 |
532 |
549 } |
537 } |
550 } else if ( newlevel > level[v] ) newlevel = level[v]; |
538 } else if ( newlevel > level[v] ) newlevel = level[v]; |
551 } //for out edges wv |
539 } //for out edges wv |
552 |
540 |
553 if ( exc > 0 ) { |
541 if ( exc > 0 ) { |
554 InEdgeIt e; |
542 for(InEdgeIt e(*g,w) ; e!=INVALID; ++e) { |
555 for(g->first(e,w); g->valid(e); g->next(e)) { |
543 |
556 |
|
557 if( (*flow)[e] <= 0 ) continue; |
544 if( (*flow)[e] <= 0 ) continue; |
558 Node v=g->tail(e); |
545 Node v=g->tail(e); |
559 |
546 |
560 if( lev > level[v] ) { //Push is allowed now |
547 if( lev > level[v] ) { //Push is allowed now |
561 |
548 |
582 } //for in edges vw |
569 } //for in edges vw |
583 |
570 |
584 } // if w still has excess after the out edge for cycle |
571 } // if w still has excess after the out edge for cycle |
585 |
572 |
586 excess.set(w, exc); |
573 excess.set(w, exc); |
587 |
574 |
588 return newlevel; |
575 return newlevel; |
589 } |
576 } |
590 |
577 |
591 |
578 |
592 |
579 |
593 void preflowPreproc(FlowEnum fe, NNMap& next, VecFirst& first, |
580 void preflowPreproc(FlowEnum fe, NNMap& next, VecFirst& first, |
594 VecNode& level_list, NNMap& left, NNMap& right) |
581 VecNode& level_list, NNMap& left, NNMap& right) |
595 { |
582 { |
596 switch (fe) { //setting excess |
583 switch (fe) { //setting excess |
597 case NO_FLOW: |
584 case NO_FLOW: |
585 for(EdgeIt e(*g); e!=INVALID; ++e) flow->set(e,0); |
|
586 for(NodeIt v(*g); v!=INVALID; ++v) excess.set(v,0); |
|
587 break; |
|
588 case ZERO_FLOW: |
|
589 for(NodeIt v(*g); v!=INVALID; ++v) excess.set(v,0); |
|
590 break; |
|
591 case GEN_FLOW: |
|
592 for(NodeIt v(*g); v!=INVALID; ++v) excess.set(v,0); |
|
598 { |
593 { |
599 EdgeIt e; |
|
600 for(g->first(e); g->valid(e); g->next(e)) flow->set(e,0); |
|
601 |
|
602 NodeIt v; |
|
603 for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0); |
|
604 break; |
|
605 } |
|
606 case ZERO_FLOW: |
|
607 { |
|
608 NodeIt v; |
|
609 for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0); |
|
610 break; |
|
611 } |
|
612 case GEN_FLOW: |
|
613 { |
|
614 NodeIt v; |
|
615 for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0); |
|
616 |
|
617 Num exc=0; |
594 Num exc=0; |
618 InEdgeIt e; |
595 for(InEdgeIt e(*g,t) ; e!=INVALID; ++e) exc+=(*flow)[e]; |
619 for(g->first(e,t); g->valid(e); g->next(e)) exc+=(*flow)[e]; |
596 for(OutEdgeIt e(*g,t) ; e!=INVALID; ++e) exc-=(*flow)[e]; |
620 OutEdgeIt f; |
|
621 for(g->first(f,t); g->valid(f); g->next(f)) exc-=(*flow)[f]; |
|
622 excess.set(t,exc); |
597 excess.set(t,exc); |
623 break; |
598 } |
624 } |
599 break; |
625 default: break; |
600 default: |
626 } |
601 break; |
627 |
602 } |
628 NodeIt v; |
603 |
629 for(g->first(v); g->valid(v); g->next(v)) level.set(v,n); |
604 for(NodeIt v(*g); v!=INVALID; ++v) level.set(v,n); |
630 //setting each node to level n |
605 //setting each node to level n |
631 |
606 |
632 std::queue<Node> bfs_queue; |
607 std::queue<Node> bfs_queue; |
633 |
608 |
634 |
609 |
635 switch (fe) { |
610 switch (fe) { |
636 case NO_FLOW: //flow is already set to const zero |
611 case NO_FLOW: //flow is already set to const zero |
637 case ZERO_FLOW: |
612 case ZERO_FLOW: |
638 { |
613 //Reverse_bfs from t, to find the starting level. |
639 //Reverse_bfs from t, to find the starting level. |
614 level.set(t,0); |
640 level.set(t,0); |
615 bfs_queue.push(t); |
641 bfs_queue.push(t); |
616 |
642 |
617 while (!bfs_queue.empty()) { |
643 while (!bfs_queue.empty()) { |
618 |
644 |
619 Node v=bfs_queue.front(); |
645 Node v=bfs_queue.front(); |
620 bfs_queue.pop(); |
646 bfs_queue.pop(); |
621 int l=level[v]+1; |
647 int l=level[v]+1; |
622 |
648 |
623 for(InEdgeIt e(*g,v) ; e!=INVALID; ++e) { |
649 InEdgeIt e; |
624 Node w=g->tail(e); |
650 for(g->first(e,v); g->valid(e); g->next(e)) { |
625 if ( level[w] == n && w != s ) { |
651 Node w=g->tail(e); |
626 bfs_queue.push(w); |
652 if ( level[w] == n && w != s ) { |
627 Node z=level_list[l]; |
653 bfs_queue.push(w); |
628 if ( z!=INVALID ) left.set(z,w); |
654 Node z=level_list[l]; |
629 right.set(w,z); |
655 if ( g->valid(z) ) left.set(z,w); |
630 level_list[l]=w; |
656 right.set(w,z); |
631 level.set(w, l); |
657 level_list[l]=w; |
632 } |
658 level.set(w, l); |
633 } |
659 } |
634 } |
660 } |
635 |
661 } |
636 //the starting flow |
662 |
637 for(OutEdgeIt e(*g,s) ; e!=INVALID; ++e) |
663 //the starting flow |
638 { |
664 OutEdgeIt e; |
639 Num c=(*capacity)[e]; |
665 for(g->first(e,s); g->valid(e); g->next(e)) |
640 if ( c <= 0 ) continue; |
641 Node w=g->head(e); |
|
642 if ( level[w] < n ) { |
|
643 if ( excess[w] <= 0 && w!=t ) //putting into the stack |
|
644 { |
|
645 next.set(w,first[level[w]]); |
|
646 first[level[w]]=w; |
|
647 } |
|
648 flow->set(e, c); |
|
649 excess.set(w, excess[w]+c); |
|
650 } |
|
651 } |
|
652 break; |
|
653 case GEN_FLOW: |
|
654 //Reverse_bfs from t in the residual graph, |
|
655 //to find the starting level. |
|
656 level.set(t,0); |
|
657 bfs_queue.push(t); |
|
658 |
|
659 while (!bfs_queue.empty()) { |
|
660 |
|
661 Node v=bfs_queue.front(); |
|
662 bfs_queue.pop(); |
|
663 int l=level[v]+1; |
|
664 |
|
665 for(InEdgeIt e(*g,v) ; e!=INVALID; ++e) { |
|
666 if ( (*capacity)[e] <= (*flow)[e] ) continue; |
|
667 Node w=g->tail(e); |
|
668 if ( level[w] == n && w != s ) { |
|
669 bfs_queue.push(w); |
|
670 Node z=level_list[l]; |
|
671 if ( z!=INVALID ) left.set(z,w); |
|
672 right.set(w,z); |
|
673 level_list[l]=w; |
|
674 level.set(w, l); |
|
675 } |
|
676 } |
|
677 |
|
678 for(OutEdgeIt e(*g,v) ; e!=INVALID; ++e) { |
|
679 if ( 0 >= (*flow)[e] ) continue; |
|
680 Node w=g->head(e); |
|
681 if ( level[w] == n && w != s ) { |
|
682 bfs_queue.push(w); |
|
683 Node z=level_list[l]; |
|
684 if ( z!=INVALID ) left.set(z,w); |
|
685 right.set(w,z); |
|
686 level_list[l]=w; |
|
687 level.set(w, l); |
|
688 } |
|
689 } |
|
690 } |
|
691 |
|
692 //the starting flow |
|
693 for(OutEdgeIt e(*g,s); e!=INVALID; ++e) |
|
694 { |
|
695 Num rem=(*capacity)[e]-(*flow)[e]; |
|
696 if ( rem <= 0 ) continue; |
|
697 Node w=g->head(e); |
|
698 if ( level[w] < n ) { |
|
699 if ( excess[w] <= 0 && w!=t ) //putting into the stack |
|
700 { |
|
701 next.set(w,first[level[w]]); |
|
702 first[level[w]]=w; |
|
703 } |
|
704 flow->set(e, (*capacity)[e]); |
|
705 excess.set(w, excess[w]+rem); |
|
706 } |
|
707 } |
|
708 |
|
709 for(InEdgeIt e(*g,s); e!=INVALID; ++e) |
|
710 { |
|
711 if ( (*flow)[e] <= 0 ) continue; |
|
712 Node w=g->tail(e); |
|
713 if ( level[w] < n ) { |
|
714 if ( excess[w] <= 0 && w!=t ) |
|
715 { |
|
716 next.set(w,first[level[w]]); |
|
717 first[level[w]]=w; |
|
718 } |
|
719 excess.set(w, excess[w]+(*flow)[e]); |
|
720 flow->set(e, 0); |
|
721 } |
|
722 } |
|
723 break; |
|
724 case PRE_FLOW: |
|
725 //Reverse_bfs from t in the residual graph, |
|
726 //to find the starting level. |
|
727 level.set(t,0); |
|
728 bfs_queue.push(t); |
|
729 |
|
730 while (!bfs_queue.empty()) { |
|
731 |
|
732 Node v=bfs_queue.front(); |
|
733 bfs_queue.pop(); |
|
734 int l=level[v]+1; |
|
735 |
|
736 for(InEdgeIt e(*g,v) ; e!=INVALID; ++e) { |
|
737 if ( (*capacity)[e] <= (*flow)[e] ) continue; |
|
738 Node w=g->tail(e); |
|
739 if ( level[w] == n && w != s ) { |
|
740 bfs_queue.push(w); |
|
741 Node z=level_list[l]; |
|
742 if ( z!=INVALID ) left.set(z,w); |
|
743 right.set(w,z); |
|
744 level_list[l]=w; |
|
745 level.set(w, l); |
|
746 } |
|
747 } |
|
748 |
|
749 for(OutEdgeIt e(*g,v) ; e!=INVALID; ++e) { |
|
750 if ( 0 >= (*flow)[e] ) continue; |
|
751 Node w=g->head(e); |
|
752 if ( level[w] == n && w != s ) { |
|
753 bfs_queue.push(w); |
|
754 Node z=level_list[l]; |
|
755 if ( z!=INVALID ) left.set(z,w); |
|
756 right.set(w,z); |
|
757 level_list[l]=w; |
|
758 level.set(w, l); |
|
759 } |
|
760 } |
|
761 } |
|
762 |
|
763 |
|
764 //the starting flow |
|
765 for(OutEdgeIt e(*g,s) ; e!=INVALID; ++e) { |
|
766 Num rem=(*capacity)[e]-(*flow)[e]; |
|
767 if ( rem <= 0 ) continue; |
|
768 Node w=g->head(e); |
|
769 if ( level[w] < n ) { |
|
770 flow->set(e, (*capacity)[e]); |
|
771 excess.set(w, excess[w]+rem); |
|
772 } |
|
773 } |
|
774 |
|
775 for(InEdgeIt e(*g,s) ; e!=INVALID; ++e) { |
|
776 if ( (*flow)[e] <= 0 ) continue; |
|
777 Node w=g->tail(e); |
|
778 if ( level[w] < n ) { |
|
779 excess.set(w, excess[w]+(*flow)[e]); |
|
780 flow->set(e, 0); |
|
781 } |
|
782 } |
|
783 |
|
784 //computing the excess |
|
785 for(NodeIt w(*g); w!=INVALID; ++w) { |
|
786 Num exc=0; |
|
787 |
|
788 for(InEdgeIt e(*g,w) ; e!=INVALID; ++e) exc+=(*flow)[e]; |
|
789 for(OutEdgeIt e(*g,w) ; e!=INVALID; ++e) exc-=(*flow)[e]; |
|
790 |
|
791 excess.set(w,exc); |
|
792 |
|
793 //putting the active nodes into the stack |
|
794 int lev=level[w]; |
|
795 if ( exc > 0 && lev < n && Node(w) != t ) |
|
796 ///\bug if ( exc > 0 && lev < n && w != t ) temporarily for working with wrappers. |
|
666 { |
797 { |
667 Num c=(*capacity)[e]; |
798 next.set(w,first[lev]); |
668 if ( c <= 0 ) continue; |
799 first[lev]=w; |
669 Node w=g->head(e); |
800 } |
670 if ( level[w] < n ) { |
801 } |
671 if ( excess[w] <= 0 && w!=t ) //putting into the stack |
802 break; |
672 { |
803 } //switch |
673 next.set(w,first[level[w]]); |
|
674 first[level[w]]=w; |
|
675 } |
|
676 flow->set(e, c); |
|
677 excess.set(w, excess[w]+c); |
|
678 } |
|
679 } |
|
680 break; |
|
681 } |
|
682 |
|
683 case GEN_FLOW: |
|
684 { |
|
685 //Reverse_bfs from t in the residual graph, |
|
686 //to find the starting level. |
|
687 level.set(t,0); |
|
688 bfs_queue.push(t); |
|
689 |
|
690 while (!bfs_queue.empty()) { |
|
691 |
|
692 Node v=bfs_queue.front(); |
|
693 bfs_queue.pop(); |
|
694 int l=level[v]+1; |
|
695 |
|
696 InEdgeIt e; |
|
697 for(g->first(e,v); g->valid(e); g->next(e)) { |
|
698 if ( (*capacity)[e] <= (*flow)[e] ) continue; |
|
699 Node w=g->tail(e); |
|
700 if ( level[w] == n && w != s ) { |
|
701 bfs_queue.push(w); |
|
702 Node z=level_list[l]; |
|
703 if ( g->valid(z) ) left.set(z,w); |
|
704 right.set(w,z); |
|
705 level_list[l]=w; |
|
706 level.set(w, l); |
|
707 } |
|
708 } |
|
709 |
|
710 OutEdgeIt f; |
|
711 for(g->first(f,v); g->valid(f); g->next(f)) { |
|
712 if ( 0 >= (*flow)[f] ) continue; |
|
713 Node w=g->head(f); |
|
714 if ( level[w] == n && w != s ) { |
|
715 bfs_queue.push(w); |
|
716 Node z=level_list[l]; |
|
717 if ( g->valid(z) ) left.set(z,w); |
|
718 right.set(w,z); |
|
719 level_list[l]=w; |
|
720 level.set(w, l); |
|
721 } |
|
722 } |
|
723 } |
|
724 |
|
725 //the starting flow |
|
726 OutEdgeIt e; |
|
727 for(g->first(e,s); g->valid(e); g->next(e)) |
|
728 { |
|
729 Num rem=(*capacity)[e]-(*flow)[e]; |
|
730 if ( rem <= 0 ) continue; |
|
731 Node w=g->head(e); |
|
732 if ( level[w] < n ) { |
|
733 if ( excess[w] <= 0 && w!=t ) //putting into the stack |
|
734 { |
|
735 next.set(w,first[level[w]]); |
|
736 first[level[w]]=w; |
|
737 } |
|
738 flow->set(e, (*capacity)[e]); |
|
739 excess.set(w, excess[w]+rem); |
|
740 } |
|
741 } |
|
742 |
|
743 InEdgeIt f; |
|
744 for(g->first(f,s); g->valid(f); g->next(f)) |
|
745 { |
|
746 if ( (*flow)[f] <= 0 ) continue; |
|
747 Node w=g->tail(f); |
|
748 if ( level[w] < n ) { |
|
749 if ( excess[w] <= 0 && w!=t ) |
|
750 { |
|
751 next.set(w,first[level[w]]); |
|
752 first[level[w]]=w; |
|
753 } |
|
754 excess.set(w, excess[w]+(*flow)[f]); |
|
755 flow->set(f, 0); |
|
756 } |
|
757 } |
|
758 break; |
|
759 } //case GEN_FLOW |
|
760 |
|
761 case PRE_FLOW: |
|
762 { |
|
763 //Reverse_bfs from t in the residual graph, |
|
764 //to find the starting level. |
|
765 level.set(t,0); |
|
766 bfs_queue.push(t); |
|
767 |
|
768 while (!bfs_queue.empty()) { |
|
769 |
|
770 Node v=bfs_queue.front(); |
|
771 bfs_queue.pop(); |
|
772 int l=level[v]+1; |
|
773 |
|
774 InEdgeIt e; |
|
775 for(g->first(e,v); g->valid(e); g->next(e)) { |
|
776 if ( (*capacity)[e] <= (*flow)[e] ) continue; |
|
777 Node w=g->tail(e); |
|
778 if ( level[w] == n && w != s ) { |
|
779 bfs_queue.push(w); |
|
780 Node z=level_list[l]; |
|
781 if ( g->valid(z) ) left.set(z,w); |
|
782 right.set(w,z); |
|
783 level_list[l]=w; |
|
784 level.set(w, l); |
|
785 } |
|
786 } |
|
787 |
|
788 OutEdgeIt f; |
|
789 for(g->first(f,v); g->valid(f); g->next(f)) { |
|
790 if ( 0 >= (*flow)[f] ) continue; |
|
791 Node w=g->head(f); |
|
792 if ( level[w] == n && w != s ) { |
|
793 bfs_queue.push(w); |
|
794 Node z=level_list[l]; |
|
795 if ( g->valid(z) ) left.set(z,w); |
|
796 right.set(w,z); |
|
797 level_list[l]=w; |
|
798 level.set(w, l); |
|
799 } |
|
800 } |
|
801 } |
|
802 |
|
803 |
|
804 //the starting flow |
|
805 OutEdgeIt e; |
|
806 for(g->first(e,s); g->valid(e); g->next(e)) |
|
807 { |
|
808 Num rem=(*capacity)[e]-(*flow)[e]; |
|
809 if ( rem <= 0 ) continue; |
|
810 Node w=g->head(e); |
|
811 if ( level[w] < n ) { |
|
812 flow->set(e, (*capacity)[e]); |
|
813 excess.set(w, excess[w]+rem); |
|
814 } |
|
815 } |
|
816 |
|
817 InEdgeIt f; |
|
818 for(g->first(f,s); g->valid(f); g->next(f)) |
|
819 { |
|
820 if ( (*flow)[f] <= 0 ) continue; |
|
821 Node w=g->tail(f); |
|
822 if ( level[w] < n ) { |
|
823 excess.set(w, excess[w]+(*flow)[f]); |
|
824 flow->set(f, 0); |
|
825 } |
|
826 } |
|
827 |
|
828 NodeIt w; //computing the excess |
|
829 for(g->first(w); g->valid(w); g->next(w)) { |
|
830 Num exc=0; |
|
831 |
|
832 InEdgeIt e; |
|
833 for(g->first(e,w); g->valid(e); g->next(e)) exc+=(*flow)[e]; |
|
834 OutEdgeIt f; |
|
835 for(g->first(f,w); g->valid(f); g->next(f)) exc-=(*flow)[f]; |
|
836 |
|
837 excess.set(w,exc); |
|
838 |
|
839 //putting the active nodes into the stack |
|
840 int lev=level[w]; |
|
841 if ( exc > 0 && lev < n && Node(w) != t ) |
|
842 ///\bug if ( exc > 0 && lev < n && w != t ) tempararily for working with wrappers. in hugo 0.2 it will work. Amugy mukodik sage_graph-fal, de smart_graph-fal nem, azt hozzatennem. |
|
843 { |
|
844 next.set(w,first[lev]); |
|
845 first[lev]=w; |
|
846 } |
|
847 } |
|
848 break; |
|
849 } //case PRE_FLOW |
|
850 } |
|
851 } //preflowPreproc |
804 } //preflowPreproc |
852 |
805 |
853 |
806 |
854 void relabel(Node w, int newlevel, NNMap& next, VecFirst& first, |
807 void relabel(Node w, int newlevel, NNMap& next, VecFirst& first, |
855 VecNode& level_list, NNMap& left, |
808 VecNode& level_list, NNMap& left, |
860 |
813 |
861 Node right_n=right[w]; |
814 Node right_n=right[w]; |
862 Node left_n=left[w]; |
815 Node left_n=left[w]; |
863 |
816 |
864 //unlacing starts |
817 //unlacing starts |
865 if ( g->valid(right_n) ) { |
818 if ( right_n!=INVALID ) { |
866 if ( g->valid(left_n) ) { |
819 if ( left_n!=INVALID ) { |
867 right.set(left_n, right_n); |
820 right.set(left_n, right_n); |
868 left.set(right_n, left_n); |
821 left.set(right_n, left_n); |
869 } else { |
822 } else { |
870 level_list[lev]=right_n; |
823 level_list[lev]=right_n; |
871 left.set(right_n, INVALID); |
824 left.set(right_n, INVALID); |
872 } |
825 } |
873 } else { |
826 } else { |
874 if ( g->valid(left_n) ) { |
827 if ( left_n!=INVALID ) { |
875 right.set(left_n, INVALID); |
828 right.set(left_n, INVALID); |
876 } else { |
829 } else { |
877 level_list[lev]=INVALID; |
830 level_list[lev]=INVALID; |
878 } |
831 } |
879 } |
832 } |
880 //unlacing ends |
833 //unlacing ends |
881 |
834 |
882 if ( !g->valid(level_list[lev]) ) { |
835 if ( level_list[lev]==INVALID ) { |
883 |
836 |
884 //gapping starts |
837 //gapping starts |
885 for (int i=lev; i!=k ; ) { |
838 for (int i=lev; i!=k ; ) { |
886 Node v=level_list[++i]; |
839 Node v=level_list[++i]; |
887 while ( g->valid(v) ) { |
840 while ( v!=INVALID ) { |
888 level.set(v,n); |
841 level.set(v,n); |
889 v=right[v]; |
842 v=right[v]; |
890 } |
843 } |
891 level_list[i]=INVALID; |
844 level_list[i]=INVALID; |
892 if ( !what_heur ) first[i]=INVALID; |
845 if ( !what_heur ) first[i]=INVALID; |
905 next.set(w,first[newlevel]); |
858 next.set(w,first[newlevel]); |
906 first[newlevel]=w; |
859 first[newlevel]=w; |
907 if ( what_heur ) b=newlevel; |
860 if ( what_heur ) b=newlevel; |
908 if ( k < newlevel ) ++k; //now k=newlevel |
861 if ( k < newlevel ) ++k; //now k=newlevel |
909 Node z=level_list[newlevel]; |
862 Node z=level_list[newlevel]; |
910 if ( g->valid(z) ) left.set(z,w); |
863 if ( z!=INVALID ) left.set(z,w); |
911 right.set(w,z); |
864 right.set(w,z); |
912 left.set(w,INVALID); |
865 left.set(w,INVALID); |
913 level_list[newlevel]=w; |
866 level_list[newlevel]=w; |
914 } |
867 } |
915 } |
868 } |
916 } //relabel |
869 } //relabel |
917 |
870 |
918 void printexcess() {//// |
871 void printexcess() {//// |
919 std::cout << "Excesses:" <<std::endl; |
872 std::cout << "Excesses:" <<std::endl; |
920 |
873 |
921 NodeIt v; |
874 for(NodeIt v(*g); v!=INVALID ; ++v) { |
922 for(g->first(v); g->valid(v); g->next(v)) { |
|
923 std::cout << 1+(g->id(v)) << ":" << excess[v]<<std::endl; |
875 std::cout << 1+(g->id(v)) << ":" << excess[v]<<std::endl; |
924 } |
876 } |
925 } |
877 } |
926 |
878 |
927 void printlevel() {//// |
879 void printlevel() {//// |
928 std::cout << "Levels:" <<std::endl; |
880 std::cout << "Levels:" <<std::endl; |
929 |
881 |
930 NodeIt v; |
882 for(NodeIt v(*g); v!=INVALID ; ++v) { |
931 for(g->first(v); g->valid(v); g->next(v)) { |
|
932 std::cout << 1+(g->id(v)) << ":" << level[v]<<std::endl; |
883 std::cout << 1+(g->id(v)) << ":" << level[v]<<std::endl; |
933 } |
884 } |
934 } |
885 } |
935 |
886 |
936 void printactive() {//// |
887 void printactive() {//// |
937 std::cout << "Levels:" <<std::endl; |
888 std::cout << "Levels:" <<std::endl; |
938 |
889 |
939 NodeIt v; |
890 for(NodeIt v(*g); v!=INVALID ; ++v) { |
940 for(g->first(v); g->valid(v); g->next(v)) { |
|
941 std::cout << 1+(g->id(v)) << ":" << level[v]<<std::endl; |
891 std::cout << 1+(g->id(v)) << ":" << level[v]<<std::endl; |
942 } |
892 } |
943 } |
893 } |
944 |
894 |
945 |
895 |