108 /// |
108 /// |
109 /// \warning Both value types must be signed and all input data must |
109 /// \warning Both value types must be signed and all input data must |
110 /// be integer. |
110 /// be integer. |
111 /// \warning This algorithm does not support negative costs for such |
111 /// \warning This algorithm does not support negative costs for such |
112 /// arcs that have infinite upper bound. |
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 #ifdef DOXYGEN |
117 #ifdef DOXYGEN |
114 template <typename GR, typename V, typename C, typename TR> |
118 template <typename GR, typename V, typename C, typename TR> |
115 #else |
119 #else |
116 template < typename GR, typename V = int, typename C = V, |
120 template < typename GR, typename V = int, typename C = V, |
117 typename TR = CostScalingDefaultTraits<GR, V, C> > |
121 typename TR = CostScalingDefaultTraits<GR, V, C> > |
157 /// over the feasible flows, but this algroithm cannot handle |
161 /// over the feasible flows, but this algroithm cannot handle |
158 /// these cases. |
162 /// these cases. |
159 UNBOUNDED |
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 private: |
193 private: |
163 |
194 |
164 TEMPLATE_DIGRAPH_TYPEDEFS(GR); |
195 TEMPLATE_DIGRAPH_TYPEDEFS(GR); |
165 |
196 |
166 typedef std::vector<int> IntVector; |
197 typedef std::vector<int> IntVector; |
503 /// |
534 /// |
504 /// This function can be called more than once. All the parameters |
535 /// This function can be called more than once. All the parameters |
505 /// that have been given are kept for the next call, unless |
536 /// that have been given are kept for the next call, unless |
506 /// \ref reset() is called, thus only the modified parameters |
537 /// \ref reset() is called, thus only the modified parameters |
507 /// have to be set again. See \ref reset() for examples. |
538 /// have to be set again. See \ref reset() for examples. |
508 /// However the underlying digraph must not be modified after this |
539 /// However, the underlying digraph must not be modified after this |
509 /// class have been constructed, since it copies the digraph. |
540 /// class have been constructed, since it copies and extends the graph. |
510 /// |
541 /// |
511 /// \param partial_augment By default the algorithm performs |
542 /// \param method The internal method that will be used in the |
512 /// partial augment and relabel operations in the cost scaling |
543 /// algorithm. For more information, see \ref Method. |
513 /// phases. Set this parameter to \c false for using local push and |
544 /// \param factor The cost scaling factor. It must be larger than one. |
514 /// relabel operations instead. |
|
515 /// |
545 /// |
516 /// \return \c INFEASIBLE if no feasible flow exists, |
546 /// \return \c INFEASIBLE if no feasible flow exists, |
517 /// \n \c OPTIMAL if the problem has optimal solution |
547 /// \n \c OPTIMAL if the problem has optimal solution |
518 /// (i.e. it is feasible and bounded), and the algorithm has found |
548 /// (i.e. it is feasible and bounded), and the algorithm has found |
519 /// optimal flow and node potentials (primal and dual solutions), |
549 /// optimal flow and node potentials (primal and dual solutions), |
521 /// and infinite upper bound. It means that the objective function |
551 /// and infinite upper bound. It means that the objective function |
522 /// is unbounded on that arc, however note that it could actually be |
552 /// is unbounded on that arc, however note that it could actually be |
523 /// bounded over the feasible flows, but this algroithm cannot handle |
553 /// bounded over the feasible flows, but this algroithm cannot handle |
524 /// these cases. |
554 /// these cases. |
525 /// |
555 /// |
526 /// \see ProblemType |
556 /// \see ProblemType, Method |
527 ProblemType run(bool partial_augment = true) { |
557 ProblemType run(Method method = PARTIAL_AUGMENT, int factor = 8) { |
|
558 _alpha = factor; |
528 ProblemType pt = init(); |
559 ProblemType pt = init(); |
529 if (pt != OPTIMAL) return pt; |
560 if (pt != OPTIMAL) return pt; |
530 start(partial_augment); |
561 start(method); |
531 return OPTIMAL; |
562 return OPTIMAL; |
532 } |
563 } |
533 |
564 |
534 /// \brief Reset all the parameters that have been given before. |
565 /// \brief Reset all the parameters that have been given before. |
535 /// |
566 /// |
678 private: |
709 private: |
679 |
710 |
680 // Initialize the algorithm |
711 // Initialize the algorithm |
681 ProblemType init() { |
712 ProblemType init() { |
682 if (_res_node_num == 0) return INFEASIBLE; |
713 if (_res_node_num == 0) return INFEASIBLE; |
683 |
|
684 // Scaling factor |
|
685 _alpha = 8; |
|
686 |
714 |
687 // Check the sum of supply values |
715 // Check the sum of supply values |
688 _sum_supply = 0; |
716 _sum_supply = 0; |
689 for (int i = 0; i != _root; ++i) { |
717 for (int i = 0; i != _root; ++i) { |
690 _sum_supply += _supply[i]; |
718 _sum_supply += _supply[i]; |
815 |
843 |
816 return OPTIMAL; |
844 return OPTIMAL; |
817 } |
845 } |
818 |
846 |
819 // Execute the algorithm and transform the results |
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 // Execute the algorithm |
852 // Execute the algorithm |
822 if (partial_augment) { |
853 switch (method) { |
823 startPartialAugment(); |
854 case PUSH: |
824 } else { |
855 startPush(); |
825 startPushRelabel(); |
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 // Compute node potentials for the original costs |
865 // Compute node potentials for the original costs |
829 _arc_vec.clear(); |
866 _arc_vec.clear(); |
830 _cost_vec.clear(); |
867 _cost_vec.clear(); |
849 if (!_forward[j]) _res_cap[j] += _lower[j]; |
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 |
891 /// Execute the algorithm performing augment and relabel operations |
855 /// relabel operations |
892 void startAugment(int max_length = std::numeric_limits<int>::max()) { |
856 void startPartialAugment() { |
|
857 // Paramters for heuristics |
893 // Paramters for heuristics |
858 const int BF_HEURISTIC_EPSILON_BOUND = 1000; |
894 const int BF_HEURISTIC_EPSILON_BOUND = 1000; |
859 const int BF_HEURISTIC_BOUND_FACTOR = 3; |
895 const int BF_HEURISTIC_BOUND_FACTOR = 3; |
860 // Maximum augment path length |
|
861 const int MAX_PATH_LENGTH = 4; |
|
862 |
896 |
863 // Perform cost scaling phases |
897 // Perform cost scaling phases |
864 IntVector pred_arc(_res_node_num); |
898 IntVector pred_arc(_res_node_num); |
865 std::vector<int> path_nodes; |
899 std::vector<int> path_nodes; |
866 for ( ; _epsilon >= 1; _epsilon = _epsilon < _alpha && _epsilon > 1 ? |
900 for ( ; _epsilon >= 1; _epsilon = _epsilon < _alpha && _epsilon > 1 ? |
923 path_nodes.push_back(start); |
957 path_nodes.push_back(start); |
924 |
958 |
925 // Find an augmenting path from the start node |
959 // Find an augmenting path from the start node |
926 int tip = start; |
960 int tip = start; |
927 while (_excess[tip] >= 0 && |
961 while (_excess[tip] >= 0 && |
928 int(path_nodes.size()) <= MAX_PATH_LENGTH) { |
962 int(path_nodes.size()) <= max_length) { |
929 int u; |
963 int u; |
930 LargeCost min_red_cost, rc; |
964 LargeCost min_red_cost, rc; |
931 int last_out = _sum_supply < 0 ? |
965 int last_out = _sum_supply < 0 ? |
932 _first_out[tip+1] : _first_out[tip+1] - 1; |
966 _first_out[tip+1] : _first_out[tip+1] - 1; |
933 for (int a = _next_out[tip]; a != last_out; ++a) { |
967 for (int a = _next_out[tip]; a != last_out; ++a) { |
982 } |
1016 } |
983 } |
1017 } |
984 } |
1018 } |
985 |
1019 |
986 /// Execute the algorithm performing push and relabel operations |
1020 /// Execute the algorithm performing push and relabel operations |
987 void startPushRelabel() { |
1021 void startPush() { |
988 // Paramters for heuristics |
1022 // Paramters for heuristics |
989 const int BF_HEURISTIC_EPSILON_BOUND = 1000; |
1023 const int BF_HEURISTIC_EPSILON_BOUND = 1000; |
990 const int BF_HEURISTIC_BOUND_FACTOR = 3; |
1024 const int BF_HEURISTIC_BOUND_FACTOR = 3; |
991 |
1025 |
992 // Perform cost scaling phases |
1026 // Perform cost scaling phases |