33 #include <lemon/static_graph.h> |
33 #include <lemon/static_graph.h> |
34 #include <lemon/adaptors.h> |
34 #include <lemon/adaptors.h> |
35 #include <lemon/circulation.h> |
35 #include <lemon/circulation.h> |
36 #include <lemon/bellman_ford.h> |
36 #include <lemon/bellman_ford.h> |
37 #include <lemon/howard_mmc.h> |
37 #include <lemon/howard_mmc.h> |
|
38 #include <lemon/hartmann_orlin_mmc.h> |
38 |
39 |
39 namespace lemon { |
40 namespace lemon { |
40 |
41 |
41 /// \addtogroup min_cost_flow_algs |
42 /// \addtogroup min_cost_flow_algs |
42 /// @{ |
43 /// @{ |
925 } |
926 } |
926 } |
927 } |
927 |
928 |
928 // Execute the "Minimum Mean Cycle Canceling" method |
929 // Execute the "Minimum Mean Cycle Canceling" method |
929 void startMinMeanCycleCanceling() { |
930 void startMinMeanCycleCanceling() { |
930 typedef SimplePath<StaticDigraph> SPath; |
931 typedef Path<StaticDigraph> SPath; |
931 typedef typename SPath::ArcIt SPathArcIt; |
932 typedef typename SPath::ArcIt SPathArcIt; |
932 typedef typename HowardMmc<StaticDigraph, CostArcMap> |
933 typedef typename HowardMmc<StaticDigraph, CostArcMap> |
933 ::template SetPath<SPath>::Create MMC; |
934 ::template SetPath<SPath>::Create HwMmc; |
|
935 typedef typename HartmannOrlinMmc<StaticDigraph, CostArcMap> |
|
936 ::template SetPath<SPath>::Create HoMmc; |
|
937 |
|
938 const double HW_ITER_LIMIT_FACTOR = 1.0; |
|
939 const int HW_ITER_LIMIT_MIN_VALUE = 5; |
|
940 |
|
941 const int hw_iter_limit = |
|
942 std::max(static_cast<int>(HW_ITER_LIMIT_FACTOR * _node_num), |
|
943 HW_ITER_LIMIT_MIN_VALUE); |
934 |
944 |
935 SPath cycle; |
945 SPath cycle; |
936 MMC mmc(_sgr, _cost_map); |
946 HwMmc hw_mmc(_sgr, _cost_map); |
937 mmc.cycle(cycle); |
947 hw_mmc.cycle(cycle); |
938 buildResidualNetwork(); |
948 buildResidualNetwork(); |
939 while (mmc.findCycleMean() && mmc.cycleCost() < 0) { |
949 while (true) { |
940 // Find the cycle |
950 |
941 mmc.findCycle(); |
951 typename HwMmc::TerminationCause hw_tc = |
942 |
952 hw_mmc.findCycleMean(hw_iter_limit); |
|
953 if (hw_tc == HwMmc::ITERATION_LIMIT) { |
|
954 // Howard's algorithm reached the iteration limit, start a |
|
955 // strongly polynomial algorithm instead |
|
956 HoMmc ho_mmc(_sgr, _cost_map); |
|
957 ho_mmc.cycle(cycle); |
|
958 // Find a minimum mean cycle (Hartmann-Orlin algorithm) |
|
959 if (!(ho_mmc.findCycleMean() && ho_mmc.cycleCost() < 0)) break; |
|
960 ho_mmc.findCycle(); |
|
961 } else { |
|
962 // Find a minimum mean cycle (Howard algorithm) |
|
963 if (!(hw_tc == HwMmc::OPTIMAL && hw_mmc.cycleCost() < 0)) break; |
|
964 hw_mmc.findCycle(); |
|
965 } |
|
966 |
943 // Compute delta value |
967 // Compute delta value |
944 Value delta = INF; |
968 Value delta = INF; |
945 for (SPathArcIt a(cycle); a != INVALID; ++a) { |
969 for (SPathArcIt a(cycle); a != INVALID; ++a) { |
946 Value d = _res_cap[_id_vec[_sgr.id(a)]]; |
970 Value d = _res_cap[_id_vec[_sgr.id(a)]]; |
947 if (d < delta) delta = d; |
971 if (d < delta) delta = d; |
962 // Execute the "Cancel-and-Tighten" method |
986 // Execute the "Cancel-and-Tighten" method |
963 void startCancelAndTighten() { |
987 void startCancelAndTighten() { |
964 // Constants for the min mean cycle computations |
988 // Constants for the min mean cycle computations |
965 const double LIMIT_FACTOR = 1.0; |
989 const double LIMIT_FACTOR = 1.0; |
966 const int MIN_LIMIT = 5; |
990 const int MIN_LIMIT = 5; |
|
991 const double HW_ITER_LIMIT_FACTOR = 1.0; |
|
992 const int HW_ITER_LIMIT_MIN_VALUE = 5; |
|
993 |
|
994 const int hw_iter_limit = |
|
995 std::max(static_cast<int>(HW_ITER_LIMIT_FACTOR * _node_num), |
|
996 HW_ITER_LIMIT_MIN_VALUE); |
967 |
997 |
968 // Contruct auxiliary data vectors |
998 // Contruct auxiliary data vectors |
969 DoubleVector pi(_res_node_num, 0.0); |
999 DoubleVector pi(_res_node_num, 0.0); |
970 IntVector level(_res_node_num); |
1000 IntVector level(_res_node_num); |
971 BoolVector reached(_res_node_num); |
1001 BoolVector reached(_res_node_num); |
1135 curr = _cost[a] + pu - pi[_target[a]]; |
1165 curr = _cost[a] + pu - pi[_target[a]]; |
1136 if (-curr > epsilon) epsilon = -curr; |
1166 if (-curr > epsilon) epsilon = -curr; |
1137 } |
1167 } |
1138 } |
1168 } |
1139 } else { |
1169 } else { |
1140 typedef HowardMmc<StaticDigraph, CostArcMap> MMC; |
1170 typedef HowardMmc<StaticDigraph, CostArcMap> HwMmc; |
|
1171 typedef HartmannOrlinMmc<StaticDigraph, CostArcMap> HoMmc; |
1141 typedef typename BellmanFord<StaticDigraph, CostArcMap> |
1172 typedef typename BellmanFord<StaticDigraph, CostArcMap> |
1142 ::template SetDistMap<CostNodeMap>::Create BF; |
1173 ::template SetDistMap<CostNodeMap>::Create BF; |
1143 |
1174 |
1144 // Set epsilon to the minimum cycle mean |
1175 // Set epsilon to the minimum cycle mean |
|
1176 Cost cycle_cost = 0; |
|
1177 int cycle_size = 1; |
1145 buildResidualNetwork(); |
1178 buildResidualNetwork(); |
1146 MMC mmc(_sgr, _cost_map); |
1179 HwMmc hw_mmc(_sgr, _cost_map); |
1147 mmc.findCycleMean(); |
1180 if (hw_mmc.findCycleMean(hw_iter_limit) == HwMmc::ITERATION_LIMIT) { |
1148 epsilon = -mmc.cycleMean(); |
1181 // Howard's algorithm reached the iteration limit, start a |
1149 Cost cycle_cost = mmc.cycleCost(); |
1182 // strongly polynomial algorithm instead |
1150 int cycle_size = mmc.cycleSize(); |
1183 HoMmc ho_mmc(_sgr, _cost_map); |
|
1184 ho_mmc.findCycleMean(); |
|
1185 epsilon = -ho_mmc.cycleMean(); |
|
1186 cycle_cost = ho_mmc.cycleCost(); |
|
1187 cycle_size = ho_mmc.cycleSize(); |
|
1188 } else { |
|
1189 // Set epsilon |
|
1190 epsilon = -hw_mmc.cycleMean(); |
|
1191 cycle_cost = hw_mmc.cycleCost(); |
|
1192 cycle_size = hw_mmc.cycleSize(); |
|
1193 } |
1151 |
1194 |
1152 // Compute feasible potentials for the current epsilon |
1195 // Compute feasible potentials for the current epsilon |
1153 for (int i = 0; i != int(_cost_vec.size()); ++i) { |
1196 for (int i = 0; i != int(_cost_vec.size()); ++i) { |
1154 _cost_vec[i] = cycle_size * _cost_vec[i] - cycle_cost; |
1197 _cost_vec[i] = cycle_size * _cost_vec[i] - cycle_cost; |
1155 } |
1198 } |