gravatar
kpeter (Peter Kovacs)
kpeter@inf.elte.hu
More options for run() in scaling MCF algorithms (#180) - Three methods can be selected and the scaling factor can be given for CostScaling. - The scaling factor can be given for CapacityScaling.
0 2 0
default
2 files changed with 71 insertions and 42 deletions:
↑ Collapse diff ↑
Ignore white space 12 line context
... ...
@@ -170,13 +170,13 @@
170 170
    CostVector _pi;
171 171
    ValueVector _excess;
172 172
    IntVector _excess_nodes;
173 173
    IntVector _deficit_nodes;
174 174

	
175 175
    Value _delta;
176
    int _phase_num;
176
    int _factor;
177 177
    IntVector _pred;
178 178

	
179 179
  public:
180 180
  
181 181
    /// \brief Constant for infinite upper bounds (capacities).
182 182
    ///
... ...
@@ -510,32 +510,32 @@
510 510
    ///
511 511
    /// This function can be called more than once. All the parameters
512 512
    /// that have been given are kept for the next call, unless
513 513
    /// \ref reset() is called, thus only the modified parameters
514 514
    /// have to be set again. See \ref reset() for examples.
515 515
    /// However the underlying digraph must not be modified after this
516
    /// class have been constructed, since it copies the digraph.
516
    /// class have been constructed, since it copies and extends the graph.
517 517
    ///
518
    /// \param scaling Enable or disable capacity scaling.
519
    /// If the maximum upper bound and/or the amount of total supply
520
    /// is rather small, the algorithm could be slightly faster without
521
    /// scaling.
518
    /// \param factor The capacity scaling factor. It must be larger than
519
    /// one to use scaling. If it is less or equal to one, then scaling
520
    /// will be disabled.
522 521
    ///
523 522
    /// \return \c INFEASIBLE if no feasible flow exists,
524 523
    /// \n \c OPTIMAL if the problem has optimal solution
525 524
    /// (i.e. it is feasible and bounded), and the algorithm has found
526 525
    /// optimal flow and node potentials (primal and dual solutions),
527 526
    /// \n \c UNBOUNDED if the digraph contains an arc of negative cost
528 527
    /// and infinite upper bound. It means that the objective function
529 528
    /// is unbounded on that arc, however note that it could actually be
530 529
    /// bounded over the feasible flows, but this algroithm cannot handle
531 530
    /// these cases.
532 531
    ///
533 532
    /// \see ProblemType
534
    ProblemType run(bool scaling = true) {
535
      ProblemType pt = init(scaling);
533
    ProblemType run(int factor = 4) {
534
      _factor = factor;
535
      ProblemType pt = init();
536 536
      if (pt != OPTIMAL) return pt;
537 537
      return start();
538 538
    }
539 539

	
540 540
    /// \brief Reset all the parameters that have been given before.
541 541
    ///
... ...
@@ -543,13 +543,13 @@
543 543
    /// before using functions \ref lowerMap(), \ref upperMap(),
544 544
    /// \ref costMap(), \ref supplyMap(), \ref stSupply().
545 545
    ///
546 546
    /// It is useful for multiple run() calls. If this function is not
547 547
    /// used, all the parameters given before are kept for the next
548 548
    /// \ref run() call.
549
    /// However the underlying digraph must not be modified after this
549
    /// However, the underlying digraph must not be modified after this
550 550
    /// class have been constructed, since it copies and extends the graph.
551 551
    ///
552 552
    /// For example,
553 553
    /// \code
554 554
    ///   CapacityScaling<ListDigraph> cs(graph);
555 555
    ///
... ...
@@ -674,13 +674,13 @@
674 674

	
675 675
    /// @}
676 676

	
677 677
  private:
678 678

	
679 679
    // Initialize the algorithm
680
    ProblemType init(bool scaling) {
680
    ProblemType init() {
681 681
      if (_node_num == 0) return INFEASIBLE;
682 682

	
683 683
      // Check the sum of supply values
684 684
      _sum_supply = 0;
685 685
      for (int i = 0; i != _root; ++i) {
686 686
        _sum_supply += _supply[i];
... ...
@@ -755,27 +755,25 @@
755 755
          _cost[a] = 0;
756 756
          _cost[_reverse[a]] = 0;
757 757
        }
758 758
      }
759 759

	
760 760
      // Initialize delta value
761
      if (scaling) {
761
      if (_factor > 1) {
762 762
        // With scaling
763 763
        Value max_sup = 0, max_dem = 0;
764 764
        for (int i = 0; i != _node_num; ++i) {
765 765
          if ( _excess[i] > max_sup) max_sup =  _excess[i];
766 766
          if (-_excess[i] > max_dem) max_dem = -_excess[i];
767 767
        }
768 768
        Value max_cap = 0;
769 769
        for (int j = 0; j != _res_arc_num; ++j) {
770 770
          if (_res_cap[j] > max_cap) max_cap = _res_cap[j];
771 771
        }
772 772
        max_sup = std::min(std::min(max_sup, max_dem), max_cap);
773
        _phase_num = 0;
774
        for (_delta = 1; 2 * _delta <= max_sup; _delta *= 2)
775
          ++_phase_num;
773
        for (_delta = 1; 2 * _delta <= max_sup; _delta *= 2) ;
776 774
      } else {
777 775
        // Without scaling
778 776
        _delta = 1;
779 777
      }
780 778

	
781 779
      return OPTIMAL;
... ...
@@ -808,14 +806,12 @@
808 806
    }
809 807

	
810 808
    // Execute the capacity scaling algorithm
811 809
    ProblemType startWithScaling() {
812 810
      // Perform capacity scaling phases
813 811
      int s, t;
814
      int phase_cnt = 0;
815
      int factor = 4;
816 812
      ResidualDijkstra _dijkstra(*this);
817 813
      while (true) {
818 814
        // Saturate all arcs not satisfying the optimality condition
819 815
        for (int u = 0; u != _node_num; ++u) {
820 816
          for (int a = _first_out[u]; a != _first_out[u+1]; ++a) {
821 817
            int v = _target[a];
... ...
@@ -884,14 +880,13 @@
884 880
          _excess[t] += d;
885 881

	
886 882
          if (_excess[s] < _delta) ++next_node;
887 883
        }
888 884

	
889 885
        if (_delta == 1) break;
890
        if (++phase_cnt == _phase_num / 4) factor = 2;
891
        _delta = _delta <= factor ? 1 : _delta / factor;
886
        _delta = _delta <= _factor ? 1 : _delta / _factor;
892 887
      }
893 888

	
894 889
      return OPTIMAL;
895 890
    }
896 891

	
897 892
    // Execute the successive shortest path algorithm
Ignore white space 6 line context
... ...
@@ -107,12 +107,16 @@
107 107
  /// algorithm. By default it is the same as \c V.
108 108
  ///
109 109
  /// \warning Both value types must be signed and all input data must
110 110
  /// be integer.
111 111
  /// \warning This algorithm does not support negative costs for such
112 112
  /// arcs that have infinite upper bound.
113
  ///
114
  /// \note %CostScaling provides three different internal methods,
115
  /// from which the most efficient one is used by default.
116
  /// For more information, see \ref Method.
113 117
#ifdef DOXYGEN
114 118
  template <typename GR, typename V, typename C, typename TR>
115 119
#else
116 120
  template < typename GR, typename V = int, typename C = V,
117 121
             typename TR = CostScalingDefaultTraits<GR, V, C> >
118 122
#endif
... ...
@@ -156,12 +160,39 @@
156 160
      /// on that arc, however note that it could actually be bounded
157 161
      /// over the feasible flows, but this algroithm cannot handle
158 162
      /// these cases.
159 163
      UNBOUNDED
160 164
    };
161 165

	
166
    /// \brief Constants for selecting the internal method.
167
    ///
168
    /// Enum type containing constants for selecting the internal method
169
    /// for the \ref run() function.
170
    ///
171
    /// \ref CostScaling provides three internal methods that differ mainly
172
    /// in their base operations, which are used in conjunction with the
173
    /// relabel operation.
174
    /// By default, the so called \ref PARTIAL_AUGMENT
175
    /// "Partial Augment-Relabel" method is used, which proved to be
176
    /// the most efficient and the most robust on various test inputs.
177
    /// However, the other methods can be selected using the \ref run()
178
    /// function with the proper parameter.
179
    enum Method {
180
      /// Local push operations are used, i.e. flow is moved only on one
181
      /// admissible arc at once.
182
      PUSH,
183
      /// Augment operations are used, i.e. flow is moved on admissible
184
      /// paths from a node with excess to a node with deficit.
185
      AUGMENT,
186
      /// Partial augment operations are used, i.e. flow is moved on 
187
      /// admissible paths started from a node with excess, but the
188
      /// lengths of these paths are limited. This method can be viewed
189
      /// as a combined version of the previous two operations.
190
      PARTIAL_AUGMENT
191
    };
192

	
162 193
  private:
163 194

	
164 195
    TEMPLATE_DIGRAPH_TYPEDEFS(GR);
165 196

	
166 197
    typedef std::vector<int> IntVector;
167 198
    typedef std::vector<char> BoolVector;
... ...
@@ -502,35 +533,35 @@
502 533
    /// \endcode
503 534
    ///
504 535
    /// This function can be called more than once. All the parameters
505 536
    /// that have been given are kept for the next call, unless
506 537
    /// \ref reset() is called, thus only the modified parameters
507 538
    /// have to be set again. See \ref reset() for examples.
508
    /// However the underlying digraph must not be modified after this
509
    /// class have been constructed, since it copies the digraph.
539
    /// However, the underlying digraph must not be modified after this
540
    /// class have been constructed, since it copies and extends the graph.
510 541
    ///
511
    /// \param partial_augment By default the algorithm performs
512
    /// partial augment and relabel operations in the cost scaling
513
    /// phases. Set this parameter to \c false for using local push and
514
    /// relabel operations instead.
542
    /// \param method The internal method that will be used in the
543
    /// algorithm. For more information, see \ref Method.
544
    /// \param factor The cost scaling factor. It must be larger than one.
515 545
    ///
516 546
    /// \return \c INFEASIBLE if no feasible flow exists,
517 547
    /// \n \c OPTIMAL if the problem has optimal solution
518 548
    /// (i.e. it is feasible and bounded), and the algorithm has found
519 549
    /// optimal flow and node potentials (primal and dual solutions),
520 550
    /// \n \c UNBOUNDED if the digraph contains an arc of negative cost
521 551
    /// and infinite upper bound. It means that the objective function
522 552
    /// is unbounded on that arc, however note that it could actually be
523 553
    /// bounded over the feasible flows, but this algroithm cannot handle
524 554
    /// these cases.
525 555
    ///
526
    /// \see ProblemType
527
    ProblemType run(bool partial_augment = true) {
556
    /// \see ProblemType, Method
557
    ProblemType run(Method method = PARTIAL_AUGMENT, int factor = 8) {
558
      _alpha = factor;
528 559
      ProblemType pt = init();
529 560
      if (pt != OPTIMAL) return pt;
530
      start(partial_augment);
561
      start(method);
531 562
      return OPTIMAL;
532 563
    }
533 564

	
534 565
    /// \brief Reset all the parameters that have been given before.
535 566
    ///
536 567
    /// This function resets all the paramaters that have been given
... ...
@@ -678,15 +709,12 @@
678 709
  private:
679 710

	
680 711
    // Initialize the algorithm
681 712
    ProblemType init() {
682 713
      if (_res_node_num == 0) return INFEASIBLE;
683 714

	
684
      // Scaling factor
685
      _alpha = 8;
686

	
687 715
      // Check the sum of supply values
688 716
      _sum_supply = 0;
689 717
      for (int i = 0; i != _root; ++i) {
690 718
        _sum_supply += _supply[i];
691 719
      }
692 720
      if (_sum_supply > 0) return INFEASIBLE;
... ...
@@ -814,18 +842,27 @@
814 842
      }
815 843
      
816 844
      return OPTIMAL;
817 845
    }
818 846

	
819 847
    // Execute the algorithm and transform the results
820
    void start(bool partial_augment) {
848
    void start(Method method) {
849
      // Maximum path length for partial augment
850
      const int MAX_PATH_LENGTH = 4;
851
      
821 852
      // Execute the algorithm
822
      if (partial_augment) {
823
        startPartialAugment();
824
      } else {
825
        startPushRelabel();
853
      switch (method) {
854
        case PUSH:
855
          startPush();
856
          break;
857
        case AUGMENT:
858
          startAugment();
859
          break;
860
        case PARTIAL_AUGMENT:
861
          startAugment(MAX_PATH_LENGTH);
862
          break;
826 863
      }
827 864

	
828 865
      // Compute node potentials for the original costs
829 866
      _arc_vec.clear();
830 867
      _cost_vec.clear();
831 868
      for (int j = 0; j != _res_arc_num; ++j) {
... ...
@@ -848,20 +885,17 @@
848 885
        for (int j = 0; j != limit; ++j) {
849 886
          if (!_forward[j]) _res_cap[j] += _lower[j];
850 887
        }
851 888
      }
852 889
    }
853 890

	
854
    /// Execute the algorithm performing partial augmentation and
855
    /// relabel operations
856
    void startPartialAugment() {
891
    /// Execute the algorithm performing augment and relabel operations
892
    void startAugment(int max_length = std::numeric_limits<int>::max()) {
857 893
      // Paramters for heuristics
858 894
      const int BF_HEURISTIC_EPSILON_BOUND = 1000;
859 895
      const int BF_HEURISTIC_BOUND_FACTOR  = 3;
860
      // Maximum augment path length
861
      const int MAX_PATH_LENGTH = 4;
862 896

	
863 897
      // Perform cost scaling phases
864 898
      IntVector pred_arc(_res_node_num);
865 899
      std::vector<int> path_nodes;
866 900
      for ( ; _epsilon >= 1; _epsilon = _epsilon < _alpha && _epsilon > 1 ?
867 901
                                        1 : _epsilon / _alpha )
... ...
@@ -922,13 +956,13 @@
922 956
          path_nodes.clear();
923 957
          path_nodes.push_back(start);
924 958

	
925 959
          // Find an augmenting path from the start node
926 960
          int tip = start;
927 961
          while (_excess[tip] >= 0 &&
928
                 int(path_nodes.size()) <= MAX_PATH_LENGTH) {
962
                 int(path_nodes.size()) <= max_length) {
929 963
            int u;
930 964
            LargeCost min_red_cost, rc;
931 965
            int last_out = _sum_supply < 0 ?
932 966
              _first_out[tip+1] : _first_out[tip+1] - 1;
933 967
            for (int a = _next_out[tip]; a != last_out; ++a) {
934 968
              if (_res_cap[a] > 0 &&
... ...
@@ -981,13 +1015,13 @@
981 1015
          }
982 1016
        }
983 1017
      }
984 1018
    }
985 1019

	
986 1020
    /// Execute the algorithm performing push and relabel operations
987
    void startPushRelabel() {
1021
    void startPush() {
988 1022
      // Paramters for heuristics
989 1023
      const int BF_HEURISTIC_EPSILON_BOUND = 1000;
990 1024
      const int BF_HEURISTIC_BOUND_FACTOR  = 3;
991 1025

	
992 1026
      // Perform cost scaling phases
993 1027
      BoolVector hyper(_res_node_num, false);
0 comments (0 inline)