Changeset 2556:74c2c81055e1 in lemon-0.x
- Timestamp:
- 01/13/08 11:32:14 (16 years ago)
- Branch:
- default
- Phase:
- public
- Convert:
- svn:c9d7d8f5-90d6-0310-b91f-818b3a526b0e/lemon/trunk@3437
- Location:
- lemon
- Files:
-
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
lemon/capacity_scaling.h
r2553 r2556 23 23 /// 24 24 /// \file 25 /// \brief The capacity scaling algorithm for finding a minimum cost 26 /// flow. 25 /// \brief The capacity scaling algorithm for finding a minimum cost flow. 27 26 28 27 #include <lemon/graph_adaptor.h> … … 50 49 /// 51 50 /// \warning 52 /// - Edge capacities and costs should be non negative integers.51 /// - Edge capacities and costs should be non-negative integers. 53 52 /// However \c CostMap::Value should be signed type. 54 53 /// - Supply values should be signed integers. … … 67 66 class CapacityScaling 68 67 { 69 typedef typename Graph::Node Node; 70 typedef typename Graph::NodeIt NodeIt; 71 typedef typename Graph::Edge Edge; 72 typedef typename Graph::EdgeIt EdgeIt; 73 typedef typename Graph::InEdgeIt InEdgeIt; 74 typedef typename Graph::OutEdgeIt OutEdgeIt; 68 GRAPH_TYPEDEFS(typename Graph); 75 69 76 70 typedef typename LowerMap::Value Lower; … … 78 72 typedef typename CostMap::Value Cost; 79 73 typedef typename SupplyMap::Value Supply; 80 typedef typename Graph::template EdgeMap<Capacity> Capacity RefMap;81 typedef typename Graph::template NodeMap<Supply> Supply RefMap;74 typedef typename Graph::template EdgeMap<Capacity> CapacityEdgeMap; 75 typedef typename Graph::template NodeMap<Supply> SupplyNodeMap; 82 76 typedef typename Graph::template NodeMap<Edge> PredMap; 83 77 84 78 public: 85 79 86 /// \briefType to enable or disable capacity scaling.80 /// Type to enable or disable capacity scaling. 87 81 enum ScalingEnum { 88 82 WITH_SCALING = 0, … … 90 84 }; 91 85 92 /// \briefThe type of the flow map.93 typedef CapacityRefMapFlowMap;94 /// \briefThe type of the potential map.86 /// The type of the flow map. 87 typedef typename Graph::template EdgeMap<Capacity> FlowMap; 88 /// The type of the potential map. 95 89 typedef typename Graph::template NodeMap<Cost> PotentialMap; 96 90 … … 111 105 protected: 112 106 113 /// \briefThe directed graph the algorithm runs on.107 /// The directed graph the algorithm runs on. 114 108 const Graph &graph; 115 109 116 /// \briefThe flow map.110 /// The flow map. 117 111 const FlowMap &flow; 118 /// \briefThe residual capacity map.119 const Capacity RefMap &res_cap;120 /// \briefThe cost map.112 /// The residual capacity map. 113 const CapacityEdgeMap &res_cap; 114 /// The cost map. 121 115 const CostMap &cost; 122 /// \briefThe excess map.123 const Supply RefMap &excess;124 125 /// \briefThe potential map.116 /// The excess map. 117 const SupplyNodeMap &excess; 118 119 /// The potential map. 126 120 PotentialMap &potential; 127 121 128 /// \briefThe distance map.122 /// The distance map. 129 123 CostNodeMap dist; 130 /// \briefThe map of predecessors edges.124 /// The map of predecessors edges. 131 125 PredMap &pred; 132 /// \briefThe processed (i.e. permanently labeled) nodes.126 /// The processed (i.e. permanently labeled) nodes. 133 127 std::vector<Node> proc_nodes; 134 128 135 129 public: 136 130 137 /// \briefThe constructor of the class.131 /// The constructor of the class. 138 132 ResidualDijkstra( const Graph &_graph, 139 133 const FlowMap &_flow, 140 const Capacity RefMap &_res_cap,134 const CapacityEdgeMap &_res_cap, 141 135 const CostMap &_cost, 142 136 const SupplyMap &_excess, … … 148 142 {} 149 143 150 /// \briefRuns the algorithm from the given source node.144 /// Runs the algorithm from the given source node. 151 145 Node run(Node s, Capacity delta) { 152 146 HeapCrossRef heap_cross_ref(graph, Heap::PRE_HEAP); … … 223 217 protected: 224 218 225 /// \briefThe directed graph the algorithm runs on.219 /// The directed graph the algorithm runs on. 226 220 const Graph &graph; 227 /// \briefThe original lower bound map.221 /// The original lower bound map. 228 222 const LowerMap *lower; 229 /// \briefThe modified capacity map.230 Capacity RefMap capacity;231 /// \briefThe cost map.223 /// The modified capacity map. 224 CapacityEdgeMap capacity; 225 /// The cost map. 232 226 const CostMap &cost; 233 /// \brief The modified supply map. 234 SupplyRefMap supply; 235 /// \brief The sum of supply values equals zero. 227 /// The modified supply map. 228 SupplyNodeMap supply; 236 229 bool valid_supply; 237 230 238 /// \briefThe edge map of the current flow.231 /// The edge map of the current flow. 239 232 FlowMap flow; 240 /// \briefThe potential node map.233 /// The potential node map. 241 234 PotentialMap potential; 242 235 243 /// \briefThe residual capacity map.244 Capacity RefMap res_cap;245 /// \briefThe excess map.246 Supply RefMap excess;247 /// \brief The excess nodes (i.e.nodes with positive excess).236 /// The residual capacity map. 237 CapacityEdgeMap res_cap; 238 /// The excess map. 239 SupplyNodeMap excess; 240 /// The excess nodes (i.e. the nodes with positive excess). 248 241 std::vector<Node> excess_nodes; 249 /// \brief The index of the next excess node.250 int next_node;251 252 /// \briefThe scaling status (enabled or disabled).242 /// The deficit nodes (i.e. the nodes with negative excess). 243 std::vector<Node> deficit_nodes; 244 245 /// The scaling status (enabled or disabled). 253 246 ScalingEnum scaling; 254 /// \brief Thedelta parameter used for capacity scaling.247 /// The \c delta parameter used for capacity scaling. 255 248 Capacity delta; 256 /// \brief The maximum number of phases. 257 Capacity phase_num; 258 /// \brief The deficit nodes. 259 std::vector<Node> deficit_nodes; 249 /// The maximum number of phases. 250 int phase_num; 260 251 261 252 /// \brief Implementation of the \ref Dijkstra algorithm for 262 253 /// finding augmenting shortest paths in the residual network. 263 254 ResidualDijkstra dijkstra; 264 /// \briefThe map of predecessors edges.255 /// The map of predecessors edges. 265 256 PredMap pred; 266 257 … … 286 277 dijkstra(graph, flow, res_cap, cost, excess, potential, pred) 287 278 { 288 // Removing non zero lower bounds279 // Removing non-zero lower bounds 289 280 capacity = subMap(_capacity, _lower); 290 281 res_cap = capacity; … … 336 327 /// \param _t The target node. 337 328 /// \param _flow_value The required amount of flow from node \c _s 338 /// to node \c _t (i.e. the supply of \c _s and the demand of 339 /// \c _t). 329 /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t). 340 330 CapacityScaling( const Graph &_graph, 341 331 const LowerMap &_lower, … … 349 339 dijkstra(graph, flow, res_cap, cost, excess, potential) 350 340 { 351 // Removing non zero lower bounds341 // Removing non-zero lower bounds 352 342 capacity = subMap(_capacity, _lower); 353 343 res_cap = capacity; … … 375 365 /// \param _t The target node. 376 366 /// \param _flow_value The required amount of flow from node \c _s 377 /// to node \c _t (i.e. the supply of \c _s and the demand of 378 /// \c _t). 367 /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t). 379 368 CapacityScaling( const Graph &_graph, 380 369 const CapacityMap &_capacity, … … 392 381 } 393 382 383 /// \brief Runs the algorithm. 384 /// 385 /// Runs the algorithm. 386 /// 387 /// \param scaling_mode The scaling mode. In case of WITH_SCALING 388 /// capacity scaling is enabled in the algorithm (this is the 389 /// default) otherwise it is disabled. 390 /// If the maximum edge capacity and/or the amount of total supply 391 /// is small, the algorithm could be slightly faster without 392 /// scaling. 393 /// 394 /// \return \c true if a feasible flow can be found. 395 bool run(int scaling_mode = WITH_SCALING) { 396 return init(scaling_mode) && start(); 397 } 398 394 399 /// \brief Returns a const reference to the flow map. 395 400 /// … … 425 430 } 426 431 427 /// \brief Runs the algorithm.428 ///429 /// Runs the algorithm.430 ///431 /// \param scaling_mode The scaling mode. In case of WITH_SCALING432 /// capacity scaling is enabled in the algorithm (this is the433 /// default value) otherwise it is disabled.434 /// If the maximum edge capacity and/or the amount of total supply435 /// is small, the algorithm could be faster without scaling.436 ///437 /// \return \c true if a feasible flow can be found.438 bool run(int scaling_mode = WITH_SCALING) {439 return init(scaling_mode) && start();440 }441 442 432 protected: 443 433 444 /// \briefInitializes the algorithm.434 /// Initializes the algorithm. 445 435 bool init(int scaling_mode) { 446 436 if (!valid_supply) return false; … … 466 456 } 467 457 468 /// \briefExecutes the algorithm.458 /// Executes the algorithm. 469 459 bool start() { 470 460 if (delta > 1) … … 507 497 if (excess[n] <= -delta) deficit_nodes.push_back(n); 508 498 } 509 next_node = 0;499 int next_node = 0; 510 500 511 501 // Finding augmenting shortest paths … … 573 563 } 574 564 575 // Handling non zero lower bounds565 // Handling non-zero lower bounds 576 566 if (lower) { 577 567 for (EdgeIt e(graph); e != INVALID; ++e) … … 585 575 bool startWithoutScaling() { 586 576 // Finding excess nodes 587 for (NodeIt n(graph); n != INVALID; ++n) {577 for (NodeIt n(graph); n != INVALID; ++n) 588 578 if (excess[n] > 0) excess_nodes.push_back(n); 589 }590 579 if (excess_nodes.size() == 0) return true; 591 next_node = 0;580 int next_node = 0; 592 581 593 582 // Finding shortest paths … … 632 621 } 633 622 634 // Handling non zero lower bounds623 // Handling non-zero lower bounds 635 624 if (lower) { 636 625 for (EdgeIt e(graph); e != INVALID; ++e) -
lemon/cycle_canceling.h
r2553 r2556 35 35 #ifdef LIMITED_CYCLE_CANCELING 36 36 #include <lemon/bellman_ford.h> 37 /// \brief The maximum number of iterations for the first execution 38 /// of the \ref lemon::BellmanFord "Bellman-Ford" algorithm. 39 /// It should be at least 2. 40 #define STARTING_LIMIT 2 41 /// \brief The iteration limit for the 42 /// \ref lemon::BellmanFord "Bellman-Ford" algorithm is multiplied by 43 /// <tt>ALPHA_MUL / ALPHA_DIV</tt> in every round. 44 /// <tt>ALPHA_MUL / ALPHA_DIV</tt> must be greater than 1. 45 #define ALPHA_MUL 3 46 /// \brief The iteration limit for the 47 /// \ref lemon::BellmanFord "Bellman-Ford" algorithm is multiplied by 48 /// <tt>ALPHA_MUL / ALPHA_DIV</tt> in every round. 49 /// <tt>ALPHA_MUL / ALPHA_DIV</tt> must be greater than 1. 50 #define ALPHA_DIV 2 37 // The maximum number of iterations for the first execution of the 38 // Bellman-Ford algorithm. It should be at least 2. 39 #define STARTING_LIMIT 2 40 // The iteration limit for the Bellman-Ford algorithm is multiplied by 41 // <tt>ALPHA_MUL / ALPHA_DIV</tt> in every round. 42 // <tt>ALPHA_MUL / ALPHA_DIV</tt> must be greater than 1. 43 #define ALPHA_MUL 3 44 #define ALPHA_DIV 2 51 45 52 46 //#define _ONLY_ONE_CYCLE_ 53 47 //#define _NO_BACK_STEP_ 54 //#define _DEBUG_ITER_55 48 #endif 56 49 … … 60 53 #endif 61 54 55 //#define _DEBUG_ITER_ 56 62 57 namespace lemon { 63 58 … … 65 60 /// @{ 66 61 67 /// \brief Implementation of a cycle-canceling algorithm for finding68 /// a minimum cost flow.62 /// \brief Implementation of a cycle-canceling algorithm for 63 /// finding a minimum cost flow. 69 64 /// 70 /// \ref lemon::CycleCanceling "CycleCanceling" implements a71 /// cycle-canceling algorithm forfinding a minimum cost flow.65 /// \ref CycleCanceling implements a cycle-canceling algorithm for 66 /// finding a minimum cost flow. 72 67 /// 73 68 /// \param Graph The directed graph type the algorithm runs on. … … 78 73 /// 79 74 /// \warning 80 /// - Edge capacities and costs should be non negative integers.81 /// 75 /// - Edge capacities and costs should be non-negative integers. 76 /// However \c CostMap::Value should be signed type. 82 77 /// - Supply values should be signed integers. 83 78 /// - \c LowerMap::Value must be convertible to 84 /// 85 /// 79 /// \c CapacityMap::Value and \c CapacityMap::Value must be 80 /// convertible to \c SupplyMap::Value. 86 81 /// 87 82 /// \author Peter Kovacs … … 95 90 class CycleCanceling 96 91 { 97 typedef typename Graph::Node Node; 98 typedef typename Graph::NodeIt NodeIt; 99 typedef typename Graph::Edge Edge; 100 typedef typename Graph::EdgeIt EdgeIt; 101 typedef typename Graph::InEdgeIt InEdgeIt; 102 typedef typename Graph::OutEdgeIt OutEdgeIt; 92 GRAPH_TYPEDEFS(typename Graph); 103 93 104 94 typedef typename LowerMap::Value Lower; … … 106 96 typedef typename CostMap::Value Cost; 107 97 typedef typename SupplyMap::Value Supply; 108 typedef typename Graph::template EdgeMap<Capacity> Capacity RefMap;109 typedef typename Graph::template NodeMap<Supply> Supply RefMap;98 typedef typename Graph::template EdgeMap<Capacity> CapacityEdgeMap; 99 typedef typename Graph::template NodeMap<Supply> SupplyNodeMap; 110 100 111 101 typedef ResGraphAdaptor< const Graph, Capacity, 112 CapacityRefMap, CapacityRefMap > ResGraph;102 CapacityEdgeMap, CapacityEdgeMap > ResGraph; 113 103 typedef typename ResGraph::Node ResNode; 114 104 typedef typename ResGraph::NodeIt ResNodeIt; … … 118 108 public: 119 109 120 /// \briefThe type of the flow map.121 typedef CapacityRefMapFlowMap;110 /// The type of the flow map. 111 typedef typename Graph::template EdgeMap<Capacity> FlowMap; 122 112 123 113 protected: 124 114 125 /// \briefMap adaptor class for handling residual edge costs.115 /// Map adaptor class for handling residual edge costs. 126 116 class ResCostMap : public MapBase<ResEdge, Cost> 127 117 { … … 135 125 136 126 Cost operator[](const ResEdge &e) const { 137 127 return ResGraph::forward(e) ? cost_map[e] : -cost_map[e]; 138 128 } 139 129 … … 142 132 protected: 143 133 144 /// \briefThe directed graph the algorithm runs on.134 /// The directed graph the algorithm runs on. 145 135 const Graph &graph; 146 /// \briefThe original lower bound map.136 /// The original lower bound map. 147 137 const LowerMap *lower; 148 /// \briefThe modified capacity map.149 Capacity RefMap capacity;150 /// \briefThe cost map.138 /// The modified capacity map. 139 CapacityEdgeMap capacity; 140 /// The cost map. 151 141 const CostMap &cost; 152 /// \brief The modified supply map. 153 SupplyRefMap supply; 154 /// \brief The sum of supply values equals zero. 142 /// The modified supply map. 143 SupplyNodeMap supply; 155 144 bool valid_supply; 156 145 157 /// \briefThe current flow.146 /// The current flow. 158 147 FlowMap flow; 159 /// \briefThe residual graph.148 /// The residual graph. 160 149 ResGraph res_graph; 161 /// \briefThe residual cost map.150 /// The residual cost map. 162 151 ResCostMap res_cost; 163 152 … … 174 163 /// \param _supply The supply values of the nodes (signed). 175 164 CycleCanceling( const Graph &_graph, 176 177 178 179 165 const LowerMap &_lower, 166 const CapacityMap &_capacity, 167 const CostMap &_cost, 168 const SupplyMap &_supply ) : 180 169 graph(_graph), lower(&_lower), capacity(_graph), cost(_cost), 181 170 supply(_graph), flow(_graph, 0), 182 171 res_graph(_graph, capacity, flow), res_cost(cost) 183 172 { 184 // Removing non zero lower bounds173 // Removing non-zero lower bounds 185 174 capacity = subMap(_capacity, _lower); 186 175 Supply sum = 0; 187 176 for (NodeIt n(graph); n != INVALID; ++n) { 188 189 190 191 192 193 177 Supply s = _supply[n]; 178 for (InEdgeIt e(graph, n); e != INVALID; ++e) 179 s += _lower[e]; 180 for (OutEdgeIt e(graph, n); e != INVALID; ++e) 181 s -= _lower[e]; 182 sum += (supply[n] = s); 194 183 } 195 184 valid_supply = sum == 0; … … 205 194 /// \param _supply The supply values of the nodes (signed). 206 195 CycleCanceling( const Graph &_graph, 207 208 209 196 const CapacityMap &_capacity, 197 const CostMap &_cost, 198 const SupplyMap &_supply ) : 210 199 graph(_graph), lower(NULL), capacity(_capacity), cost(_cost), 211 200 supply(_supply), flow(_graph, 0), … … 218 207 } 219 208 220 221 209 /// \brief Simple constructor of the class (with lower bounds). 222 210 /// … … 230 218 /// \param _t The target node. 231 219 /// \param _flow_value The required amount of flow from node \c _s 232 /// to node \c _t (i.e. the supply of \c _s and the demand of 233 /// \c _t). 220 /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t). 234 221 CycleCanceling( const Graph &_graph, 235 236 237 238 239 222 const LowerMap &_lower, 223 const CapacityMap &_capacity, 224 const CostMap &_cost, 225 Node _s, Node _t, 226 Supply _flow_value ) : 240 227 graph(_graph), lower(&_lower), capacity(_graph), cost(_cost), 241 228 supply(_graph), flow(_graph, 0), 242 229 res_graph(_graph, capacity, flow), res_cost(cost) 243 230 { 244 // Removing non zero lower bounds231 // Removing non-zero lower bounds 245 232 capacity = subMap(_capacity, _lower); 246 233 for (NodeIt n(graph); n != INVALID; ++n) { 247 248 249 250 251 252 253 254 234 Supply s = 0; 235 if (n == _s) s = _flow_value; 236 if (n == _t) s = -_flow_value; 237 for (InEdgeIt e(graph, n); e != INVALID; ++e) 238 s += _lower[e]; 239 for (OutEdgeIt e(graph, n); e != INVALID; ++e) 240 s -= _lower[e]; 241 supply[n] = s; 255 242 } 256 243 valid_supply = true; … … 267 254 /// \param _t The target node. 268 255 /// \param _flow_value The required amount of flow from node \c _s 269 /// to node \c _t (i.e. the supply of \c _s and the demand of 270 /// \c _t). 256 /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t). 271 257 CycleCanceling( const Graph &_graph, 272 273 274 275 258 const CapacityMap &_capacity, 259 const CostMap &_cost, 260 Node _s, Node _t, 261 Supply _flow_value ) : 276 262 graph(_graph), lower(NULL), capacity(_capacity), cost(_cost), 277 263 supply(_graph, 0), flow(_graph, 0), … … 283 269 } 284 270 271 /// \brief Runs the algorithm. 272 /// 273 /// Runs the algorithm. 274 /// 275 /// \return \c true if a feasible flow can be found. 276 bool run() { 277 return init() && start(); 278 } 279 285 280 /// \brief Returns a const reference to the flow map. 286 281 /// … … 301 296 Cost c = 0; 302 297 for (EdgeIt e(graph); e != INVALID; ++e) 303 298 c += flow[e] * cost[e]; 304 299 return c; 305 300 } 306 301 307 /// \brief Runs the algorithm.308 ///309 /// Runs the algorithm.310 ///311 /// \return \c true if a feasible flow can be found.312 bool run() {313 return init() && start();314 }315 316 302 protected: 317 303 318 /// \briefInitializes the algorithm.304 /// Initializes the algorithm. 319 305 bool init() { 320 306 // Checking the sum of supply values … … 324 310 325 311 // Finding a feasible flow 326 Circulation< Graph, ConstMap<Edge, Capacity>, Capacity RefMap,327 328 329 312 Circulation< Graph, ConstMap<Edge, Capacity>, CapacityEdgeMap, 313 SupplyMap > 314 circulation( graph, constMap<Edge>((Capacity)0), capacity, 315 supply ); 330 316 circulation.flowMap(flow); 331 317 return circulation.run(); … … 334 320 #ifdef LIMITED_CYCLE_CANCELING 335 321 /// \brief Executes a cycle-canceling algorithm using 336 /// \ref lemon::BellmanFord "Bellman-Ford" algorithm with limited 337 /// iteration count. 322 /// \ref Bellman-Ford algorithm with limited iteration count. 338 323 bool start() { 339 324 typename BellmanFord<ResGraph, ResCostMap>::PredMap pred(res_graph); … … 348 333 bool optimal = false; 349 334 while (!optimal) { 350 351 352 353 354 355 335 BellmanFord<ResGraph, ResCostMap> bf(res_graph, res_cost); 336 bf.predMap(pred); 337 bf.init(0); 338 int iter_num = 0; 339 bool cycle_found = false; 340 while (!cycle_found) { 356 341 #ifdef _NO_BACK_STEP_ 357 358 342 int curr_iter_num = length_bound <= node_num ? 343 length_bound - iter_num : node_num - iter_num; 359 344 #else 360 361 362 #endif 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 #ifdef _DEBUG_ITER_ 402 403 #endif 404 405 406 345 int curr_iter_num = iter_num + length_bound <= node_num ? 346 length_bound : node_num - iter_num; 347 #endif 348 iter_num += curr_iter_num; 349 int real_iter_num = curr_iter_num; 350 for (int i = 0; i < curr_iter_num; ++i) { 351 if (bf.processNextWeakRound()) { 352 real_iter_num = i; 353 break; 354 } 355 } 356 if (real_iter_num < curr_iter_num) { 357 optimal = true; 358 break; 359 } else { 360 // Searching for node disjoint negative cycles 361 for (ResNodeIt n(res_graph); n != INVALID; ++n) 362 visited[n] = 0; 363 int id = 0; 364 for (ResNodeIt n(res_graph); n != INVALID; ++n) { 365 if (visited[n] > 0) continue; 366 visited[n] = ++id; 367 ResNode u = pred[n] == INVALID ? 368 INVALID : res_graph.source(pred[n]); 369 while (u != INVALID && visited[u] == 0) { 370 visited[u] = id; 371 u = pred[u] == INVALID ? 372 INVALID : res_graph.source(pred[u]); 373 } 374 if (u != INVALID && visited[u] == id) { 375 // Finding the negative cycle 376 cycle_found = true; 377 cycle.clear(); 378 ResEdge e = pred[u]; 379 cycle.push_back(e); 380 Capacity d = res_graph.rescap(e); 381 while (res_graph.source(e) != u) { 382 cycle.push_back(e = pred[res_graph.source(e)]); 383 if (res_graph.rescap(e) < d) 384 d = res_graph.rescap(e); 385 } 386 #ifdef _DEBUG_ITER_ 387 ++cycle_num; 388 #endif 389 // Augmenting along the cycle 390 for (int i = 0; i < cycle.size(); ++i) 391 res_graph.augment(cycle[i], d); 407 392 #ifdef _ONLY_ONE_CYCLE_ 408 409 #endif 410 411 412 413 414 415 416 393 break; 394 #endif 395 } 396 } 397 } 398 399 if (!cycle_found) 400 length_bound = length_bound * ALPHA_MUL / ALPHA_DIV; 401 } 417 402 } 418 403 419 404 #ifdef _DEBUG_ITER_ 420 405 std::cout << "Limited cycle-canceling algorithm finished. " 421 422 423 #endif 424 425 // Handling non zero lower bounds406 << "Found " << cycle_num << " negative cycles." 407 << std::endl; 408 #endif 409 410 // Handling non-zero lower bounds 426 411 if (lower) { 427 428 412 for (EdgeIt e(graph); e != INVALID; ++e) 413 flow[e] += (*lower)[e]; 429 414 } 430 415 return true; … … 434 419 #ifdef MIN_MEAN_CYCLE_CANCELING 435 420 /// \brief Executes the minimum mean cycle-canceling algorithm 436 /// using \ref lemon::MinMeanCycle "MinMeanCycle" class.421 /// using \ref MinMeanCycle. 437 422 bool start() { 438 423 typedef Path<ResGraph> ResPath; … … 445 430 mmc.cyclePath(cycle).init(); 446 431 if (mmc.findMinMean()) { 447 448 #ifdef _DEBUG_ITER_ 449 ++iter;450 #endif 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 432 while (mmc.cycleLength() < 0) { 433 #ifdef _DEBUG_ITER_ 434 ++cycle_num; 435 #endif 436 // Finding the cycle 437 mmc.findCycle(); 438 439 // Finding the largest flow amount that can be augmented 440 // along the cycle 441 Capacity delta = 0; 442 for (typename ResPath::EdgeIt e(cycle); e != INVALID; ++e) { 443 if (delta == 0 || res_graph.rescap(e) < delta) 444 delta = res_graph.rescap(e); 445 } 446 447 // Augmenting along the cycle 448 for (typename ResPath::EdgeIt e(cycle); e != INVALID; ++e) 449 res_graph.augment(e, delta); 450 451 // Finding the minimum cycle mean for the modified residual 452 // graph 453 mmc.reset(); 454 if (!mmc.findMinMean()) break; 455 } 471 456 } 472 457 473 458 #ifdef _DEBUG_ITER_ 474 459 std::cout << "Minimum mean cycle-canceling algorithm finished. " 475 476 477 #endif 478 479 // Handling non zero lower bounds460 << "Found " << cycle_num << " negative cycles." 461 << std::endl; 462 #endif 463 464 // Handling non-zero lower bounds 480 465 if (lower) { 481 482 466 for (EdgeIt e(graph); e != INVALID; ++e) 467 flow[e] += (*lower)[e]; 483 468 } 484 469 return true; -
lemon/min_cost_flow.h
r2553 r2556 34 34 /// \brief An efficient algorithm for finding a minimum cost flow. 35 35 /// 36 /// \ref lemon::MinCostFlow "MinCostFlow" implements an efficient37 /// a lgorithm for finding aminimum cost flow.36 /// \ref MinCostFlow provides an efficient algorithm for finding 37 /// a minimum cost flow. 38 38 /// 39 /// \note \ref lemon::MinCostFlow "MinCostFlow" is just an alias for40 /// \ref lemon::NetworkSimplex "NetworkSimplex", wich is the most41 /// efficient implementation for the minimum cost flow problem in the42 /// LEMON library according to our benchmarktests.39 /// \note \ref MinCostFlow is just an alias for \ref NetworkSimplex, 40 /// which is the most efficient implementation for the minimum cost 41 /// flow problem in the LEMON library according to our benchmark 42 /// tests. 43 43 /// 44 44 /// \note There are three implementations for the minimum cost flow 45 /// problem, that can be used in the same way. 46 /// - \ref lemon::CapacityScaling "CapacityScaling" 47 /// - \ref lemon::CycleCanceling "CycleCanceling" 48 /// - \ref lemon::NetworkSimplex "NetworkSimplex" 45 /// problem, that can be used exactly the same way. 46 /// - \ref CapacityScaling 47 /// - \ref CycleCanceling 48 /// - \ref NetworkSimplex 49 /// 50 /// \note For the detailed documentation of this class see 51 /// \ref NetworkSimplex. 49 52 /// 50 53 /// \param Graph The directed graph type the algorithm runs on. … … 55 58 /// 56 59 /// \warning 57 /// - Edge capacities and costs should be non negative integers.58 /// 60 /// - Edge capacities and costs should be non-negative integers. 61 /// However \c CostMap::Value should be signed type. 59 62 /// - Supply values should be signed integers. 60 63 /// - \c LowerMap::Value must be convertible to 61 /// 62 /// 64 /// \c CapacityMap::Value and \c CapacityMap::Value must be 65 /// convertible to \c SupplyMap::Value. 63 66 /// 64 67 /// \author Peter Kovacs … … 71 74 <typename CapacityMap::Value> > 72 75 class MinCostFlow : 73 public NetworkSimplex< Graph, 74 LowerMap, 75 CapacityMap, 76 CostMap, 77 SupplyMap > 76 public NetworkSimplex< Graph, LowerMap, CapacityMap, 77 CostMap, SupplyMap > 78 78 { 79 79 typedef NetworkSimplex< Graph, LowerMap, CapacityMap, 80 CostMap, SupplyMap > 81 MinCostFlowImpl; 80 CostMap, SupplyMap > MinCostFlowImpl; 82 81 typedef typename Graph::Node Node; 83 82 typedef typename SupplyMap::Value Supply; … … 85 84 public: 86 85 87 /// \briefGeneral constructor of the class (with lower bounds).88 MinCostFlow( const Graph & _graph,89 const LowerMap &_lower,90 const CapacityMap &_capacity,91 const CostMap &_cost,92 const SupplyMap &_supply ) :93 MinCostFlowImpl( _graph, _lower, _capacity, _cost, _supply) {}86 /// General constructor of the class (with lower bounds). 87 MinCostFlow( const Graph &graph, 88 const LowerMap &lower, 89 const CapacityMap &capacity, 90 const CostMap &cost, 91 const SupplyMap &supply ) : 92 MinCostFlowImpl(graph, lower, capacity, cost, supply) {} 94 93 95 /// \briefGeneral constructor of the class (without lower bounds).96 MinCostFlow( const Graph & _graph,97 const CapacityMap &_capacity,98 const CostMap &_cost,99 const SupplyMap &_supply ) :100 MinCostFlowImpl( _graph, _capacity, _cost, _supply) {}94 /// General constructor of the class (without lower bounds). 95 MinCostFlow( const Graph &graph, 96 const CapacityMap &capacity, 97 const CostMap &_ost, 98 const SupplyMap &supply ) : 99 MinCostFlowImpl(graph, capacity, cost, supply) {} 101 100 102 /// \briefSimple constructor of the class (with lower bounds).103 MinCostFlow( const Graph & _graph,104 const LowerMap &_lower,105 const CapacityMap &_capacity,106 const CostMap &_cost,107 Node _s, Node _t,108 Supply _flow_value ) :109 MinCostFlowImpl( _graph, _lower, _capacity, _cost,110 _s, _t, _flow_value ) {}101 /// Simple constructor of the class (with lower bounds). 102 MinCostFlow( const Graph &graph, 103 const LowerMap &lower, 104 const CapacityMap &capacity, 105 const CostMap &_ost, 106 Node s, Node t, 107 Supply flow_value ) : 108 MinCostFlowImpl( graph, lower, capacity, cost, 109 s, t, flow_value ) {} 111 110 112 /// \briefSimple constructor of the class (without lower bounds).113 MinCostFlow( const Graph & _graph,114 const CapacityMap &_capacity,115 const CostMap &_cost,116 Node _s, Node _t,117 Supply _flow_value ) :118 MinCostFlowImpl( _graph, _capacity, _cost,119 _s, _t, _flow_value ) {}111 /// Simple constructor of the class (without lower bounds). 112 MinCostFlow( const Graph &graph, 113 const CapacityMap &capacity, 114 const CostMap &cost, 115 Node s, Node t, 116 Supply flow_value ) : 117 MinCostFlowImpl( graph, capacity, cost, 118 s, t, flow_value ) {} 120 119 121 120 }; //class MinCostFlow -
lemon/min_cost_max_flow.h
r2553 r2556 23 23 /// 24 24 /// \file 25 /// \brief An efficient algorithm for finding a minimum cost maximum 26 /// flow. 25 /// \brief An efficient algorithm for finding a minimum cost maximum flow. 27 26 28 27 #include <lemon/preflow.h> … … 35 34 /// @{ 36 35 37 /// \brief An efficient algorithm for finding a minimum cost maximum38 /// flow.36 /// \brief An efficient algorithm for finding a minimum cost 37 /// maximum flow. 39 38 /// 40 /// \ref lemon::MinCostFlow "MinCostMaxFlow" implements an efficient 41 /// algorithm for finding a maximum flow having minimal total cost 42 /// from a given source node to a given target node in a directed 43 /// graph. 39 /// \ref MinCostMaxFlow implements an efficient algorithm for 40 /// finding a maximum flow having minimal total cost from a given 41 /// source node to a given target node in a directed graph. 44 42 /// 45 /// \note \ref lemon::MinCostMaxFlow "MinCostMaxFlow" uses 46 /// \ref lemon::Preflow "Preflow" algorithm for finding the maximum 47 /// flow value and \ref lemon::NetworkSimplex "NetworkSimplex" 48 /// algorithm for finding a minimum cost flow of that value. 43 /// \note \ref MinCostMaxFlow uses \ref Preflow algorithm for finding 44 /// the maximum flow value and \ref NetworkSimplex algorithm for 45 /// finding a minimum cost flow of that value. 49 46 /// 50 47 /// \param Graph The directed graph type the algorithm runs on. … … 53 50 /// 54 51 /// \warning 55 /// - Edge capacities and costs should be non negative integers.56 /// 52 /// - Edge capacities and costs should be non-negative integers. 53 /// However \c CostMap::Value should be signed type. 57 54 /// 58 55 /// \author Peter Kovacs 59 56 60 57 template < typename Graph, 61 62 58 typename CapacityMap = typename Graph::template EdgeMap<int>, 59 typename CostMap = typename Graph::template EdgeMap<int> > 63 60 class MinCostMaxFlow 64 61 { … … 67 64 68 65 typedef typename CapacityMap::Value Capacity; 66 typedef typename CostMap::Value Cost; 69 67 typedef typename Graph::template NodeMap<Capacity> SupplyMap; 70 68 typedef NetworkSimplex< Graph, CapacityMap, CapacityMap, 71 69 CostMap, SupplyMap > 72 70 MinCostFlowImpl; 73 71 74 72 public: 75 73 76 /// \briefThe type of the flow map.74 /// The type of the flow map. 77 75 typedef typename Graph::template EdgeMap<Capacity> FlowMap; 78 typedef typename CostMap::Value Cost;79 76 80 77 private: 81 78 82 /// \briefThe directed graph the algorithm runs on.79 /// The directed graph the algorithm runs on. 83 80 const Graph &graph; 84 /// \briefThe modified capacity map.81 /// The modified capacity map. 85 82 const CapacityMap &capacity; 86 /// \briefThe cost map.83 /// The cost map. 87 84 const CostMap &cost; 88 /// \brief The source node. 85 /// The edge map of the found flow. 86 FlowMap flow; 87 /// The source node. 89 88 Node source; 90 /// \briefThe target node.89 /// The target node. 91 90 Node target; 92 /// \brief The edge map of the found flow.93 FlowMap flow;94 91 95 typedef Preflow<Graph, CapacityMap> PreflowImpl; 96 /// \brief \ref lemon::Preflow "Preflow" class for finding the 97 /// maximum flow value. 98 PreflowImpl preflow; 92 typedef Preflow<Graph, CapacityMap> MaxFlowImpl; 93 /// \brief \ref Preflow class for finding the maximum flow value. 94 MaxFlowImpl preflow; 99 95 100 96 public: … … 110 106 /// \param _t The target node. 111 107 MinCostMaxFlow( const Graph &_graph, 112 113 114 108 const CapacityMap &_capacity, 109 const CostMap &_cost, 110 Node _s, Node _t ) : 115 111 graph(_graph), capacity(_capacity), cost(_cost), 116 112 source(_s), target(_t), flow(_graph), … … 136 132 Cost c = 0; 137 133 for (typename Graph::EdgeIt e(graph); e != INVALID; ++e) 138 134 c += flow[e] * cost[e]; 139 135 return c; 140 136 } … … 145 141 preflow.runMinCut(); 146 142 MinCostFlowImpl mcf_impl( graph, capacity, cost, 147 143 source, target, preflow.flowValue() ); 148 144 mcf_impl.run(); 149 145 flow = mcf_impl.flowMap(); -
lemon/network_simplex.h
r2553 r2556 23 23 /// 24 24 /// \file 25 /// \brief The network simplex algorithm for finding a minimum cost 26 /// flow. 25 /// \brief The network simplex algorithm for finding a minimum cost flow. 27 26 28 27 #include <limits> … … 41 40 42 41 43 // / \briefState constant for edges at their lower bounds.44 #define LOWER 45 // / \briefState constant for edges in the spanning tree.46 #define TREE 47 // / \briefState constant for edges at their upper bounds.48 #define UPPER 42 // State constant for edges at their lower bounds. 43 #define LOWER 1 44 // State constant for edges in the spanning tree. 45 #define TREE 0 46 // State constant for edges at their upper bounds. 47 #define UPPER -1 49 48 50 49 #ifdef EDGE_BLOCK_PIVOT 51 50 #include <cmath> 52 /// \brief Lower bound for the size of blocks. 53 #define MIN_BLOCK_SIZE 10 51 #define MIN_BLOCK_SIZE 10 54 52 #endif 55 53 56 54 #ifdef CANDIDATE_LIST_PIVOT 57 55 #include <vector> 58 #define LIST_LENGTH_DIV 59 #define MINOR_LIMIT_DIV 56 #define LIST_LENGTH_DIV 50 57 #define MINOR_LIMIT_DIV 200 60 58 #endif 61 59 … … 64 62 #include <algorithm> 65 63 #define LIST_LENGTH_DIV 100 66 #define LOWER_DIV 64 #define LOWER_DIV 4 67 65 #endif 68 66 … … 75 73 /// finding a minimum cost flow. 76 74 /// 77 /// \ref lemon::NetworkSimplex "NetworkSimplex" implements the78 /// network simplex algorithm forfinding a minimum cost flow.75 /// \ref NetworkSimplex implements the network simplex algorithm for 76 /// finding a minimum cost flow. 79 77 /// 80 78 /// \param Graph The directed graph type the algorithm runs on. … … 85 83 /// 86 84 /// \warning 87 /// - Edge capacities and costs should be non negative integers.88 /// 85 /// - Edge capacities and costs should be non-negative integers. 86 /// However \c CostMap::Value should be signed type. 89 87 /// - Supply values should be signed integers. 90 88 /// - \c LowerMap::Value must be convertible to 91 /// 92 /// 89 /// \c CapacityMap::Value and \c CapacityMap::Value must be 90 /// convertible to \c SupplyMap::Value. 93 91 /// 94 92 /// \author Peter Kovacs … … 108 106 109 107 typedef SmartGraph SGraph; 110 typedef typename SGraph::Node Node; 111 typedef typename SGraph::NodeIt NodeIt; 112 typedef typename SGraph::Edge Edge; 113 typedef typename SGraph::EdgeIt EdgeIt; 114 typedef typename SGraph::InEdgeIt InEdgeIt; 115 typedef typename SGraph::OutEdgeIt OutEdgeIt; 108 GRAPH_TYPEDEFS(typename SGraph); 116 109 117 110 typedef typename SGraph::template EdgeMap<Lower> SLowerMap; … … 132 125 public: 133 126 134 /// \briefThe type of the flow map.127 /// The type of the flow map. 135 128 typedef typename Graph::template EdgeMap<Capacity> FlowMap; 136 /// \briefThe type of the potential map.129 /// The type of the potential map. 137 130 typedef typename Graph::template NodeMap<Cost> PotentialMap; 138 131 139 132 protected: 140 133 141 /// \briefMap adaptor class for handling reduced edge costs.134 /// Map adaptor class for handling reduced edge costs. 142 135 class ReducedCostMap : public MapBase<Edge, Cost> 143 136 { … … 151 144 152 145 ReducedCostMap( const SGraph &_gr, 153 154 155 146 const SCostMap &_cm, 147 const SPotentialMap &_pm ) : 148 gr(_gr), cost_map(_cm), pot_map(_pm) {} 156 149 157 150 Cost operator[](const Edge &e) const { 158 151 return cost_map[e] - pot_map[gr.source(e)] + pot_map[gr.target(e)]; 159 152 } 160 153 … … 163 156 protected: 164 157 165 /// \briefThe directed graph the algorithm runs on.158 /// The directed graph the algorithm runs on. 166 159 SGraph graph; 167 /// \briefThe original graph.160 /// The original graph. 168 161 const Graph &graph_ref; 169 /// \briefThe original lower bound map.162 /// The original lower bound map. 170 163 const LowerMap *lower; 171 /// \briefThe capacity map.164 /// The capacity map. 172 165 SCapacityMap capacity; 173 /// \briefThe cost map.166 /// The cost map. 174 167 SCostMap cost; 175 /// \briefThe supply map.168 /// The supply map. 176 169 SSupplyMap supply; 177 /// \briefThe reduced cost map.170 /// The reduced cost map. 178 171 ReducedCostMap red_cost; 179 /// \brief The sum of supply values equals zero.180 172 bool valid_supply; 181 173 182 /// \briefThe edge map of the current flow.174 /// The edge map of the current flow. 183 175 SCapacityMap flow; 184 /// \brief The edge map of the found flow on the original graph. 176 /// The potential node map. 177 SPotentialMap potential; 185 178 FlowMap flow_result; 186 /// \brief The potential node map.187 SPotentialMap potential;188 /// \brief The potential node map on the original graph.189 179 PotentialMap potential_result; 190 180 191 /// \briefNode reference for the original graph.181 /// Node reference for the original graph. 192 182 NodeRefMap node_ref; 193 /// \briefEdge reference for the original graph.183 /// Edge reference for the original graph. 194 184 EdgeRefMap edge_ref; 195 185 196 /// \brief Thedepth node map of the spanning tree structure.186 /// The \c depth node map of the spanning tree structure. 197 187 IntNodeMap depth; 198 /// \brief Theparent node map of the spanning tree structure.188 /// The \c parent node map of the spanning tree structure. 199 189 NodeNodeMap parent; 200 /// \brief Thepred_edge node map of the spanning tree structure.190 /// The \c pred_edge node map of the spanning tree structure. 201 191 EdgeNodeMap pred_edge; 202 /// \brief Thethread node map of the spanning tree structure.192 /// The \c thread node map of the spanning tree structure. 203 193 NodeNodeMap thread; 204 /// \brief Theforward node map of the spanning tree structure.194 /// The \c forward node map of the spanning tree structure. 205 195 BoolNodeMap forward; 206 /// \briefThe state edge map.196 /// The state edge map. 207 197 IntEdgeMap state; 208 198 /// The root node of the starting spanning tree. 199 Node root; 200 201 // The entering edge of the current pivot iteration. 202 Edge in_edge; 203 // Temporary nodes used in the current pivot iteration. 204 Node join, u_in, v_in, u_out, v_out; 205 Node right, first, second, last; 206 Node stem, par_stem, new_stem; 207 // The maximum augment amount along the found cycle in the current 208 // pivot iteration. 209 Capacity delta; 209 210 210 211 #ifdef EDGE_BLOCK_PIVOT 211 /// \briefThe size of blocks for the "Edge Block" pivot rule.212 /// The size of blocks for the "Edge Block" pivot rule. 212 213 int block_size; 213 214 #endif … … 235 236 #endif 236 237 237 // Root node of the starting spanning tree.238 Node root;239 // The entering edge of the current pivot iteration.240 Edge in_edge;241 // Temporary nodes used in the current pivot iteration.242 Node join, u_in, v_in, u_out, v_out;243 Node right, first, second, last;244 Node stem, par_stem, new_stem;245 // The maximum augment amount along the cycle in the current pivot246 // iteration.247 Capacity delta;248 249 238 public : 250 239 … … 259 248 /// \param _supply The supply values of the nodes (signed). 260 249 NetworkSimplex( const Graph &_graph, 261 262 263 264 250 const LowerMap &_lower, 251 const CapacityMap &_capacity, 252 const CostMap &_cost, 253 const SupplyMap &_supply ) : 265 254 graph_ref(_graph), lower(&_lower), capacity(graph), cost(graph), 266 255 supply(graph), flow(graph), flow_result(_graph), potential(graph), … … 273 262 Supply sum = 0; 274 263 for (typename Graph::NodeIt n(graph_ref); n != INVALID; ++n) 275 264 sum += _supply[n]; 276 265 if (!(valid_supply = sum == 0)) return; 277 266 … … 280 269 graph.reserveEdge(countEdges(graph_ref) + countNodes(graph_ref)); 281 270 copyGraph(graph, graph_ref) 282 283 284 285 286 287 // Removing non zero lower bounds271 .edgeMap(cost, _cost) 272 .nodeRef(node_ref) 273 .edgeRef(edge_ref) 274 .run(); 275 276 // Removing non-zero lower bounds 288 277 for (typename Graph::EdgeIt e(graph_ref); e != INVALID; ++e) { 289 278 capacity[edge_ref[e]] = _capacity[e] - _lower[e]; 290 279 } 291 280 for (typename Graph::NodeIt n(graph_ref); n != INVALID; ++n) { 292 293 294 295 296 297 281 Supply s = _supply[n]; 282 for (typename Graph::InEdgeIt e(graph_ref, n); e != INVALID; ++e) 283 s += _lower[e]; 284 for (typename Graph::OutEdgeIt e(graph_ref, n); e != INVALID; ++e) 285 s -= _lower[e]; 286 supply[node_ref[n]] = s; 298 287 } 299 288 } … … 308 297 /// \param _supply The supply values of the nodes (signed). 309 298 NetworkSimplex( const Graph &_graph, 310 311 312 299 const CapacityMap &_capacity, 300 const CostMap &_cost, 301 const SupplyMap &_supply ) : 313 302 graph_ref(_graph), lower(NULL), capacity(graph), cost(graph), 314 303 supply(graph), flow(graph), flow_result(_graph), potential(graph), … … 321 310 Supply sum = 0; 322 311 for (typename Graph::NodeIt n(graph_ref); n != INVALID; ++n) 323 312 sum += _supply[n]; 324 313 if (!(valid_supply = sum == 0)) return; 325 314 326 315 // Copying graph_ref to graph 327 316 copyGraph(graph, graph_ref) 328 329 330 331 332 333 317 .edgeMap(capacity, _capacity) 318 .edgeMap(cost, _cost) 319 .nodeMap(supply, _supply) 320 .nodeRef(node_ref) 321 .edgeRef(edge_ref) 322 .run(); 334 323 } 335 324 … … 345 334 /// \param _t The target node. 346 335 /// \param _flow_value The required amount of flow from node \c _s 347 /// to node \c _t (i.e. the supply of \c _s and the demand of 348 /// \c _t). 336 /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t). 349 337 NetworkSimplex( const Graph &_graph, 350 351 352 353 354 355 338 const LowerMap &_lower, 339 const CapacityMap &_capacity, 340 const CostMap &_cost, 341 typename Graph::Node _s, 342 typename Graph::Node _t, 343 typename SupplyMap::Value _flow_value ) : 356 344 graph_ref(_graph), lower(&_lower), capacity(graph), cost(graph), 357 345 supply(graph), flow(graph), flow_result(_graph), potential(graph), … … 363 351 // Copying graph_ref to graph 364 352 copyGraph(graph, graph_ref) 365 366 367 368 369 370 // Removing non zero lower bounds353 .edgeMap(cost, _cost) 354 .nodeRef(node_ref) 355 .edgeRef(edge_ref) 356 .run(); 357 358 // Removing non-zero lower bounds 371 359 for (typename Graph::EdgeIt e(graph_ref); e != INVALID; ++e) { 372 360 capacity[edge_ref[e]] = _capacity[e] - _lower[e]; 373 361 } 374 362 for (typename Graph::NodeIt n(graph_ref); n != INVALID; ++n) { 375 376 377 378 379 380 381 382 363 Supply s = 0; 364 if (n == _s) s = _flow_value; 365 if (n == _t) s = -_flow_value; 366 for (typename Graph::InEdgeIt e(graph_ref, n); e != INVALID; ++e) 367 s += _lower[e]; 368 for (typename Graph::OutEdgeIt e(graph_ref, n); e != INVALID; ++e) 369 s -= _lower[e]; 370 supply[node_ref[n]] = s; 383 371 } 384 372 valid_supply = true; … … 395 383 /// \param _t The target node. 396 384 /// \param _flow_value The required amount of flow from node \c _s 397 /// to node \c _t (i.e. the supply of \c _s and the demand of 398 /// \c _t). 385 /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t). 399 386 NetworkSimplex( const Graph &_graph, 400 401 402 403 404 387 const CapacityMap &_capacity, 388 const CostMap &_cost, 389 typename Graph::Node _s, 390 typename Graph::Node _t, 391 typename SupplyMap::Value _flow_value ) : 405 392 graph_ref(_graph), lower(NULL), capacity(graph), cost(graph), 406 393 supply(graph, 0), flow(graph), flow_result(_graph), potential(graph), … … 412 399 // Copying graph_ref to graph 413 400 copyGraph(graph, graph_ref) 414 415 416 417 418 401 .edgeMap(capacity, _capacity) 402 .edgeMap(cost, _cost) 403 .nodeRef(node_ref) 404 .edgeRef(edge_ref) 405 .run(); 419 406 supply[node_ref[_s]] = _flow_value; 420 407 supply[node_ref[_t]] = -_flow_value; 421 408 valid_supply = true; 409 } 410 411 /// \brief Runs the algorithm. 412 /// 413 /// Runs the algorithm. 414 /// 415 /// \return \c true if a feasible flow can be found. 416 bool run() { 417 return init() && start(); 422 418 } 423 419 … … 451 447 Cost c = 0; 452 448 for (typename Graph::EdgeIt e(graph_ref); e != INVALID; ++e) 453 449 c += flow_result[e] * cost[edge_ref[e]]; 454 450 return c; 455 }456 457 /// \brief Runs the algorithm.458 ///459 /// Runs the algorithm.460 ///461 /// \return \c true if a feasible flow can be found.462 bool run() {463 return init() && start();464 451 } 465 452 … … 473 460 // Initializing state and flow maps 474 461 for (EdgeIt e(graph); e != INVALID; ++e) { 475 476 462 flow[e] = 0; 463 state[e] = LOWER; 477 464 } 478 465 … … 492 479 Cost max_cost = std::numeric_limits<Cost>::max() / 4; 493 480 for (NodeIt u(graph); u != INVALID; ++u) { 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 481 if (u == root) continue; 482 thread[last] = u; 483 last = u; 484 parent[u] = root; 485 depth[u] = 1; 486 sum += supply[u]; 487 if (supply[u] >= 0) { 488 e = graph.addEdge(u, root); 489 flow[e] = supply[u]; 490 forward[u] = true; 491 potential[u] = max_cost; 492 } else { 493 e = graph.addEdge(root, u); 494 flow[e] = -supply[u]; 495 forward[u] = false; 496 potential[u] = -max_cost; 497 } 498 cost[e] = max_cost; 499 capacity[e] = std::numeric_limits<Capacity>::max(); 500 state[e] = TREE; 501 pred_edge[u] = e; 515 502 } 516 503 thread[last] = root; … … 521 508 block_size = 2 * int(sqrt(countEdges(graph))); 522 509 if (block_size < MIN_BLOCK_SIZE) block_size = MIN_BLOCK_SIZE; 523 // block_size = edge_num >= BLOCK_NUM * MIN_BLOCK_SIZE ?524 // edge_num / BLOCK_NUM : MIN_BLOCK_SIZE;525 510 #endif 526 511 #ifdef CANDIDATE_LIST_PIVOT … … 544 529 bool findEnteringEdge(EdgeIt &next_edge) { 545 530 for (EdgeIt e = next_edge; e != INVALID; ++e) { 546 547 548 549 550 531 if (state[e] * red_cost[e] < 0) { 532 in_edge = e; 533 next_edge = ++e; 534 return true; 535 } 551 536 } 552 537 for (EdgeIt e(graph); e != next_edge; ++e) { 553 554 555 556 557 538 if (state[e] * red_cost[e] < 0) { 539 in_edge = e; 540 next_edge = ++e; 541 return true; 542 } 558 543 } 559 544 return false; … … 567 552 Cost min = 0; 568 553 for (EdgeIt e(graph); e != INVALID; ++e) { 569 570 571 572 554 if (state[e] * red_cost[e] < min) { 555 min = state[e] * red_cost[e]; 556 in_edge = e; 557 } 573 558 } 574 559 return min < 0; … … 585 570 int cnt = 0; 586 571 for (EdgeIt e = next_edge; e != INVALID; ++e) { 587 588 589 590 591 592 593 594 572 if ((curr = state[e] * red_cost[e]) < min) { 573 min = curr; 574 min_edge = e; 575 } 576 if (++cnt == block_size) { 577 if (min < 0) break; 578 cnt = 0; 579 } 595 580 } 596 581 if (!(min < 0)) { 597 598 599 600 601 602 603 604 605 606 582 for (EdgeIt e(graph); e != next_edge; ++e) { 583 if ((curr = state[e] * red_cost[e]) < min) { 584 min = curr; 585 min_edge = e; 586 } 587 if (++cnt == block_size) { 588 if (min < 0) break; 589 cnt = 0; 590 } 591 } 607 592 } 608 593 in_edge = min_edge; 609 594 if ((next_edge = ++min_edge) == INVALID) 610 595 next_edge = EdgeIt(graph); 611 596 return min < 0; 612 597 } … … 620 605 621 606 if (minor_count >= minor_limit || candidates.size() == 0) { 622 623 624 625 626 627 628 629 630 607 // Major iteration 608 candidates.clear(); 609 for (EdgeIt e(graph); e != INVALID; ++e) { 610 if (state[e] * red_cost[e] < 0) { 611 candidates.push_back(e); 612 if (candidates.size() == list_length) break; 613 } 614 } 615 if (candidates.size() == 0) return false; 631 616 } 632 617 … … 637 622 for (int i = 0; i < candidates.size(); ++i) { 638 623 e = candidates[i]; 639 640 641 642 624 if (state[e] * red_cost[e] < min) { 625 min = state[e] * red_cost[e]; 626 in_edge = e; 627 } 643 628 } 644 629 return true; … … 656 641 public: 657 642 SortFunc(const IntEdgeMap &_st, const ReducedCostMap &_rc) : 658 643 st(_st), rc(_rc) {} 659 644 bool operator()(const Edge &e1, const Edge &e2) { 660 645 return st[e1] * rc[e1] < st[e2] * rc[e2]; 661 646 } 662 647 }; … … 669 654 // Minor iteration 670 655 while (list_index < candidates.size()) { 671 672 656 in_edge = candidates[list_index++]; 657 if (state[in_edge] * red_cost[in_edge] < 0) return true; 673 658 } 674 659 … … 677 662 Cost curr, min = 0; 678 663 for (EdgeIt e(graph); e != INVALID; ++e) { 679 680 681 682 683 664 if ((curr = state[e] * red_cost[e]) < min / LOWER_DIV) { 665 candidates.push_back(e); 666 if (curr < min) min = curr; 667 if (candidates.size() == list_length) break; 668 } 684 669 } 685 670 if (candidates.size() == 0) return false; … … 697 682 Node v = graph.target(in_edge); 698 683 while (u != v) { 699 700 701 702 703 704 684 if (depth[u] == depth[v]) { 685 u = parent[u]; 686 v = parent[v]; 687 } 688 else if (depth[u] > depth[v]) u = parent[u]; 689 else v = parent[v]; 705 690 } 706 691 return u; … … 713 698 // of the cycle 714 699 if (state[in_edge] == LOWER) { 715 716 second= graph.target(in_edge);700 first = graph.source(in_edge); 701 second = graph.target(in_edge); 717 702 } else { 718 719 second= graph.source(in_edge);703 first = graph.target(in_edge); 704 second = graph.source(in_edge); 720 705 } 721 706 delta = capacity[in_edge]; … … 727 712 // root node 728 713 for (Node u = first; u != join; u = parent[u]) { 729 730 731 732 733 734 735 736 737 714 e = pred_edge[u]; 715 d = forward[u] ? flow[e] : capacity[e] - flow[e]; 716 if (d < delta) { 717 delta = d; 718 u_out = u; 719 u_in = first; 720 v_in = second; 721 result = true; 722 } 738 723 } 739 724 // Searching the cycle along the path form the second node to the 740 725 // root node 741 726 for (Node u = second; u != join; u = parent[u]) { 742 743 744 745 746 747 748 749 750 727 e = pred_edge[u]; 728 d = forward[u] ? capacity[e] - flow[e] : flow[e]; 729 if (d <= delta) { 730 delta = d; 731 u_out = u; 732 u_in = second; 733 v_in = first; 734 result = true; 735 } 751 736 } 752 737 return result; 753 738 } 754 739 755 /// \brief Changes flow andstate edge maps.740 /// \brief Changes \c flow and \c state edge maps. 756 741 void changeFlows(bool change) { 757 742 // Augmenting along the cycle 758 743 if (delta > 0) { 759 760 761 762 763 764 765 766 744 Capacity val = state[in_edge] * delta; 745 flow[in_edge] += val; 746 for (Node u = graph.source(in_edge); u != join; u = parent[u]) { 747 flow[pred_edge[u]] += forward[u] ? -val : val; 748 } 749 for (Node u = graph.target(in_edge); u != join; u = parent[u]) { 750 flow[pred_edge[u]] += forward[u] ? val : -val; 751 } 767 752 } 768 753 // Updating the state of the entering and leaving edges 769 754 if (change) { 770 771 772 755 state[in_edge] = TREE; 756 state[pred_edge[u_out]] = 757 (flow[pred_edge[u_out]] == 0) ? LOWER : UPPER; 773 758 } else { 774 775 } 776 } 777 778 /// \brief Updates thread andparent node maps.759 state[in_edge] = -state[in_edge]; 760 } 761 } 762 763 /// \brief Updates \c thread and \c parent node maps. 779 764 void updateThreadParent() { 780 765 Node u; … … 784 769 bool par_first = false; 785 770 if (join == v_out) { 786 787 788 789 790 791 771 for (u = join; u != u_in && u != v_in; u = thread[u]) ; 772 if (u == v_in) { 773 par_first = true; 774 while (thread[u] != u_out) u = thread[u]; 775 first = u; 776 } 792 777 } 793 778 … … 797 782 right = thread[u]; 798 783 if (thread[v_in] == u_out) { 799 800 784 for (last = u; depth[last] > depth[u_out]; last = thread[last]) ; 785 if (last == u_out) last = thread[last]; 801 786 } 802 787 else last = thread[v_in]; … … 806 791 par_stem = v_in; 807 792 while (stem != u_out) { 808 809 810 811 812 813 814 815 816 817 818 819 820 821 822 823 824 793 thread[u] = new_stem = parent[stem]; 794 795 // Finding the node just before the stem node (u) according to 796 // the original thread index 797 for (u = new_stem; thread[u] != stem; u = thread[u]) ; 798 thread[u] = right; 799 800 // Changing the parent node of stem and shifting stem and 801 // par_stem nodes 802 parent[stem] = par_stem; 803 par_stem = stem; 804 stem = new_stem; 805 806 // Finding the last successor of stem (u) and the node after it 807 // (right) according to the thread index 808 for (u = stem; depth[thread[u]] > depth[stem]; u = thread[u]) ; 809 right = thread[u]; 825 810 } 826 811 parent[u_out] = par_stem; … … 828 813 829 814 if (join == v_out && par_first) { 830 815 if (first != v_in) thread[first] = right; 831 816 } else { 832 833 834 } 835 } 836 837 /// \brief Updates pred_edge andforward node maps.817 for (u = v_out; thread[u] != u_out; u = thread[u]) ; 818 thread[u] = right; 819 } 820 } 821 822 /// \brief Updates \c pred_edge and \c forward node maps. 838 823 void updatePredEdge() { 839 824 Node u = u_out, v; 840 825 while (u != u_in) { 841 842 843 844 826 v = parent[u]; 827 pred_edge[u] = pred_edge[v]; 828 forward[u] = !forward[v]; 829 u = v; 845 830 } 846 831 pred_edge[u_in] = in_edge; … … 848 833 } 849 834 850 /// \brief Updates depth andpotential node maps.835 /// \brief Updates \c depth and \c potential node maps. 851 836 void updateDepthPotential() { 852 837 depth[u_in] = depth[v_in] + 1; 853 838 potential[u_in] = forward[u_in] ? 854 855 839 potential[v_in] + cost[pred_edge[u_in]] : 840 potential[v_in] - cost[pred_edge[u_in]]; 856 841 857 842 Node u = thread[u_in], v; 858 843 while (true) { 859 860 861 862 863 864 865 866 844 v = parent[u]; 845 if (v == INVALID) break; 846 depth[u] = depth[v] + 1; 847 potential[u] = forward[u] ? 848 potential[v] + cost[pred_edge[u]] : 849 potential[v] - cost[pred_edge[u]]; 850 if (depth[u] <= depth[v_in]) break; 851 u = thread[u]; 867 852 } 868 853 } … … 881 866 #endif 882 867 { 883 884 885 886 887 888 889 890 868 join = findJoinNode(); 869 bool change = findLeavingEdge(); 870 changeFlows(change); 871 if (change) { 872 updateThreadParent(); 873 updatePredEdge(); 874 updateDepthPotential(); 875 } 891 876 #ifdef _DEBUG_ITER_ 892 877 ++iter_num; 893 878 #endif 894 879 } … … 896 881 #ifdef _DEBUG_ITER_ 897 882 std::cout << "Network Simplex algorithm finished. " << iter_num 898 883 << " pivot iterations performed." << std::endl; 899 884 #endif 900 885 … … 902 887 // artificial edges 903 888 for (InEdgeIt e(graph, root); e != INVALID; ++e) 904 889 if (flow[e] > 0) return false; 905 890 for (OutEdgeIt e(graph, root); e != INVALID; ++e) 906 891 if (flow[e] > 0) return false; 907 892 908 893 // Copying flow values to flow_result 909 894 if (lower) { 910 911 895 for (typename Graph::EdgeIt e(graph_ref); e != INVALID; ++e) 896 flow_result[e] = (*lower)[e] + flow[edge_ref[e]]; 912 897 } else { 913 914 898 for (typename Graph::EdgeIt e(graph_ref); e != INVALID; ++e) 899 flow_result[e] = flow[edge_ref[e]]; 915 900 } 916 901 // Copying potential values to potential_result 917 902 for (typename Graph::NodeIt n(graph_ref); n != INVALID; ++n) 918 903 potential_result[n] = potential[node_ref[n]]; 919 904 920 905 return true;
Note: See TracChangeset
for help on using the changeset viewer.