gravatar
alpar (Alpar Juttner)
alpar@cs.elte.hu
Merge bugfix #372 to branch 1.1
0 1 0
merge 1.1
0 files changed with 24 insertions and 20 deletions:
↑ Collapse diff ↑
Show white space 384 line context
... ...
@@ -369,541 +369,545 @@
369 369
    /// using this function.
370 370
    const Elevator& elevator() const {
371 371
      return *_level;
372 372
    }
373 373

	
374 374
    /// \brief Sets the tolerance used by algorithm.
375 375
    ///
376 376
    /// Sets the tolerance used by algorithm.
377 377
    Preflow& tolerance(const Tolerance& tolerance) {
378 378
      _tolerance = tolerance;
379 379
      return *this;
380 380
    }
381 381

	
382 382
    /// \brief Returns a const reference to the tolerance.
383 383
    ///
384 384
    /// Returns a const reference to the tolerance.
385 385
    const Tolerance& tolerance() const {
386 386
      return _tolerance;
387 387
    }
388 388

	
389 389
    /// \name Execution Control
390 390
    /// The simplest way to execute the preflow algorithm is to use
391 391
    /// \ref run() or \ref runMinCut().\n
392 392
    /// If you need more control on the initial solution or the execution,
393 393
    /// first you have to call one of the \ref init() functions, then
394 394
    /// \ref startFirstPhase() and if you need it \ref startSecondPhase().
395 395

	
396 396
    ///@{
397 397

	
398 398
    /// \brief Initializes the internal data structures.
399 399
    ///
400 400
    /// Initializes the internal data structures and sets the initial
401 401
    /// flow to zero on each arc.
402 402
    void init() {
403 403
      createStructures();
404 404

	
405 405
      _phase = true;
406 406
      for (NodeIt n(_graph); n != INVALID; ++n) {
407 407
        (*_excess)[n] = 0;
408 408
      }
409 409

	
410 410
      for (ArcIt e(_graph); e != INVALID; ++e) {
411 411
        _flow->set(e, 0);
412 412
      }
413 413

	
414 414
      typename Digraph::template NodeMap<bool> reached(_graph, false);
415 415

	
416 416
      _level->initStart();
417 417
      _level->initAddItem(_target);
418 418

	
419 419
      std::vector<Node> queue;
420 420
      reached[_source] = true;
421 421

	
422 422
      queue.push_back(_target);
423 423
      reached[_target] = true;
424 424
      while (!queue.empty()) {
425 425
        _level->initNewLevel();
426 426
        std::vector<Node> nqueue;
427 427
        for (int i = 0; i < int(queue.size()); ++i) {
428 428
          Node n = queue[i];
429 429
          for (InArcIt e(_graph, n); e != INVALID; ++e) {
430 430
            Node u = _graph.source(e);
431 431
            if (!reached[u] && _tolerance.positive((*_capacity)[e])) {
432 432
              reached[u] = true;
433 433
              _level->initAddItem(u);
434 434
              nqueue.push_back(u);
435 435
            }
436 436
          }
437 437
        }
438 438
        queue.swap(nqueue);
439 439
      }
440 440
      _level->initFinish();
441 441

	
442 442
      for (OutArcIt e(_graph, _source); e != INVALID; ++e) {
443 443
        if (_tolerance.positive((*_capacity)[e])) {
444 444
          Node u = _graph.target(e);
445 445
          if ((*_level)[u] == _level->maxLevel()) continue;
446 446
          _flow->set(e, (*_capacity)[e]);
447 447
          (*_excess)[u] += (*_capacity)[e];
448 448
          if (u != _target && !_level->active(u)) {
449 449
            _level->activate(u);
450 450
          }
451 451
        }
452 452
      }
453 453
    }
454 454

	
455 455
    /// \brief Initializes the internal data structures using the
456 456
    /// given flow map.
457 457
    ///
458 458
    /// Initializes the internal data structures and sets the initial
459 459
    /// flow to the given \c flowMap. The \c flowMap should contain a
460 460
    /// flow or at least a preflow, i.e. at each node excluding the
461 461
    /// source node the incoming flow should greater or equal to the
462 462
    /// outgoing flow.
463 463
    /// \return \c false if the given \c flowMap is not a preflow.
464 464
    template <typename FlowMap>
