48 /// \tparam SupplyMap The type of the supply map. |
48 /// \tparam SupplyMap The type of the supply map. |
49 /// |
49 /// |
50 /// \warning |
50 /// \warning |
51 /// - Edge capacities and costs should be \e non-negative \e integers. |
51 /// - Edge capacities and costs should be \e non-negative \e integers. |
52 /// - Supply values should be \e signed \e integers. |
52 /// - Supply values should be \e signed \e integers. |
53 /// - \c LowerMap::Value must be convertible to \c CapacityMap::Value. |
53 /// - The value types of the maps should be convertible to each other. |
54 /// - \c CapacityMap::Value and \c SupplyMap::Value must be |
54 /// - \c CostMap::Value must be signed type. |
55 /// convertible to each other. |
|
56 /// - All value types must be convertible to \c CostMap::Value, which |
|
57 /// must be signed type. |
|
58 /// |
55 /// |
59 /// \author Peter Kovacs |
56 /// \author Peter Kovacs |
60 |
57 |
61 template < typename Graph, |
58 template < typename Graph, |
62 typename LowerMap = typename Graph::template EdgeMap<int>, |
59 typename LowerMap = typename Graph::template EdgeMap<int>, |
118 // The processed (i.e. permanently labeled) nodes |
115 // The processed (i.e. permanently labeled) nodes |
119 std::vector<Node> _proc_nodes; |
116 std::vector<Node> _proc_nodes; |
120 |
117 |
121 public: |
118 public: |
122 |
119 |
123 /// The constructor of the class. |
120 /// Constructor. |
124 ResidualDijkstra( const Graph &graph, |
121 ResidualDijkstra( const Graph &graph, |
125 const FlowMap &flow, |
122 const FlowMap &flow, |
126 const CapacityEdgeMap &res_cap, |
123 const CapacityEdgeMap &res_cap, |
127 const CostMap &cost, |
124 const CostMap &cost, |
128 const SupplyMap &excess, |
125 const SupplyMap &excess, |
219 // The modified supply map |
216 // The modified supply map |
220 SupplyNodeMap _supply; |
217 SupplyNodeMap _supply; |
221 bool _valid_supply; |
218 bool _valid_supply; |
222 |
219 |
223 // Edge map of the current flow |
220 // Edge map of the current flow |
224 FlowMap _flow; |
221 FlowMap *_flow; |
|
222 bool _local_flow; |
225 // Node map of the current potentials |
223 // Node map of the current potentials |
226 PotentialMap _potential; |
224 PotentialMap *_potential; |
|
225 bool _local_potential; |
227 |
226 |
228 // The residual capacity map |
227 // The residual capacity map |
229 CapacityEdgeMap _res_cap; |
228 CapacityEdgeMap _res_cap; |
230 // The excess map |
229 // The excess map |
231 SupplyNodeMap _excess; |
230 SupplyNodeMap _excess; |
241 |
240 |
242 // The pred edge map |
241 // The pred edge map |
243 PredMap _pred; |
242 PredMap _pred; |
244 // Implementation of the Dijkstra algorithm for finding augmenting |
243 // Implementation of the Dijkstra algorithm for finding augmenting |
245 // shortest paths in the residual network |
244 // shortest paths in the residual network |
246 ResidualDijkstra _dijkstra; |
245 ResidualDijkstra *_dijkstra; |
247 |
246 |
248 public : |
247 public: |
249 |
248 |
250 /// \brief General constructor of the class (with lower bounds). |
249 /// \brief General constructor (with lower bounds). |
251 /// |
250 /// |
252 /// General constructor of the class (with lower bounds). |
251 /// General constructor (with lower bounds). |
253 /// |
252 /// |
254 /// \param graph The directed graph the algorithm runs on. |
253 /// \param graph The directed graph the algorithm runs on. |
255 /// \param lower The lower bounds of the edges. |
254 /// \param lower The lower bounds of the edges. |
256 /// \param capacity The capacities (upper bounds) of the edges. |
255 /// \param capacity The capacities (upper bounds) of the edges. |
257 /// \param cost The cost (length) values of the edges. |
256 /// \param cost The cost (length) values of the edges. |
260 const LowerMap &lower, |
259 const LowerMap &lower, |
261 const CapacityMap &capacity, |
260 const CapacityMap &capacity, |
262 const CostMap &cost, |
261 const CostMap &cost, |
263 const SupplyMap &supply ) : |
262 const SupplyMap &supply ) : |
264 _graph(graph), _lower(&lower), _capacity(graph), _cost(cost), |
263 _graph(graph), _lower(&lower), _capacity(graph), _cost(cost), |
265 _supply(graph), _flow(graph, 0), _potential(graph, 0), |
264 _supply(graph), _flow(0), _local_flow(false), |
266 _res_cap(graph), _excess(graph), _pred(graph), |
265 _potential(0), _local_potential(false), |
267 _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred) |
266 _res_cap(graph), _excess(graph), _pred(graph) |
268 { |
267 { |
269 // Removing non-zero lower bounds |
268 // Removing non-zero lower bounds |
270 _capacity = subMap(capacity, lower); |
269 _capacity = subMap(capacity, lower); |
271 _res_cap = _capacity; |
270 _res_cap = _capacity; |
272 Supply sum = 0; |
271 Supply sum = 0; |
280 sum += s; |
279 sum += s; |
281 } |
280 } |
282 _valid_supply = sum == 0; |
281 _valid_supply = sum == 0; |
283 } |
282 } |
284 |
283 |
285 /// \brief General constructor of the class (without lower bounds). |
284 /// \brief General constructor (without lower bounds). |
286 /// |
285 /// |
287 /// General constructor of the class (without lower bounds). |
286 /// General constructor (without lower bounds). |
288 /// |
287 /// |
289 /// \param graph The directed graph the algorithm runs on. |
288 /// \param graph The directed graph the algorithm runs on. |
290 /// \param capacity The capacities (upper bounds) of the edges. |
289 /// \param capacity The capacities (upper bounds) of the edges. |
291 /// \param cost The cost (length) values of the edges. |
290 /// \param cost The cost (length) values of the edges. |
292 /// \param supply The supply values of the nodes (signed). |
291 /// \param supply The supply values of the nodes (signed). |
293 CapacityScaling( const Graph &graph, |
292 CapacityScaling( const Graph &graph, |
294 const CapacityMap &capacity, |
293 const CapacityMap &capacity, |
295 const CostMap &cost, |
294 const CostMap &cost, |
296 const SupplyMap &supply ) : |
295 const SupplyMap &supply ) : |
297 _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost), |
296 _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost), |
298 _supply(supply), _flow(graph, 0), _potential(graph, 0), |
297 _supply(supply), _flow(0), _local_flow(false), |
299 _res_cap(capacity), _excess(graph), _pred(graph), |
298 _potential(0), _local_potential(false), |
300 _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred) |
299 _res_cap(capacity), _excess(graph), _pred(graph) |
301 { |
300 { |
302 // Checking the sum of supply values |
301 // Checking the sum of supply values |
303 Supply sum = 0; |
302 Supply sum = 0; |
304 for (NodeIt n(_graph); n != INVALID; ++n) sum += _supply[n]; |
303 for (NodeIt n(_graph); n != INVALID; ++n) sum += _supply[n]; |
305 _valid_supply = sum == 0; |
304 _valid_supply = sum == 0; |
306 } |
305 } |
307 |
306 |
308 /// \brief Simple constructor of the class (with lower bounds). |
307 /// \brief Simple constructor (with lower bounds). |
309 /// |
308 /// |
310 /// Simple constructor of the class (with lower bounds). |
309 /// Simple constructor (with lower bounds). |
311 /// |
310 /// |
312 /// \param graph The directed graph the algorithm runs on. |
311 /// \param graph The directed graph the algorithm runs on. |
313 /// \param lower The lower bounds of the edges. |
312 /// \param lower The lower bounds of the edges. |
314 /// \param capacity The capacities (upper bounds) of the edges. |
313 /// \param capacity The capacities (upper bounds) of the edges. |
315 /// \param cost The cost (length) values of the edges. |
314 /// \param cost The cost (length) values of the edges. |
322 const CapacityMap &capacity, |
321 const CapacityMap &capacity, |
323 const CostMap &cost, |
322 const CostMap &cost, |
324 Node s, Node t, |
323 Node s, Node t, |
325 Supply flow_value ) : |
324 Supply flow_value ) : |
326 _graph(graph), _lower(&lower), _capacity(graph), _cost(cost), |
325 _graph(graph), _lower(&lower), _capacity(graph), _cost(cost), |
327 _supply(graph), _flow(graph, 0), _potential(graph, 0), |
326 _supply(graph), _flow(0), _local_flow(false), |
328 _res_cap(graph), _excess(graph), _pred(graph), |
327 _potential(0), _local_potential(false), |
329 _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred) |
328 _res_cap(graph), _excess(graph), _pred(graph) |
330 { |
329 { |
331 // Removing non-zero lower bounds |
330 // Removing non-zero lower bounds |
332 _capacity = subMap(capacity, lower); |
331 _capacity = subMap(capacity, lower); |
333 _res_cap = _capacity; |
332 _res_cap = _capacity; |
334 for (NodeIt n(_graph); n != INVALID; ++n) { |
333 for (NodeIt n(_graph); n != INVALID; ++n) { |
342 _supply[n] = sum; |
341 _supply[n] = sum; |
343 } |
342 } |
344 _valid_supply = true; |
343 _valid_supply = true; |
345 } |
344 } |
346 |
345 |
347 /// \brief Simple constructor of the class (without lower bounds). |
346 /// \brief Simple constructor (without lower bounds). |
348 /// |
347 /// |
349 /// Simple constructor of the class (without lower bounds). |
348 /// Simple constructor (without lower bounds). |
350 /// |
349 /// |
351 /// \param graph The directed graph the algorithm runs on. |
350 /// \param graph The directed graph the algorithm runs on. |
352 /// \param capacity The capacities (upper bounds) of the edges. |
351 /// \param capacity The capacities (upper bounds) of the edges. |
353 /// \param cost The cost (length) values of the edges. |
352 /// \param cost The cost (length) values of the edges. |
354 /// \param s The source node. |
353 /// \param s The source node. |
359 const CapacityMap &capacity, |
358 const CapacityMap &capacity, |
360 const CostMap &cost, |
359 const CostMap &cost, |
361 Node s, Node t, |
360 Node s, Node t, |
362 Supply flow_value ) : |
361 Supply flow_value ) : |
363 _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost), |
362 _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost), |
364 _supply(graph, 0), _flow(graph, 0), _potential(graph, 0), |
363 _supply(graph, 0), _flow(0), _local_flow(false), |
365 _res_cap(capacity), _excess(graph), _pred(graph), |
364 _potential(0), _local_potential(false), |
366 _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred) |
365 _res_cap(capacity), _excess(graph), _pred(graph) |
367 { |
366 { |
368 _supply[s] = flow_value; |
367 _supply[s] = flow_value; |
369 _supply[t] = -flow_value; |
368 _supply[t] = -flow_value; |
370 _valid_supply = true; |
369 _valid_supply = true; |
371 } |
370 } |
372 |
371 |
|
372 /// Destructor. |
|
373 ~CapacityScaling() { |
|
374 if (_local_flow) delete _flow; |
|
375 if (_local_potential) delete _potential; |
|
376 delete _dijkstra; |
|
377 } |
|
378 |
|
379 /// \brief Sets the flow map. |
|
380 /// |
|
381 /// Sets the flow map. |
|
382 /// |
|
383 /// \return \c (*this) |
|
384 CapacityScaling& flowMap(FlowMap &map) { |
|
385 if (_local_flow) { |
|
386 delete _flow; |
|
387 _local_flow = false; |
|
388 } |
|
389 _flow = ↦ |
|
390 return *this; |
|
391 } |
|
392 |
|
393 /// \brief Sets the potential map. |
|
394 /// |
|
395 /// Sets the potential map. |
|
396 /// |
|
397 /// \return \c (*this) |
|
398 CapacityScaling& potentialMap(PotentialMap &map) { |
|
399 if (_local_potential) { |
|
400 delete _potential; |
|
401 _local_potential = false; |
|
402 } |
|
403 _potential = ↦ |
|
404 return *this; |
|
405 } |
|
406 |
|
407 /// \name Execution control |
|
408 /// The only way to execute the algorithm is to call the run() |
|
409 /// function. |
|
410 |
|
411 /// @{ |
|
412 |
373 /// \brief Runs the algorithm. |
413 /// \brief Runs the algorithm. |
374 /// |
414 /// |
375 /// Runs the algorithm. |
415 /// Runs the algorithm. |
376 /// |
416 /// |
377 /// \param scaling Enable or disable capacity scaling. |
417 /// \param scaling Enable or disable capacity scaling. |
382 /// \return \c true if a feasible flow can be found. |
422 /// \return \c true if a feasible flow can be found. |
383 bool run(bool scaling = true) { |
423 bool run(bool scaling = true) { |
384 return init(scaling) && start(); |
424 return init(scaling) && start(); |
385 } |
425 } |
386 |
426 |
|
427 /// @} |
|
428 |
|
429 /// \name Query Functions |
|
430 /// The result of the algorithm can be obtained using these |
|
431 /// functions. |
|
432 /// \n run() must be called before using them. |
|
433 |
|
434 /// @{ |
|
435 |
387 /// \brief Returns a const reference to the edge map storing the |
436 /// \brief Returns a const reference to the edge map storing the |
388 /// found flow. |
437 /// found flow. |
389 /// |
438 /// |
390 /// Returns a const reference to the edge map storing the found flow. |
439 /// Returns a const reference to the edge map storing the found flow. |
391 /// |
440 /// |
392 /// \pre \ref run() must be called before using this function. |
441 /// \pre \ref run() must be called before using this function. |
393 const FlowMap& flowMap() const { |
442 const FlowMap& flowMap() const { |
394 return _flow; |
443 return *_flow; |
395 } |
444 } |
396 |
445 |
397 /// \brief Returns a const reference to the node map storing the |
446 /// \brief Returns a const reference to the node map storing the |
398 /// found potentials (the dual solution). |
447 /// found potentials (the dual solution). |
399 /// |
448 /// |
400 /// Returns a const reference to the node map storing the found |
449 /// Returns a const reference to the node map storing the found |
401 /// potentials (the dual solution). |
450 /// potentials (the dual solution). |
402 /// |
451 /// |
403 /// \pre \ref run() must be called before using this function. |
452 /// \pre \ref run() must be called before using this function. |
404 const PotentialMap& potentialMap() const { |
453 const PotentialMap& potentialMap() const { |
405 return _potential; |
454 return *_potential; |
|
455 } |
|
456 |
|
457 /// \brief Returns the flow on the edge. |
|
458 /// |
|
459 /// Returns the flow on the edge. |
|
460 /// |
|
461 /// \pre \ref run() must be called before using this function. |
|
462 Capacity flow(const Edge& edge) const { |
|
463 return (*_flow)[edge]; |
|
464 } |
|
465 |
|
466 /// \brief Returns the potential of the node. |
|
467 /// |
|
468 /// Returns the potential of the node. |
|
469 /// |
|
470 /// \pre \ref run() must be called before using this function. |
|
471 Cost potential(const Node& node) const { |
|
472 return (*_potential)[node]; |
406 } |
473 } |
407 |
474 |
408 /// \brief Returns the total cost of the found flow. |
475 /// \brief Returns the total cost of the found flow. |
409 /// |
476 /// |
410 /// Returns the total cost of the found flow. The complexity of the |
477 /// Returns the total cost of the found flow. The complexity of the |
412 /// |
479 /// |
413 /// \pre \ref run() must be called before using this function. |
480 /// \pre \ref run() must be called before using this function. |
414 Cost totalCost() const { |
481 Cost totalCost() const { |
415 Cost c = 0; |
482 Cost c = 0; |
416 for (EdgeIt e(_graph); e != INVALID; ++e) |
483 for (EdgeIt e(_graph); e != INVALID; ++e) |
417 c += _flow[e] * _cost[e]; |
484 c += (*_flow)[e] * _cost[e]; |
418 return c; |
485 return c; |
419 } |
486 } |
|
487 |
|
488 /// @} |
420 |
489 |
421 private: |
490 private: |
422 |
491 |
423 /// Initializes the algorithm. |
492 /// Initializes the algorithm. |
424 bool init(bool scaling) { |
493 bool init(bool scaling) { |
425 if (!_valid_supply) return false; |
494 if (!_valid_supply) return false; |
|
495 |
|
496 // Initializing maps |
|
497 if (!_flow) { |
|
498 _flow = new FlowMap(_graph); |
|
499 _local_flow = true; |
|
500 } |
|
501 if (!_potential) { |
|
502 _potential = new PotentialMap(_graph); |
|
503 _local_potential = true; |
|
504 } |
|
505 for (EdgeIt e(_graph); e != INVALID; ++e) (*_flow)[e] = 0; |
|
506 for (NodeIt n(_graph); n != INVALID; ++n) (*_potential)[n] = 0; |
426 _excess = _supply; |
507 _excess = _supply; |
427 |
508 |
428 // Initilaizing delta value |
509 _dijkstra = new ResidualDijkstra( _graph, *_flow, _res_cap, _cost, |
|
510 _excess, *_potential, _pred ); |
|
511 |
|
512 // Initializing delta value |
429 if (scaling) { |
513 if (scaling) { |
430 // With scaling |
514 // With scaling |
431 Supply max_sup = 0, max_dem = 0; |
515 Supply max_sup = 0, max_dem = 0; |
432 for (NodeIt n(_graph); n != INVALID; ++n) { |
516 for (NodeIt n(_graph); n != INVALID; ++n) { |
433 if ( _supply[n] > max_sup) max_sup = _supply[n]; |
517 if ( _supply[n] > max_sup) max_sup = _supply[n]; |
459 int factor = 4; |
544 int factor = 4; |
460 while (true) { |
545 while (true) { |
461 // Saturating all edges not satisfying the optimality condition |
546 // Saturating all edges not satisfying the optimality condition |
462 for (EdgeIt e(_graph); e != INVALID; ++e) { |
547 for (EdgeIt e(_graph); e != INVALID; ++e) { |
463 Node u = _graph.source(e), v = _graph.target(e); |
548 Node u = _graph.source(e), v = _graph.target(e); |
464 Cost c = _cost[e] + _potential[u] - _potential[v]; |
549 Cost c = _cost[e] + (*_potential)[u] - (*_potential)[v]; |
465 if (c < 0 && _res_cap[e] >= _delta) { |
550 if (c < 0 && _res_cap[e] >= _delta) { |
466 _excess[u] -= _res_cap[e]; |
551 _excess[u] -= _res_cap[e]; |
467 _excess[v] += _res_cap[e]; |
552 _excess[v] += _res_cap[e]; |
468 _flow[e] = _capacity[e]; |
553 (*_flow)[e] = _capacity[e]; |
469 _res_cap[e] = 0; |
554 _res_cap[e] = 0; |
470 } |
555 } |
471 else if (c > 0 && _flow[e] >= _delta) { |
556 else if (c > 0 && (*_flow)[e] >= _delta) { |
472 _excess[u] += _flow[e]; |
557 _excess[u] += (*_flow)[e]; |
473 _excess[v] -= _flow[e]; |
558 _excess[v] -= (*_flow)[e]; |
474 _flow[e] = 0; |
559 (*_flow)[e] = 0; |
475 _res_cap[e] = _capacity[e]; |
560 _res_cap[e] = _capacity[e]; |
476 } |
561 } |
477 } |
562 } |
478 |
563 |
479 // Finding excess nodes and deficit nodes |
564 // Finding excess nodes and deficit nodes |
570 while ( _excess[_excess_nodes[next_node]] > 0 || |
655 while ( _excess[_excess_nodes[next_node]] > 0 || |
571 ++next_node < int(_excess_nodes.size()) ) |
656 ++next_node < int(_excess_nodes.size()) ) |
572 { |
657 { |
573 // Running Dijkstra |
658 // Running Dijkstra |
574 s = _excess_nodes[next_node]; |
659 s = _excess_nodes[next_node]; |
575 if ((t = _dijkstra.run(s, 1)) == INVALID) |
660 if ((t = _dijkstra->run(s, 1)) == INVALID) |
576 return false; |
661 return false; |
577 |
662 |
578 // Augmenting along a shortest path from s to t |
663 // Augmenting along a shortest path from s to t |
579 Capacity d = _excess[s] < -_excess[t] ? _excess[s] : -_excess[t]; |
664 Capacity d = _excess[s] < -_excess[t] ? _excess[s] : -_excess[t]; |
580 Node u = t; |
665 Node u = t; |