Changeset 615:b6b31b75b522 in lemon0.x for src/work/jacint
 Timestamp:
 05/11/04 21:50:21 (17 years ago)
 Branch:
 default
 Phase:
 public
 Convert:
 svn:c9d7d8f590d60310b91f818b3a526b0e/lemon/trunk@800
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

src/work/jacint/max_flow.h
r602 r615 2 2 3 3 /* 4 Heuristics: 4 Heuristics: 5 5 2 phase 6 6 gap … … 9 9 runs heuristic 'highest label' for H1*n relabels 10 10 runs heuristic 'bound decrease' for H0*n relabels, starts with 'highest label' 11 11 12 12 Parameters H0 and H1 are initialized to 20 and 1. 13 13 14 14 Constructors: 15 15 16 Preflow(Graph, Node, Node, CapMap, FlowMap, bool) : bool must be false if 16 Preflow(Graph, Node, Node, CapMap, FlowMap, bool) : bool must be false if 17 17 FlowMap is not constant zero, and should be true if it is 18 18 … … 23 23 Num flowValue() : returns the value of a maximum flow 24 24 25 void minMinCut(CutMap& M) : sets M to the characteristic vector of the 25 void minMinCut(CutMap& M) : sets M to the characteristic vector of the 26 26 minimum min cut. M should be a map of bools initialized to false. ??Is it OK? 27 27 28 void maxMinCut(CutMap& M) : sets M to the characteristic vector of the 28 void maxMinCut(CutMap& M) : sets M to the characteristic vector of the 29 29 maximum min cut. M should be a map of bools initialized to false. 30 30 31 void minCut(CutMap& M) : sets M to the characteristic vector of 31 void minCut(CutMap& M) : sets M to the characteristic vector of 32 32 a min cut. M should be a map of bools initialized to false. 33 33 … … 36 36 #ifndef HUGO_MAX_FLOW_H 37 37 #define HUGO_MAX_FLOW_H 38 39 #define H0 2040 #define H1 141 38 42 39 #include <vector> … … 51 48 52 49 /// \file 53 /// \brief Dimacs file format reader. 50 /// \brief Maximum flows. 51 /// \ingroup galgs 54 52 55 53 namespace hugo { … … 58 56 // ///\author Marton Makai, Jacint Szabo 59 57 /// A class for computing max flows and related quantities. 60 template <typename Graph, typename Num, 61 typename CapMap=typename Graph::template EdgeMap<Num>, 58 /// \ingroup galgs 59 template <typename Graph, typename Num, 60 typename CapMap=typename Graph::template EdgeMap<Num>, 62 61 typename FlowMap=typename Graph::template EdgeMap<Num> > 63 62 class MaxFlow { 64 63 protected: 65 64 typedef typename Graph::Node Node; 66 65 typedef typename Graph::NodeIt NodeIt; … … 75 74 Node s; 76 75 Node t; 77 const CapMap* capacity; 76 const CapMap* capacity; 78 77 FlowMap* flow; 79 78 int n; //the number of nodes of G … … 84 83 typedef typename Graph::template NodeMap<int> ReachedMap; 85 84 ReachedMap level; 86 //level works as a bool map in augmenting path algorithms 85 //level works as a bool map in augmenting path algorithms 87 86 //and is used by bfs for storing reached information. 88 87 //In preflow, it shows levels of nodes. 89 //typename Graph::template NodeMap<int> level; 90 typename Graph::template NodeMap<Num> excess; 88 //typename Graph::template NodeMap<int> level; 89 typename Graph::template NodeMap<Num> excess; 91 90 // protected: 92 91 // MaxFlow() { } 93 // void set(const Graph& _G, Node _s, Node _t, const CapMap& _capacity, 94 // FlowMap& _flow) 92 // void set(const Graph& _G, Node _s, Node _t, const CapMap& _capacity, 93 // FlowMap& _flow) 95 94 // { 96 // g=&_G; 97 // s=_s; 98 // t=_t; 95 // g=&_G; 96 // s=_s; 97 // t=_t; 99 98 // capacity=&_capacity; 100 99 // flow=&_flow; 101 100 // n=_G.nodeNum; 102 // level.set (_G); //kellene vmi ilyesmi fv 101 // level.set (_G); //kellene vmi ilyesmi fv 103 102 // excess(_G,0); //itt is 104 103 // } 105 104 105 // constants used for heuristics 106 static const int H0=20; 107 static const int H1=1; 108 106 109 public: 107 110 108 111 ///\todo Document this. 109 112 ///\todo Maybe, it should be PRE_FLOW instead. 113 /// \c NO_FLOW means nothing, 110 114 /// \c ZERO_FLOW means something, 111 115 /// \c GEN_FLOW means something else, 112 /// \c PRE FLOW is something different.116 /// \c PRE_FLOW is something different. 113 117 enum flowEnum{ 114 ZERO_FLOW=0, 115 GEN_FLOW=1, 116 PREFLOW=2 118 ZERO_FLOW, 119 GEN_FLOW, 120 PRE_FLOW, 121 NO_FLOW 117 122 }; 118 123 119 MaxFlow(const Graph& _G, Node _s, Node _t, const CapMap& _capacity, 124 MaxFlow(const Graph& _G, Node _s, Node _t, const CapMap& _capacity, 120 125 FlowMap& _flow) : 121 g(&_G), s(_s), t(_t), capacity(&_capacity), 126 g(&_G), s(_s), t(_t), capacity(&_capacity), 122 127 flow(&_flow), n(_G.nodeNum()), level(_G), excess(_G,0) {} 123 128 124 129 /// A max flow algorithm is run. 125 ///\pre the flow have to be 0 at the beginning. 126 void run() { 127 preflow(ZERO_FLOW); 128 } 129 130 /// A preflow algorithm is run. 131 ///\pre The initial edgemap have to be a 130 /// \pre The flow have to satisfy the requirements 131 /// stated in fe. 132 void run(flowEnum fe=ZERO_FLOW) { 133 preflow(fe); 134 } 135 136 /// A preflow algorithm is run. 137 /// \pre The initial edgemap have to be a 132 138 /// zero flow if \c fe is \c ZERO_FLOW, 133 /// a flow if \c fe is \c GEN_FLOW, 134 /// and a preflow it is \c PREFLOW. 139 /// a flow if \c fe is \c GEN_FLOW, 140 /// a preflow if fe is \c PRE_FLOW and 141 /// anything if fe is NO_FLOW. 135 142 void preflow(flowEnum fe) { 136 143 preflowPhase0(fe); … … 138 145 } 139 146 140 /// Run the first phase of preflow, starting from a 0 flow, from a flow, 141 /// or from a preflow, according to \c fe.142 void preflowPhase0( flowEnum fe);147 /// Run the first phase of preflow, starting from a 0 flow, from a flow, 148 /// or from a preflow, of from undefined value according to \c fe. 149 void preflowPhase0(flowEnum fe); 143 150 144 151 /// Second phase of preflow. 145 152 void preflowPhase1(); 146 153 147 /// Starting from a flow, this method searches for an augmenting path 148 /// according to the EdmondsKarp algorithm 149 /// and augments the flow on if any. 154 /// Starting from a flow, this method searches for an augmenting path 155 /// according to the EdmondsKarp algorithm 156 /// and augments the flow on if any. 150 157 /// The return value shows if the augmentation was succesful. 151 158 bool augmentOnShortestPath(); 152 159 153 /// Starting from a flow, this method searches for an augmenting blockin 154 /// flow according to Dinits' algorithm and augments the flow on if any. 155 /// The blocking flow is computed in a physically constructed 160 /// Starting from a flow, this method searches for an augmenting blocking 161 /// flow according to Dinits' algorithm and augments the flow on if any. 162 /// The blocking flow is computed in a physically constructed 156 163 /// residual graph of type \c Mutablegraph. 157 164 /// The return value show sif the augmentation was succesful. 158 165 template<typename MutableGraph> bool augmentOnBlockingFlow(); 159 166 160 /// The same as \c augmentOnBlockingFlow<MutableGraph> but the 167 /// The same as \c augmentOnBlockingFlow<MutableGraph> but the 161 168 /// residual graph is not constructed physically. 162 169 /// The return value shows if the augmentation was succesful. … … 164 171 165 172 /// Returns the actual flow value. 166 /// More precisely, it returns the negative excess of s, thus 173 /// More precisely, it returns the negative excess of s, thus 167 174 /// this works also for preflows. 168 Num flowValue() { 175 Num flowValue() { 169 176 Num a=0; 170 FOR_EACH_INC_LOC( OutEdgeIt, e, *g, s) a+=(*flow)[e];171 FOR_EACH_INC_LOC( InEdgeIt, e, *g, s) a=(*flow)[e];177 FOR_EACH_INC_LOC(InEdgeIt, e, *g, t) a+=(*flow)[e]; 178 FOR_EACH_INC_LOC(OutEdgeIt, e, *g, t) a=(*flow)[e]; 172 179 return a; 173 180 } 174 181 175 182 /// Should be used between preflowPhase0 and preflowPhase1. 176 ///\todo We have to make some status variable which shows the actual state 177 /// of the class. This enables us to determine which methods are valid 183 /// \todo We have to make some status variable which shows the 184 /// actual state 185 /// of the class. This enables us to determine which methods are valid 178 186 /// for MinCut computation 179 187 template<typename _CutMap> … … 189 197 } 190 198 191 /// The unique inclusionwise minimum cut is computed by 199 /// The unique inclusionwise minimum cut is computed by 192 200 /// processing a bfs from s in the residual graph. 193 /// \pre flow have to be a max flow otherwise it will the whole nodeset.201 /// \pre flow have to be a max flow otherwise it will the whole nodeset. 194 202 template<typename _CutMap> 195 203 void minMinCut(_CutMap& M) { 196 204 197 205 std::queue<Node> queue; 198 199 M.set(s,true); 206 207 M.set(s,true); 200 208 queue.push(s); 201 209 … … 211 219 M.set(v, true); 212 220 } 213 } 221 } 214 222 215 223 InEdgeIt f; … … 220 228 M.set(v, true); 221 229 } 222 } 223 } 224 } 225 226 227 /// The unique inclusionwise maximum cut is computed by 230 } 231 } 232 } 233 234 235 /// The unique inclusionwise maximum cut is computed by 228 236 /// processing a reverse bfs from t in the residual graph. 229 /// \pre flow have to be a max flow otherwise it will be empty.237 /// \pre flow have to be a max flow otherwise it will be empty. 230 238 template<typename _CutMap> 231 239 void maxMinCut(_CutMap& M) { … … 237 245 238 246 std::queue<Node> queue; 239 240 M.set(t,false); 247 248 M.set(t,false); 241 249 queue.push(t); 242 250 … … 244 252 Node w=queue.front(); 245 253 queue.pop(); 246 247 254 248 255 InEdgeIt e; … … 254 261 } 255 262 } 256 263 257 264 OutEdgeIt f; 258 265 for(g>first(f,w) ; g>valid(f); g>next(f)) { … … 275 282 /// 276 283 void resetTarget(Node _t) { t=_t; } 277 284 278 285 /// capacitymap is changed. 279 286 void resetCap(const CapMap& _cap) { capacity=&_cap; } 280 281 /// flowmap is changed. 287 288 /// flowmap is changed. 282 289 void resetFlow(FlowMap& _flow) { flow=&_flow; } 283 290 … … 286 293 287 294 int push(Node w, VecStack& active) { 288 295 289 296 int lev=level[w]; 290 297 Num exc=excess[w]; 291 298 int newlevel=n; //bound on the next level of w 292 299 293 300 OutEdgeIt e; 294 301 for(g>first(e,w); g>valid(e); g>next(e)) { 295 296 if ( (*flow)[e] >= (*capacity)[e] ) continue; 297 Node v=g>head(e); 298 302 303 if ( (*flow)[e] >= (*capacity)[e] ) continue; 304 Node v=g>head(e); 305 299 306 if( lev > level[v] ) { //Push is allowed now 300 307 301 308 if ( excess[v]<=0 && v!=t && v!=s ) { 302 309 int lev_v=level[v]; 303 310 active[lev_v].push(v); 304 311 } 305 312 306 313 Num cap=(*capacity)[e]; 307 314 Num flo=(*flow)[e]; 308 315 Num remcap=capflo; 309 316 310 317 if ( remcap >= exc ) { //A nonsaturating push. 311 318 312 319 flow>set(e, flo+exc); 313 320 excess.set(v, excess[v]+exc); 314 321 exc=0; 315 break; 316 322 break; 323 317 324 } else { //A saturating push. 318 325 flow>set(e, cap); … … 321 328 } 322 329 } else if ( newlevel > level[v] ) newlevel = level[v]; 323 } //for out edges wv 324 325 if ( exc > 0 ) { 330 } //for out edges wv 331 332 if ( exc > 0 ) { 326 333 InEdgeIt e; 327 334 for(g>first(e,w); g>valid(e); g>next(e)) { 328 329 if( (*flow)[e] <= 0 ) continue; 330 Node v=g>tail(e); 331 335 336 if( (*flow)[e] <= 0 ) continue; 337 Node v=g>tail(e); 338 332 339 if( lev > level[v] ) { //Push is allowed now 333 340 334 341 if ( excess[v]<=0 && v!=t && v!=s ) { 335 342 int lev_v=level[v]; 336 343 active[lev_v].push(v); 337 344 } 338 345 339 346 Num flo=(*flow)[e]; 340 347 341 348 if ( flo >= exc ) { //A nonsaturating push. 342 349 343 350 flow>set(e, floexc); 344 351 excess.set(v, excess[v]+exc); 345 352 exc=0; 346 break; 353 break; 347 354 } else { //A saturating push. 348 355 349 356 excess.set(v, excess[v]+flo); 350 357 exc=flo; 351 358 flow>set(e,0); 352 } 359 } 353 360 } else if ( newlevel > level[v] ) newlevel = level[v]; 354 361 } //for in edges vw 355 362 356 363 } // if w still has excess after the out edge for cycle 357 364 358 365 excess.set(w, exc); 359 366 360 367 return newlevel; 361 368 } 362 369 363 370 364 void preflowPreproc ( flowEnum fe, VecStack& active,365 VecNode& level_list, NNMap& left, NNMap& right )371 void preflowPreproc(flowEnum fe, VecStack& active, 372 VecNode& level_list, NNMap& left, NNMap& right) 366 373 { 367 368 374 std::queue<Node> bfs_queue; 369 370 switch ( fe) {371 case ZERO_FLOW: 375 376 switch (fe) { 377 case ZERO_FLOW: 372 378 { 373 379 //Reverse_bfs from t, to find the starting level. 374 380 level.set(t,0); 375 381 bfs_queue.push(t); 376 382 377 383 while (!bfs_queue.empty()) { 378 379 Node v=bfs_queue.front(); 384 385 Node v=bfs_queue.front(); 380 386 bfs_queue.pop(); 381 387 int l=level[v]+1; 382 388 383 389 InEdgeIt e; 384 390 for(g>first(e,v); g>valid(e); g>next(e)) { … … 394 400 } 395 401 } 396 402 397 403 //the starting flow 398 404 OutEdgeIt e; 399 for(g>first(e,s); g>valid(e); g>next(e)) 405 for(g>first(e,s); g>valid(e); g>next(e)) 400 406 { 401 407 Num c=(*capacity)[e]; 402 408 if ( c <= 0 ) continue; 403 409 Node w=g>head(e); 404 if ( level[w] < n ) { 410 if ( level[w] < n ) { 405 411 if ( excess[w] <= 0 && w!=t ) active[level[w]].push(w); 406 flow>set(e, c); 412 flow>set(e, c); 407 413 excess.set(w, excess[w]+c); 408 414 } … … 410 416 break; 411 417 } 412 418 413 419 case GEN_FLOW: 414 case PRE FLOW:420 case PRE_FLOW: 415 421 { 416 //Reverse_bfs from t in the residual graph, 422 //Reverse_bfs from t in the residual graph, 417 423 //to find the starting level. 418 424 level.set(t,0); 419 425 bfs_queue.push(t); 420 426 421 427 while (!bfs_queue.empty()) { 422 423 Node v=bfs_queue.front(); 428 429 Node v=bfs_queue.front(); 424 430 bfs_queue.pop(); 425 431 int l=level[v]+1; 426 432 427 433 InEdgeIt e; 428 434 for(g>first(e,v); g>valid(e); g>next(e)) { … … 438 444 } 439 445 } 440 446 441 447 OutEdgeIt f; 442 448 for(g>first(f,v); g>valid(f); g>next(f)) { … … 453 459 } 454 460 } 455 456 461 462 457 463 //the starting flow 458 464 OutEdgeIt e; 459 for(g>first(e,s); g>valid(e); g>next(e)) 465 for(g>first(e,s); g>valid(e); g>next(e)) 460 466 { 461 467 Num rem=(*capacity)[e](*flow)[e]; 462 468 if ( rem <= 0 ) continue; 463 469 Node w=g>head(e); 464 if ( level[w] < n ) { 470 if ( level[w] < n ) { 465 471 if ( excess[w] <= 0 && w!=t ) active[level[w]].push(w); 466 flow>set(e, (*capacity)[e]); 472 flow>set(e, (*capacity)[e]); 467 473 excess.set(w, excess[w]+rem); 468 474 } 469 475 } 470 476 471 477 InEdgeIt f; 472 for(g>first(f,s); g>valid(f); g>next(f)) 478 for(g>first(f,s); g>valid(f); g>next(f)) 473 479 { 474 480 if ( (*flow)[f] <= 0 ) continue; 475 481 Node w=g>tail(f); 476 if ( level[w] < n ) { 482 if ( level[w] < n ) { 477 483 if ( excess[w] <= 0 && w!=t ) active[level[w]].push(w); 478 484 excess.set(w, excess[w]+(*flow)[f]); 479 flow>set(f, 0); 485 flow>set(f, 0); 480 486 } 481 } 487 } 482 488 break; 483 } //case PRE FLOW489 } //case PRE_FLOW 484 490 } 485 491 } //preflowPreproc … … 487 493 488 494 489 void relabel(Node w, int newlevel, VecStack& active, 490 VecNode& level_list, NNMap& left, 491 NNMap& right, int& b, int& k, bool what_heur ) 495 void relabel(Node w, int newlevel, VecStack& active, 496 VecNode& level_list, NNMap& left, 497 NNMap& right, int& b, int& k, bool what_heur ) 492 498 { 493 499 494 Num lev=level[w]; 495 500 Num lev=level[w]; 501 496 502 Node right_n=right[w]; 497 503 Node left_n=left[w]; 498 504 499 505 //unlacing starts 500 506 if ( g>valid(right_n) ) { … … 503 509 left.set(right_n, left_n); 504 510 } else { 505 level_list[lev]=right_n; 511 level_list[lev]=right_n; 506 512 left.set(right_n, INVALID); 507 } 513 } 508 514 } else { 509 515 if ( g>valid(left_n) ) { 510 516 right.set(left_n, INVALID); 511 } else { 512 level_list[lev]=INVALID; 513 } 514 } 517 } else { 518 level_list[lev]=INVALID; 519 } 520 } 515 521 //unlacing ends 516 522 517 523 if ( !g>valid(level_list[lev]) ) { 518 524 519 525 //gapping starts 520 526 for (int i=lev; i!=k ; ) { … … 529 535 active[i].pop(); //FIXME: ezt szebben kene 530 536 } 531 } 532 } 533 537 } 538 } 539 534 540 level.set(w,n); 535 541 b=lev1; 536 542 k=b; 537 543 //gapping ends 538 544 539 545 } else { 540 541 if ( newlevel == n ) level.set(w,n); 546 547 if ( newlevel == n ) level.set(w,n); 542 548 else { 543 549 level.set(w,++newlevel); … … 552 558 } 553 559 } 554 560 555 561 } //relabel 556 562 557 563 558 template<typename MapGraphWrapper> 564 template<typename MapGraphWrapper> 559 565 class DistanceMap { 560 566 protected: 561 567 const MapGraphWrapper* g; 562 typename MapGraphWrapper::template NodeMap<int> dist; 568 typename MapGraphWrapper::template NodeMap<int> dist; 563 569 public: 564 570 DistanceMap(MapGraphWrapper& _g) : g(&_g), dist(*g, g>nodeNum()) { } 565 void set(const typename MapGraphWrapper::Node& n, int a) { 566 dist.set(n, a); 567 } 568 int operator[](const typename MapGraphWrapper::Node& n) 571 void set(const typename MapGraphWrapper::Node& n, int a) { 572 dist.set(n, a); 573 } 574 int operator[](const typename MapGraphWrapper::Node& n) 569 575 { return dist[n]; } 570 // int get(const typename MapGraphWrapper::Node& n) const { 576 // int get(const typename MapGraphWrapper::Node& n) const { 571 577 // return dist[n]; } 572 // bool get(const typename MapGraphWrapper::Edge& e) const { 578 // bool get(const typename MapGraphWrapper::Edge& e) const { 573 579 // return (dist.get(g>tail(e))<dist.get(g>head(e))); } 574 bool operator[](const typename MapGraphWrapper::Edge& e) const { 575 return (dist[g>tail(e)]<dist[g>head(e)]); 580 bool operator[](const typename MapGraphWrapper::Edge& e) const { 581 return (dist[g>tail(e)]<dist[g>head(e)]); 576 582 } 577 583 }; 578 584 579 585 }; 580 586 581 587 582 588 template <typename Graph, typename Num, typename CapMap, typename FlowMap> 583 void MaxFlow<Graph, Num, CapMap, FlowMap>::preflowPhase0( flowEnum fe ) 589 void MaxFlow<Graph, Num, CapMap, FlowMap>::preflowPhase0( flowEnum fe ) 584 590 { 585 586 int heur0=(int)(H0*n); //time while running 'bound decrease' 591 592 int heur0=(int)(H0*n); //time while running 'bound decrease' 587 593 int heur1=(int)(H1*n); //time while running 'highest label' 588 594 int heur=heur1; //starting time interval (#of relabels) 589 595 int numrelabel=0; 590 591 bool what_heur=1; 596 597 bool what_heur=1; 592 598 //It is 0 in case 'bound decrease' and 1 in case 'highest label' 593 599 594 bool end=false; 595 //Needed for 'bound decrease', true means no active nodes are above bound b. 600 bool end=false; 601 //Needed for 'bound decrease', true means no active nodes are above bound 602 //b. 596 603 597 604 int k=n2; //bound on the highest level under n containing a node 598 605 int b=k; //bound on the highest level under n of an active node 599 606 600 607 VecStack active(n); 601 608 602 609 NNMap left(*g, INVALID); 603 610 NNMap right(*g, INVALID); … … 608 615 for(g>first(v); g>valid(v); g>next(v)) level.set(v,n); 609 616 //setting each node to level n 610 611 switch ( fe) {612 case PRE FLOW:617 618 switch (fe) { 619 case PRE_FLOW: 613 620 { 614 //co unting the excess621 //computing the excess 615 622 NodeIt v; 616 623 for(g>first(v); g>valid(v); g>next(v)) { 617 624 Num exc=0; 618 625 619 626 InEdgeIt e; 620 627 for(g>first(e,v); g>valid(e); g>next(e)) exc+=(*flow)[e]; 621 628 OutEdgeIt f; 622 629 for(g>first(f,v); g>valid(f); g>next(f)) exc=(*flow)[f]; 623 624 excess.set(v,exc); 625 630 631 excess.set(v,exc); 632 626 633 //putting the active nodes into the stack 627 634 int lev=level[v]; … … 632 639 case GEN_FLOW: 633 640 { 634 // Counting the excess of t641 //computing the excess of t 635 642 Num exc=0; 636 643 637 644 InEdgeIt e; 638 645 for(g>first(e,t); g>valid(e); g>next(e)) exc+=(*flow)[e]; 639 646 OutEdgeIt f; 640 647 for(g>first(f,t); g>valid(f); g>next(f)) exc=(*flow)[f]; 641 642 excess.set(t,exc); 643 648 649 excess.set(t,exc); 650 644 651 break; 645 652 } 646 default: 647 break; 648 } 649 650 preflowPreproc( fe, active, level_list, left, right ); 651 //End of preprocessing 652 653 653 default:; 654 } 655 656 preflowPreproc(fe, active, level_list, left, right); 657 //End of preprocessing 658 659 654 660 //Push/relabel on the highest level active nodes. 655 661 while ( true ) { … … 660 666 } else break; 661 667 } 662 663 if ( active[b].empty() ) b; 668 669 if ( active[b].empty() ) b; 664 670 else { 665 end=false; 671 end=false; 666 672 Node w=active[b].top(); 667 673 active[b].pop(); 668 674 int newlevel=push(w,active); 669 if ( excess[w] > 0 ) relabel(w, newlevel, active, level_list, 675 if ( excess[w] > 0 ) relabel(w, newlevel, active, level_list, 670 676 left, right, b, k, what_heur); 671 672 ++numrelabel; 677 678 ++numrelabel; 673 679 if ( numrelabel >= heur ) { 674 680 numrelabel=0; … … 680 686 what_heur=1; 681 687 heur=heur1; 682 b=k; 683 } 684 } 685 } 686 } 688 b=k; 689 } 690 } 691 } 692 } 687 693 } 688 694 … … 690 696 691 697 template <typename Graph, typename Num, typename CapMap, typename FlowMap> 692 void MaxFlow<Graph, Num, CapMap, FlowMap>::preflowPhase1() 698 void MaxFlow<Graph, Num, CapMap, FlowMap>::preflowPhase1() 693 699 { 694 700 695 701 int k=n2; //bound on the highest level under n containing a node 696 702 int b=k; //bound on the highest level under n of an active node 697 703 698 704 VecStack active(n); 699 705 level.set(s,0); 700 706 std::queue<Node> bfs_queue; 701 707 bfs_queue.push(s); 702 708 703 709 while (!bfs_queue.empty()) { 704 705 Node v=bfs_queue.front(); 710 711 Node v=bfs_queue.front(); 706 712 bfs_queue.pop(); 707 713 int l=level[v]+1; 708 714 709 715 InEdgeIt e; 710 716 for(g>first(e,v); g>valid(e); g>next(e)) { 711 717 if ( (*capacity)[e] <= (*flow)[e] ) continue; 712 718 Node u=g>tail(e); 713 if ( level[u] >= n ) { 719 if ( level[u] >= n ) { 714 720 bfs_queue.push(u); 715 721 level.set(u, l); … … 717 723 } 718 724 } 719 725 720 726 OutEdgeIt f; 721 727 for(g>first(f,v); g>valid(f); g>next(f)) { 722 728 if ( 0 >= (*flow)[f] ) continue; 723 729 Node u=g>head(f); 724 if ( level[u] >= n ) { 730 if ( level[u] >= n ) { 725 731 bfs_queue.push(u); 726 732 level.set(u, l); … … 732 738 733 739 while ( true ) { 734 740 735 741 if ( b == 0 ) break; 736 742 737 if ( active[b].empty() ) b; 743 if ( active[b].empty() ) b; 738 744 else { 739 745 Node w=active[b].top(); 740 746 active[b].pop(); 741 int newlevel=push(w,active); 747 int newlevel=push(w,active); 742 748 743 749 //relabel … … 754 760 755 761 template <typename Graph, typename Num, typename CapMap, typename FlowMap> 756 bool MaxFlow<Graph, Num, CapMap, FlowMap>::augmentOnShortestPath() 762 bool MaxFlow<Graph, Num, CapMap, FlowMap>::augmentOnShortestPath() 757 763 { 758 764 ResGW res_graph(*g, *capacity, *flow); 759 765 bool _augment=false; 760 766 761 767 //ReachedMap level(res_graph); 762 768 FOR_EACH_LOC(typename Graph::NodeIt, e, *g) level.set(e, 0); 763 769 BfsIterator<ResGW, ReachedMap> bfs(res_graph, level); 764 770 bfs.pushAndSetReached(s); 765 766 typename ResGW::template NodeMap<ResGWEdge> pred(res_graph); 771 772 typename ResGW::template NodeMap<ResGWEdge> pred(res_graph); 767 773 pred.set(s, INVALID); 768 774 769 775 typename ResGW::template NodeMap<Num> free(res_graph); 770 776 771 777 //searching for augmenting path 772 while ( !bfs.finished() ) { 778 while ( !bfs.finished() ) { 773 779 ResGWOutEdgeIt e=bfs; 774 780 if (res_graph.valid(e) && bfs.isBNodeNewlyReached()) { … … 779 785 free.set(w, std::min(free[v], res_graph.resCap(e))); 780 786 } else { 781 free.set(w, res_graph.resCap(e)); 787 free.set(w, res_graph.resCap(e)); 782 788 } 783 789 if (res_graph.head(e)==t) { _augment=true; break; } 784 790 } 785 791 786 792 ++bfs; 787 793 } //end of searching augmenting path … … 790 796 Node n=t; 791 797 Num augment_value=free[t]; 792 while (res_graph.valid(pred[n])) { 798 while (res_graph.valid(pred[n])) { 793 799 ResGWEdge e=pred[n]; 794 res_graph.augment(e, augment_value); 800 res_graph.augment(e, augment_value); 795 801 n=res_graph.tail(e); 796 802 } … … 806 812 807 813 808 809 810 814 template <typename Graph, typename Num, typename CapMap, typename FlowMap> 811 template<typename MutableGraph> 812 bool MaxFlow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow() 813 { 815 template<typename MutableGraph> 816 bool MaxFlow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow() 817 { 814 818 typedef MutableGraph MG; 815 819 bool _augment=false; … … 822 826 BfsIterator<ResGW, ReachedMap> bfs(res_graph, level); 823 827 bfs.pushAndSetReached(s); 824 typename ResGW::template NodeMap<int> 828 typename ResGW::template NodeMap<int> 825 829 dist(res_graph); //filled up with 0's 826 830 … … 828 832 //with the set of edges which are on shortest paths 829 833 MG F; 830 typename ResGW::template NodeMap<typename MG::Node> 834 typename ResGW::template NodeMap<typename MG::Node> 831 835 res_graph_to_F(res_graph); 832 836 { … … 842 846 typename MG::template EdgeMap<Num> residual_capacity(F); 843 847 844 while ( !bfs.finished() ) { 848 while ( !bfs.finished() ) { 845 849 ResGWOutEdgeIt e=bfs; 846 850 if (res_graph.valid(e)) { 847 851 if (bfs.isBNodeNewlyReached()) { 848 852 dist.set(res_graph.head(e), dist[res_graph.tail(e)]+1); 849 typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.tail(e)], res_graph_to_F[res_graph.head(e)]); 853 typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.tail(e)], 854 res_graph_to_F[res_graph.head(e)]); 850 855 original_edge.update(); 851 856 original_edge.set(f, e); … … 854 859 } else { 855 860 if (dist[res_graph.head(e)]==(dist[res_graph.tail(e)]+1)) { 856 typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.tail(e)], res_graph_to_F[res_graph.head(e)]); 861 typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.tail(e)], 862 res_graph_to_F[res_graph.head(e)]); 857 863 original_edge.update(); 858 864 original_edge.set(f, e); … … 877 883 typename MG::template NodeMap<Num> free(F); 878 884 879 dfs.pushAndSetReached(sF); 885 dfs.pushAndSetReached(sF); 880 886 while (!dfs.finished()) { 881 887 ++dfs; … … 888 894 free.set(w, std::min(free[v], residual_capacity[dfs])); 889 895 } else { 890 free.set(w, residual_capacity[dfs]); 891 } 892 if (w==tF) { 893 __augment=true; 896 free.set(w, residual_capacity[dfs]); 897 } 898 if (w==tF) { 899 __augment=true; 894 900 _augment=true; 895 break; 896 } 897 901 break; 902 } 903 898 904 } else { 899 905 F.erase(/*typename MG::OutEdgeIt*/(dfs)); 900 906 } 901 } 907 } 902 908 } 903 909 … … 905 911 typename MG::Node n=tF; 906 912 Num augment_value=free[tF]; 907 while (F.valid(pred[n])) { 913 while (F.valid(pred[n])) { 908 914 typename MG::Edge e=pred[n]; 909 res_graph.augment(original_edge[e], augment_value); 915 res_graph.augment(original_edge[e], augment_value); 910 916 n=F.tail(e); 911 if (residual_capacity[e]==augment_value) 912 F.erase(e); 913 else 917 if (residual_capacity[e]==augment_value) 918 F.erase(e); 919 else 914 920 residual_capacity.set(e, residual_capacity[e]augment_value); 915 921 } 916 922 } 917 918 } 919 923 924 } 925 920 926 return _augment; 921 927 } … … 924 930 925 931 926 927 928 932 template <typename Graph, typename Num, typename CapMap, typename FlowMap> 929 bool MaxFlow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow2() 933 bool MaxFlow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow2() 930 934 { 931 935 bool _augment=false; 932 936 933 937 ResGW res_graph(*g, *capacity, *flow); 934 938 935 939 //ReachedMap level(res_graph); 936 940 FOR_EACH_LOC(typename Graph::NodeIt, e, *g) level.set(e, 0); … … 939 943 bfs.pushAndSetReached(s); 940 944 DistanceMap<ResGW> dist(res_graph); 941 while ( !bfs.finished() ) { 945 while ( !bfs.finished() ) { 942 946 ResGWOutEdgeIt e=bfs; 943 947 if (res_graph.valid(e) && bfs.isBNodeNewlyReached()) { … … 949 953 //Subgraph containing the edges on some shortest paths 950 954 ConstMap<typename ResGW::Node, bool> true_map(true); 951 typedef SubGraphWrapper<ResGW, ConstMap<typename ResGW::Node, bool>, 955 typedef SubGraphWrapper<ResGW, ConstMap<typename ResGW::Node, bool>, 952 956 DistanceMap<ResGW> > FilterResGW; 953 957 FilterResGW filter_res_graph(res_graph, true_map, dist); 954 958 955 //Subgraph, which is able to delete edges which are already 959 //Subgraph, which is able to delete edges which are already 956 960 //met by the dfs 957 typename FilterResGW::template NodeMap<typename FilterResGW::OutEdgeIt> 961 typename FilterResGW::template NodeMap<typename FilterResGW::OutEdgeIt> 958 962 first_out_edges(filter_res_graph); 959 963 typename FilterResGW::NodeIt v; 960 for(filter_res_graph.first(v); filter_res_graph.valid(v); 961 filter_res_graph.next(v)) 964 for(filter_res_graph.first(v); filter_res_graph.valid(v); 965 filter_res_graph.next(v)) 962 966 { 963 967 typename FilterResGW::OutEdgeIt e; … … 975 979 __augment=false; 976 980 //computing blocking flow with dfs 977 DfsIterator< ErasingResGW, 978 typename ErasingResGW::template NodeMap<bool> > 981 DfsIterator< ErasingResGW, 982 typename ErasingResGW::template NodeMap<bool> > 979 983 dfs(erasing_res_graph); 980 984 typename ErasingResGW:: 981 template NodeMap<typename ErasingResGW::OutEdgeIt> 982 pred(erasing_res_graph); 985 template NodeMap<typename ErasingResGW::OutEdgeIt> 986 pred(erasing_res_graph); 983 987 pred.set(s, INVALID); 984 988 //invalid iterators for sources 985 989 986 typename ErasingResGW::template NodeMap<Num> 990 typename ErasingResGW::template NodeMap<Num> 987 991 free1(erasing_res_graph); 988 992 989 dfs.pushAndSetReached( 990 typename ErasingResGW::Node( 991 typename FilterResGW::Node( 992 typename ResGW::Node(s) 993 ) 994 ) 995 ); 993 dfs.pushAndSetReached 994 ///\bug hugo 0.2 995 (typename ErasingResGW::Node 996 (typename FilterResGW::Node 997 (typename ResGW::Node(s) 998 ) 999 ) 1000 ); 996 1001 while (!dfs.finished()) { 997 1002 ++dfs; 998 if (erasing_res_graph.valid( 999 typename ErasingResGW::OutEdgeIt(dfs))) 1000 { 1003 if (erasing_res_graph.valid(typename ErasingResGW::OutEdgeIt(dfs))) 1004 { 1001 1005 if (dfs.isBNodeNewlyReached()) { 1002 1006 1003 1007 typename ErasingResGW::Node v=erasing_res_graph.aNode(dfs); 1004 1008 typename ErasingResGW::Node w=erasing_res_graph.bNode(dfs); … … 1006 1010 pred.set(w, /*typename ErasingResGW::OutEdgeIt*/(dfs)); 1007 1011 if (erasing_res_graph.valid(pred[v])) { 1008 free1.set(w, std::min(free1[v], res_graph.resCap( 1009 typename ErasingResGW::OutEdgeIt(dfs)))); 1012 free1.set 1013 (w, std::min(free1[v], res_graph.resCap 1014 (typename ErasingResGW::OutEdgeIt(dfs)))); 1010 1015 } else { 1011 free1.set(w, res_graph.resCap( 1012 typename ErasingResGW::OutEdgeIt(dfs))); 1016 free1.set 1017 (w, res_graph.resCap 1018 (typename ErasingResGW::OutEdgeIt(dfs))); 1013 1019 } 1014 1015 if (w==t) { 1016 __augment=true; 1020 1021 if (w==t) { 1022 __augment=true; 1017 1023 _augment=true; 1018 break; 1024 break; 1019 1025 } 1020 1026 } else { … … 1022 1028 } 1023 1029 } 1024 } 1030 } 1025 1031 1026 1032 if (__augment) { 1027 typename ErasingResGW::Node n=typename FilterResGW::Node(typename ResGW::Node(t)); 1033 typename ErasingResGW::Node 1034 n=typename FilterResGW::Node(typename ResGW::Node(t)); 1028 1035 // typename ResGW::NodeMap<Num> a(res_graph); 1029 1036 // typename ResGW::Node b; … … 1036 1043 // Num j2=a2[b2]; 1037 1044 Num augment_value=free1[n]; 1038 while (erasing_res_graph.valid(pred[n])) { 1045 while (erasing_res_graph.valid(pred[n])) { 1039 1046 typename ErasingResGW::OutEdgeIt e=pred[n]; 1040 1047 res_graph.augment(e, augment_value); … … 1044 1051 } 1045 1052 } 1046 1047 } //while (__augment) 1048 1053 1054 } //while (__augment) 1055 1049 1056 return _augment; 1050 1057 } 1051 1058 1052 1059 1053 1054 1055 1060 } //namespace hugo 1056 1061
Note: See TracChangeset
for help on using the changeset viewer.