465 465
    bool init(const FlowMap& flowMap) {
466 466
      createStructures();
467 467

	
468 468
      for (ArcIt e(_graph); e != INVALID; ++e) {
469 469
        _flow->set(e, flowMap[e]);
470 470
      }
471 471

	
472 472
      for (NodeIt n(_graph); n != INVALID; ++n) {
473 473
        Value excess = 0;
474 474
        for (InArcIt e(_graph, n); e != INVALID; ++e) {
475 475
          excess += (*_flow)[e];
476 476
        }
477 477
        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
478 478
          excess -= (*_flow)[e];
479 479
        }
480 480
        if (excess < 0 && n != _source) return false;
481 481
        (*_excess)[n] = excess;
482 482
      }
483 483

	
484 484
      typename Digraph::template NodeMap<bool> reached(_graph, false);
485 485

	
486 486
      _level->initStart();
487 487
      _level->initAddItem(_target);
488 488

	
489 489
      std::vector<Node> queue;
490 490
      reached[_source] = true;
491 491

	
492 492
      queue.push_back(_target);
493 493
      reached[_target] = true;
494 494
      while (!queue.empty()) {
495 495
        _level->initNewLevel();
496 496
        std::vector<Node> nqueue;
497 497
        for (int i = 0; i < int(queue.size()); ++i) {
498 498
          Node n = queue[i];
499 499
          for (InArcIt e(_graph, n); e != INVALID; ++e) {
500 500
            Node u = _graph.source(e);
501 501
            if (!reached[u] &&
502 502
                _tolerance.positive((*_capacity)[e] - (*_flow)[e])) {
503 503
              reached[u] = true;
504 504
              _level->initAddItem(u);
505 505
              nqueue.push_back(u);
506 506
            }
507 507
          }
508 508
          for (OutArcIt e(_graph, n); e != INVALID; ++e) {
509 509
            Node v = _graph.target(e);
510 510
            if (!reached[v] && _tolerance.positive((*_flow)[e])) {
511 511
              reached[v] = true;
512 512
              _level->initAddItem(v);
513 513
              nqueue.push_back(v);
514 514
            }
515 515
          }
516 516
        }
517 517
        queue.swap(nqueue);
518 518
      }
519 519
      _level->initFinish();
520 520

	
521 521
      for (OutArcIt e(_graph, _source); e != INVALID; ++e) {
522 522
        Value rem = (*_capacity)[e] - (*_flow)[e];
523 523
        if (_tolerance.positive(rem)) {
524 524
          Node u = _graph.target(e);
525 525
          if ((*_level)[u] == _level->maxLevel()) continue;
526 526
          _flow->set(e, (*_capacity)[e]);
527 527
          (*_excess)[u] += rem;
528 528
          if (u != _target && !_level->active(u)) {
529 529
            _level->activate(u);
530 530
          }
531 531
        }
532 532
      }
533 533
      for (InArcIt e(_graph, _source); e != INVALID; ++e) {
534 534
        Value rem = (*_flow)[e];
535 535
        if (_tolerance.positive(rem)) {
536 536
          Node v = _graph.source(e);
537 537
          if ((*_level)[v] == _level->maxLevel()) continue;
538 538
          _flow->set(e, 0);
539 539
          (*_excess)[v] += rem;
540 540
          if (v != _target && !_level->active(v)) {
541 541
            _level->activate(v);
542 542
          }
543 543
        }
544 544
      }
545 545
      return true;
546 546
    }
547 547

	
548 548
    /// \brief Starts the first phase of the preflow algorithm.
549 549
    ///
550 550
    /// The preflow algorithm consists of two phases, this method runs
551 551
    /// the first phase. After the first phase the maximum flow value
552 552
    /// and a minimum value cut can already be computed, although a
553 553
    /// maximum flow is not yet obtained. So after calling this method
554 554
    /// \ref flowValue() returns the value of a maximum flow and \ref
555 555
    /// minCut() returns a minimum cut.
556 556
    /// \pre One of the \ref init() functions must be called before
557 557
    /// using this function.
558 558
    void startFirstPhase() {
559 559
      _phase = true;
560 560

	
561
      Node n = _level->highestActive();
562
      int level = _level->highestActiveLevel();
563
      while (n != INVALID) {
561
      while (true) {
564 562
        int num = _node_num;
565 563

	
566
        while (num > 0 && n != INVALID) {
564
        Node n = INVALID;
565
        int level = -1;
566

	
567
        while (num > 0) {
568
          n = _level->highestActive();
569
          if (n == INVALID) goto first_phase_done;
570
          level = _level->highestActiveLevel();
571
          --num;
572
          
567 573
          Value excess = (*_excess)[n];
568 574
          int new_level = _level->maxLevel();
569 575

	
570 576
          for (OutArcIt e(_graph, n); e != INVALID; ++e) {
571 577
            Value rem = (*_capacity)[e] - (*_flow)[e];
572 578
            if (!_tolerance.positive(rem)) continue;
573 579
            Node v = _graph.target(e);
574 580
            if ((*_level)[v] < level) {
575 581
              if (!_level->active(v) && v != _target) {
576 582
                _level->activate(v);
577 583
              }
578 584
              if (!_tolerance.less(rem, excess)) {
579 585
                _flow->set(e, (*_flow)[e] + excess);
580 586
                (*_excess)[v] += excess;
581 587
                excess = 0;
582 588
                goto no_more_push_1;
583 589
              } else {
584 590
                excess -= rem;
585 591
                (*_excess)[v] += rem;
586 592
                _flow->set(e, (*_capacity)[e]);
587 593
              }
588 594
            } else if (new_level > (*_level)[v]) {
589 595
              new_level = (*_level)[v];
590 596
            }
591 597
          }
592 598

	
593 599
          for (InArcIt e(_graph, n); e != INVALID; ++e) {
594 600
            Value rem = (*_flow)[e];
595 601
            if (!_tolerance.positive(rem)) continue;
596 602
            Node v = _graph.source(e);
597 603
            if ((*_level)[v] < level) {
598 604
              if (!_level->active(v) && v != _target) {
599 605
                _level->activate(v);
600 606
              }
601 607
              if (!_tolerance.less(rem, excess)) {
602 608
                _flow->set(e, (*_flow)[e] - excess);
603 609
                (*_excess)[v] += excess;
604 610
                excess = 0;
605 611
                goto no_more_push_1;
606 612
              } else {
607 613
                excess -= rem;
608 614
                (*_excess)[v] += rem;
609 615
                _flow->set(e, 0);
610 616
              }
611 617
            } else if (new_level > (*_level)[v]) {
612 618
              new_level = (*_level)[v];
613 619
            }
614 620
          }
615 621

	
616 622
        no_more_push_1:
617 623

	
618 624
          (*_excess)[n] = excess;
619 625

	
620 626
          if (excess != 0) {
621 627
            if (new_level + 1 < _level->maxLevel()) {
622 628
              _level->liftHighestActive(new_level + 1);
623 629
            } else {
624 630
              _level->liftHighestActiveToTop();
625 631
            }
626 632
            if (_level->emptyLevel(level)) {
627 633
              _level->liftToTop(level);
628 634
            }
629 635
          } else {
630 636
            _level->deactivate(n);
631 637
          }
632

	
633
          n = _level->highestActive();
634
          level = _level->highestActiveLevel();
635
          --num;
636 638
        }
637 639

	
638 640
        num = _node_num * 20;
639
        while (num > 0 && n != INVALID) {
641
        while (num > 0) {
642
          while (level >= 0 && _level->activeFree(level)) {
643
            --level;
644
          }
645
          if (level == -1) {
646
            n = _level->highestActive();
647
            level = _level->highestActiveLevel();
648
            if (n == INVALID) goto first_phase_done;
649
          } else {
650
            n = _level->activeOn(level);
651
          }
652
          --num;
653

	
640 654
          Value excess = (*_excess)[n];
641 655
          int new_level = _level->maxLevel();
642 656

	
643 657
          for (OutArcIt e(_graph, n); e != INVALID; ++e) {
644 658
            Value rem = (*_capacity)[e] - (*_flow)[e];
645 659
            if (!_tolerance.positive(rem)) continue;
646 660
            Node v = _graph.target(e);
647 661
            if ((*_level)[v] < level) {
648 662
              if (!_level->active(v) && v != _target) {
649 663
                _level->activate(v);
650 664
              }
651 665
              if (!_tolerance.less(rem, excess)) {
652 666
                _flow->set(e, (*_flow)[e] + excess);
653 667
                (*_excess)[v] += excess;
654 668
                excess = 0;
655 669
                goto no_more_push_2;
656 670
              } else {
657 671
                excess -= rem;
658 672
                (*_excess)[v] += rem;
659 673
                _flow->set(e, (*_capacity)[e]);
660 674
              }
661 675
            } else if (new_level > (*_level)[v]) {
662 676
              new_level = (*_level)[v];
663 677
            }
664 678
          }
665 679

	
666 680
          for (InArcIt e(_graph, n); e != INVALID; ++e) {
667 681
            Value rem = (*_flow)[e];
668 682
            if (!_tolerance.positive(rem)) continue;
669 683
            Node v = _graph.source(e);
670 684
            if ((*_level)[v] < level) {
671 685
              if (!_level->active(v) && v != _target) {
672 686
                _level->activate(v);
673 687
              }
674 688
              if (!_tolerance.less(rem, excess)) {
675 689
                _flow->set(e, (*_flow)[e] - excess);
676 690
                (*_excess)[v] += excess;
677 691
                excess = 0;
678 692
                goto no_more_push_2;
679 693
              } else {
680 694
                excess -= rem;
681 695
                (*_excess)[v] += rem;
682 696
                _flow->set(e, 0);
683 697
              }
684 698
            } else if (new_level > (*_level)[v]) {
685 699
              new_level = (*_level)[v];
686 700
            }
687 701
          }
688 702

	
689 703
        no_more_push_2:
690 704

	
691 705
          (*_excess)[n] = excess;
692 706

	
693 707
          if (excess != 0) {
694 708
            if (new_level + 1 < _level->maxLevel()) {
695 709
              _level->liftActiveOn(level, new_level + 1);
696 710
            } else {
697 711
              _level->liftActiveToTop(level);
698 712
            }
699 713
            if (_level->emptyLevel(level)) {
700 714
              _level->liftToTop(level);
701 715
            }
702 716
          } else {
703 717
            _level->deactivate(n);
704 718
          }
705

	
706
          while (level >= 0 && _level->activeFree(level)) {
707
            --level;
708 719
          }
709
          if (level == -1) {
710
            n = _level->highestActive();
711
            level = _level->highestActiveLevel();
712
          } else {
713
            n = _level->activeOn(level);
714 720
          }
715
          --num;
716
        }
717
      }
721
    first_phase_done:;
718 722
    }
719 723

	
720 724
    /// \brief Starts the second phase of the preflow algorithm.
721 725
    ///
722 726
    /// The preflow algorithm consists of two phases, this method runs
723 727
    /// the second phase. After calling one of the \ref init() functions
724 728
    /// and \ref startFirstPhase() and then \ref startSecondPhase(),
725 729
    /// \ref flowMap() returns a maximum flow, \ref flowValue() returns the
726 730
    /// value of a maximum flow, \ref minCut() returns a minimum cut
727 731
    /// \pre One of the \ref init() functions and \ref startFirstPhase()
728 732
    /// must be called before using this function.
729 733
    void startSecondPhase() {
730 734
      _phase = false;
731 735

	
732 736
      typename Digraph::template NodeMap<bool> reached(_graph);
733 737
      for (NodeIt n(_graph); n != INVALID; ++n) {
734 738
        reached[n] = (*_level)[n] < _level->maxLevel();
735 739
      }
736 740

	
737 741
      _level->initStart();
738 742
      _level->initAddItem(_source);
739 743

	
740 744
      std::vector<Node> queue;
741 745
      queue.push_back(_source);
742 746
      reached[_source] = true;
743 747

	
744 748
      while (!queue.empty()) {
745 749
        _level->initNewLevel();
746 750
        std::vector<Node> nqueue;
747 751
        for (int i = 0; i < int(queue.size()); ++i) {
748 752
          Node n = queue[i];
749 753
          for (OutArcIt e(_graph, n); e != INVALID; ++e) {
750 754
            Node v = _graph.target(e);
751 755
            if (!reached[v] && _tolerance.positive((*_flow)[e])) {
752 756
              reached[v] = true;
753 757
              _level->initAddItem(v);
754 758
              nqueue.push_back(v);
755 759
            }
756 760
          }
757 761
          for (InArcIt e(_graph, n); e != INVALID; ++e) {
758 762
            Node u = _graph.source(e);
759 763
            if (!reached[u] &&
760 764
                _tolerance.positive((*_capacity)[e] - (*_flow)[e])) {
761 765
              reached[u] = true;
762 766
              _level->initAddItem(u);
763 767
              nqueue.push_back(u);
764 768
            }
765 769
          }
766 770
        }
767 771
        queue.swap(nqueue);
768 772
      }
769 773
      _level->initFinish();
770 774

	
771 775
      for (NodeIt n(_graph); n != INVALID; ++n) {
772 776
        if (!reached[n]) {
773 777
          _level->dirtyTopButOne(n);
774 778
        } else if ((*_excess)[n] > 0 && _target != n) {
775 779
          _level->activate(n);
776 780
        }
777 781
      }
778 782

	
779 783
      Node n;
780 784
      while ((n = _level->highestActive()) != INVALID) {
781 785
        Value excess = (*_excess)[n];
782 786
        int level = _level->highestActiveLevel();
783 787
        int new_level = _level->maxLevel();
784 788

	
785 789
        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
786 790
          Value rem = (*_capacity)[e] - (*_flow)[e];
787 791
          if (!_tolerance.positive(rem)) continue;
788 792
          Node v = _graph.target(e);
789 793
          if ((*_level)[v] < level) {
790 794
            if (!_level->active(v) && v != _source) {
791 795
              _level->activate(v);
792 796
            }
793 797
            if (!_tolerance.less(rem, excess)) {
794 798
              _flow->set(e, (*_flow)[e] + excess);
795 799
              (*_excess)[v] += excess;
796 800
              excess = 0;
797 801
              goto no_more_push;
798 802
            } else {
799 803
              excess -= rem;
800 804
              (*_excess)[v] += rem;
801 805
              _flow->set(e, (*_capacity)[e]);
802 806
            }
803 807
          } else if (new_level > (*_level)[v]) {
804 808
            new_level = (*_level)[v];
805 809
          }
806 810
        }
807 811

	
808 812
        for (InArcIt e(_graph, n); e != INVALID; ++e) {
809 813
          Value rem = (*_flow)[e];
810 814
          if (!_tolerance.positive(rem)) continue;
811 815
          Node v = _graph.source(e);
812 816
          if ((*_level)[v] < level) {
813 817
            if (!_level->active(v) && v != _source) {
814 818
              _level->activate(v);
815 819
            }
816 820
            if (!_tolerance.less(rem, excess)) {
817 821
              _flow->set(e, (*_flow)[e] - excess);
818 822
              (*_excess)[v] += excess;
819 823
              excess = 0;
820 824
              goto no_more_push;
821 825
            } else {
822 826
              excess -= rem;
823 827
              (*_excess)[v] += rem;
824 828
              _flow->set(e, 0);
825 829
            }
826 830
          } else if (new_level > (*_level)[v]) {
827 831
            new_level = (*_level)[v];
828 832
          }
829 833
        }
830 834

	
831 835
      no_more_push:
832 836

	
833 837
        (*_excess)[n] = excess;
834 838

	
835 839
        if (excess != 0) {
836 840
          if (new_level + 1 < _level->maxLevel()) {
837 841
            _level->liftHighestActive(new_level + 1);
838 842
          } else {
839 843
            // Calculation error
840 844
            _level->liftHighestActiveToTop();
841 845
          }
842 846
          if (_level->emptyLevel(level)) {
843 847
            // Calculation error
844 848
            _level->liftToTop(level);
845 849
          }
846 850
        } else {
847 851
          _level->deactivate(n);
848 852
        }
849 853

	
850 854
      }
851 855
    }
852 856

	
853 857
    /// \brief Runs the preflow algorithm.
854 858
    ///
855 859
    /// Runs the preflow algorithm.
856 860
    /// \note pf.run() is just a shortcut of the following code.
857 861
    /// \code
858 862
    ///   pf.init();
859 863
    ///   pf.startFirstPhase();
860 864
    ///   pf.startSecondPhase();
861 865
    /// \endcode
862 866
    void run() {
863 867
      init();
864 868
      startFirstPhase();
865 869
      startSecondPhase();
866 870
    }
867 871

	
868 872
    /// \brief Runs the preflow algorithm to compute the minimum cut.
869 873
    ///
870 874
    /// Runs the preflow algorithm to compute the minimum cut.
871 875
    /// \note pf.runMinCut() is just a shortcut of the following code.
872 876
    /// \code
873 877
    ///   pf.init();
874 878
    ///   pf.startFirstPhase();
875 879
    /// \endcode
876 880
    void runMinCut() {
877 881
      init();
878 882
      startFirstPhase();
879 883
    }
880 884

	
881 885
    /// @}
882 886

	
883 887
    /// \name Query Functions
884 888
    /// The results of the preflow algorithm can be obtained using these
885 889
    /// functions.\n
886 890
    /// Either one of the \ref run() "run*()" functions or one of the
887 891
    /// \ref startFirstPhase() "start*()" functions should be called
888 892
    /// before using them.
889 893

	
890 894
    ///@{
891 895

	
892 896
    /// \brief Returns the value of the maximum flow.
893 897
    ///
894 898
    /// Returns the value of the maximum flow by returning the excess
895 899
    /// of the target node. This value equals to the value of
896 900
    /// the maximum flow already after the first phase of the algorithm.
897 901
    ///
898 902
    /// \pre Either \ref run() or \ref init() must be called before
899 903
    /// using this function.
900 904
    Value flowValue() const {
901 905
      return (*_excess)[_target];
902 906
    }
903 907

	
904 908
    /// \brief Returns the flow value on the given arc.
905 909
    ///
906 910
    /// Returns the flow value on the given arc. This method can
907 911
    /// be called after the second phase of the algorithm.
908 912
    ///
909 913
    /// \pre Either \ref run() or \ref init() must be called before
0 comments (0 inline)