1 /* -*- mode: C++; indent-tabs-mode: nil; -*-
3 * This file is a part of LEMON, a generic C++ optimization library.
5 * Copyright (C) 2003-2010
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
9 * Permission to use, modify and distribute this software is granted
10 * provided that this copyright notice appears in all copies. For
11 * precise terms see the accompanying LICENSE file.
13 * This software is provided "AS IS" with no warranty of any kind,
14 * express or implied, and with no claim as to its suitability for any
19 #ifndef LEMON_CAPACITY_SCALING_H
20 #define LEMON_CAPACITY_SCALING_H
22 /// \ingroup min_cost_flow_algs
25 /// \brief Capacity Scaling algorithm for finding a minimum cost flow.
29 #include <lemon/core.h>
30 #include <lemon/bin_heap.h>
34 /// \brief Default traits class of CapacityScaling algorithm.
36 /// Default traits class of CapacityScaling algorithm.
37 /// \tparam GR Digraph type.
38 /// \tparam V The number type used for flow amounts, capacity bounds
39 /// and supply values. By default it is \c int.
40 /// \tparam C The number type used for costs and potentials.
41 /// By default it is the same as \c V.
42 template <typename GR, typename V = int, typename C = V>
43 struct CapacityScalingDefaultTraits
45 /// The type of the digraph
47 /// The type of the flow amounts, capacity bounds and supply values
49 /// The type of the arc costs
52 /// \brief The type of the heap used for internal Dijkstra computations.
54 /// The type of the heap used for internal Dijkstra computations.
55 /// It must conform to the \ref lemon::concepts::Heap "Heap" concept,
56 /// its priority type must be \c Cost and its cross reference type
57 /// must be \ref RangeMap "RangeMap<int>".
58 typedef BinHeap<Cost, RangeMap<int> > Heap;
61 /// \addtogroup min_cost_flow_algs
64 /// \brief Implementation of the Capacity Scaling algorithm for
65 /// finding a \ref min_cost_flow "minimum cost flow".
67 /// \ref CapacityScaling implements the capacity scaling version
68 /// of the successive shortest path algorithm for finding a
69 /// \ref min_cost_flow "minimum cost flow" \ref amo93networkflows,
70 /// \ref edmondskarp72theoretical. It is an efficient dual
73 /// This algorithm is typically slower than \ref CostScaling and
74 /// \ref NetworkSimplex, but in special cases, it can be more
75 /// efficient than them.
76 /// (For more information, see \ref min_cost_flow_algs "the module page".)
78 /// Most of the parameters of the problem (except for the digraph)
79 /// can be given using separate functions, and the algorithm can be
80 /// executed using the \ref run() function. If some parameters are not
81 /// specified, then default values will be used.
83 /// \tparam GR The digraph type the algorithm runs on.
84 /// \tparam V The number type used for flow amounts, capacity bounds
85 /// and supply values in the algorithm. By default, it is \c int.
86 /// \tparam C The number type used for costs and potentials in the
87 /// algorithm. By default, it is the same as \c V.
88 /// \tparam TR The traits class that defines various types used by the
89 /// algorithm. By default, it is \ref CapacityScalingDefaultTraits
90 /// "CapacityScalingDefaultTraits<GR, V, C>".
91 /// In most cases, this parameter should not be set directly,
92 /// consider to use the named template parameters instead.
94 /// \warning Both \c V and \c C must be signed number types.
95 /// \warning Capacity bounds and supply values must be integer, but
96 /// arc costs can be arbitrary real numbers.
97 /// \warning This algorithm does not support negative costs for
98 /// arcs having infinite upper bound.
100 template <typename GR, typename V, typename C, typename TR>
102 template < typename GR, typename V = int, typename C = V,
103 typename TR = CapacityScalingDefaultTraits<GR, V, C> >
105 class CapacityScaling
109 /// The type of the digraph
110 typedef typename TR::Digraph Digraph;
111 /// The type of the flow amounts, capacity bounds and supply values
112 typedef typename TR::Value Value;
113 /// The type of the arc costs
114 typedef typename TR::Cost Cost;
116 /// The type of the heap used for internal Dijkstra computations
117 typedef typename TR::Heap Heap;
119 /// The \ref CapacityScalingDefaultTraits "traits class" of the algorithm
124 /// \brief Problem type constants for the \c run() function.
126 /// Enum type containing the problem type constants that can be
127 /// returned by the \ref run() function of the algorithm.
129 /// The problem has no feasible solution (flow).
131 /// The problem has optimal solution (i.e. it is feasible and
132 /// bounded), and the algorithm has found optimal flow and node
133 /// potentials (primal and dual solutions).
135 /// The digraph contains an arc of negative cost and infinite
136 /// upper bound. It means that the objective function is unbounded
137 /// on that arc, however, note that it could actually be bounded
138 /// over the feasible flows, but this algroithm cannot handle
145 TEMPLATE_DIGRAPH_TYPEDEFS(GR);
147 typedef std::vector<int> IntVector;
148 typedef std::vector<Value> ValueVector;
149 typedef std::vector<Cost> CostVector;
150 typedef std::vector<char> BoolVector;
151 // Note: vector<char> is used instead of vector<bool> for efficiency reasons
155 // Data related to the underlying digraph
162 // Parameters of the problem
166 // Data structures for storing the digraph
170 IntVector _first_out;
182 ValueVector _res_cap;
185 IntVector _excess_nodes;
186 IntVector _deficit_nodes;
194 /// \brief Constant for infinite upper bounds (capacities).
196 /// Constant for infinite upper bounds (capacities).
197 /// It is \c std::numeric_limits<Value>::infinity() if available,
198 /// \c std::numeric_limits<Value>::max() otherwise.
203 // Special implementation of the Dijkstra algorithm for finding
204 // shortest paths in the residual network of the digraph with
205 // respect to the reduced arc costs and modifying the node
206 // potentials according to the found distance labels.
207 class ResidualDijkstra
213 const IntVector &_first_out;
214 const IntVector &_target;
215 const CostVector &_cost;
216 const ValueVector &_res_cap;
217 const ValueVector &_excess;
221 IntVector _proc_nodes;
226 ResidualDijkstra(CapacityScaling& cs) :
227 _node_num(cs._node_num), _geq(cs._sum_supply < 0),
228 _first_out(cs._first_out), _target(cs._target), _cost(cs._cost),
229 _res_cap(cs._res_cap), _excess(cs._excess), _pi(cs._pi),
230 _pred(cs._pred), _dist(cs._node_num)
233 int run(int s, Value delta = 1) {
234 RangeMap<int> heap_cross_ref(_node_num, Heap::PRE_HEAP);
235 Heap heap(heap_cross_ref);
241 while (!heap.empty() && _excess[heap.top()] > -delta) {
242 int u = heap.top(), v;
243 Cost d = heap.prio() + _pi[u], dn;
244 _dist[u] = heap.prio();
245 _proc_nodes.push_back(u);
248 // Traverse outgoing residual arcs
249 int last_out = _geq ? _first_out[u+1] : _first_out[u+1] - 1;
250 for (int a = _first_out[u]; a != last_out; ++a) {
251 if (_res_cap[a] < delta) continue;
253 switch (heap.state(v)) {
255 heap.push(v, d + _cost[a] - _pi[v]);
259 dn = d + _cost[a] - _pi[v];
261 heap.decrease(v, dn);
265 case Heap::POST_HEAP:
270 if (heap.empty()) return -1;
272 // Update potentials of processed nodes
274 Cost dt = heap.prio();
275 for (int i = 0; i < int(_proc_nodes.size()); ++i) {
276 _pi[_proc_nodes[i]] += _dist[_proc_nodes[i]] - dt;
282 }; //class ResidualDijkstra
286 /// \name Named Template Parameters
289 template <typename T>
290 struct SetHeapTraits : public Traits {
294 /// \brief \ref named-templ-param "Named parameter" for setting
297 /// \ref named-templ-param "Named parameter" for setting \c Heap
298 /// type, which is used for internal Dijkstra computations.
299 /// It must conform to the \ref lemon::concepts::Heap "Heap" concept,
300 /// its priority type must be \c Cost and its cross reference type
301 /// must be \ref RangeMap "RangeMap<int>".
302 template <typename T>
304 : public CapacityScaling<GR, V, C, SetHeapTraits<T> > {
305 typedef CapacityScaling<GR, V, C, SetHeapTraits<T> > Create;
316 /// \brief Constructor.
318 /// The constructor of the class.
320 /// \param graph The digraph the algorithm runs on.
321 CapacityScaling(const GR& graph) :
322 _graph(graph), _node_id(graph), _arc_idf(graph), _arc_idb(graph),
323 INF(std::numeric_limits<Value>::has_infinity ?
324 std::numeric_limits<Value>::infinity() :
325 std::numeric_limits<Value>::max())
327 // Check the number types
328 LEMON_ASSERT(std::numeric_limits<Value>::is_signed,
329 "The flow type of CapacityScaling must be signed");
330 LEMON_ASSERT(std::numeric_limits<Cost>::is_signed,
331 "The cost type of CapacityScaling must be signed");
333 // Reset data structures
338 /// The parameters of the algorithm can be specified using these
343 /// \brief Set the lower bounds on the arcs.
345 /// This function sets the lower bounds on the arcs.
346 /// If it is not used before calling \ref run(), the lower bounds
347 /// will be set to zero on all arcs.
349 /// \param map An arc map storing the lower bounds.
350 /// Its \c Value type must be convertible to the \c Value type
351 /// of the algorithm.
353 /// \return <tt>(*this)</tt>
354 template <typename LowerMap>
355 CapacityScaling& lowerMap(const LowerMap& map) {
357 for (ArcIt a(_graph); a != INVALID; ++a) {
358 _lower[_arc_idf[a]] = map[a];
363 /// \brief Set the upper bounds (capacities) on the arcs.
365 /// This function sets the upper bounds (capacities) on the arcs.
366 /// If it is not used before calling \ref run(), the upper bounds
367 /// will be set to \ref INF on all arcs (i.e. the flow value will be
368 /// unbounded from above).
370 /// \param map An arc map storing the upper bounds.
371 /// Its \c Value type must be convertible to the \c Value type
372 /// of the algorithm.
374 /// \return <tt>(*this)</tt>
375 template<typename UpperMap>
376 CapacityScaling& upperMap(const UpperMap& map) {
377 for (ArcIt a(_graph); a != INVALID; ++a) {
378 _upper[_arc_idf[a]] = map[a];
383 /// \brief Set the costs of the arcs.
385 /// This function sets the costs of the arcs.
386 /// If it is not used before calling \ref run(), the costs
387 /// will be set to \c 1 on all arcs.
389 /// \param map An arc map storing the costs.
390 /// Its \c Value type must be convertible to the \c Cost type
391 /// of the algorithm.
393 /// \return <tt>(*this)</tt>
394 template<typename CostMap>
395 CapacityScaling& costMap(const CostMap& map) {
396 for (ArcIt a(_graph); a != INVALID; ++a) {
397 _cost[_arc_idf[a]] = map[a];
398 _cost[_arc_idb[a]] = -map[a];
403 /// \brief Set the supply values of the nodes.
405 /// This function sets the supply values of the nodes.
406 /// If neither this function nor \ref stSupply() is used before
407 /// calling \ref run(), the supply of each node will be set to zero.
409 /// \param map A node map storing the supply values.
410 /// Its \c Value type must be convertible to the \c Value type
411 /// of the algorithm.
413 /// \return <tt>(*this)</tt>
414 template<typename SupplyMap>
415 CapacityScaling& supplyMap(const SupplyMap& map) {
416 for (NodeIt n(_graph); n != INVALID; ++n) {
417 _supply[_node_id[n]] = map[n];
422 /// \brief Set single source and target nodes and a supply value.
424 /// This function sets a single source node and a single target node
425 /// and the required flow value.
426 /// If neither this function nor \ref supplyMap() is used before
427 /// calling \ref run(), the supply of each node will be set to zero.
429 /// Using this function has the same effect as using \ref supplyMap()
430 /// with a map in which \c k is assigned to \c s, \c -k is
431 /// assigned to \c t and all other nodes have zero supply value.
433 /// \param s The source node.
434 /// \param t The target node.
435 /// \param k The required amount of flow from node \c s to node \c t
436 /// (i.e. the supply of \c s and the demand of \c t).
438 /// \return <tt>(*this)</tt>
439 CapacityScaling& stSupply(const Node& s, const Node& t, Value k) {
440 for (int i = 0; i != _node_num; ++i) {
443 _supply[_node_id[s]] = k;
444 _supply[_node_id[t]] = -k;
450 /// \name Execution control
451 /// The algorithm can be executed using \ref run().
455 /// \brief Run the algorithm.
457 /// This function runs the algorithm.
458 /// The paramters can be specified using functions \ref lowerMap(),
459 /// \ref upperMap(), \ref costMap(), \ref supplyMap(), \ref stSupply().
462 /// CapacityScaling<ListDigraph> cs(graph);
463 /// cs.lowerMap(lower).upperMap(upper).costMap(cost)
464 /// .supplyMap(sup).run();
467 /// This function can be called more than once. All the given parameters
468 /// are kept for the next call, unless \ref resetParams() or \ref reset()
469 /// is used, thus only the modified parameters have to be set again.
470 /// If the underlying digraph was also modified after the construction
471 /// of the class (or the last \ref reset() call), then the \ref reset()
472 /// function must be called.
474 /// \param factor The capacity scaling factor. It must be larger than
475 /// one to use scaling. If it is less or equal to one, then scaling
476 /// will be disabled.
478 /// \return \c INFEASIBLE if no feasible flow exists,
479 /// \n \c OPTIMAL if the problem has optimal solution
480 /// (i.e. it is feasible and bounded), and the algorithm has found
481 /// optimal flow and node potentials (primal and dual solutions),
482 /// \n \c UNBOUNDED if the digraph contains an arc of negative cost
483 /// and infinite upper bound. It means that the objective function
484 /// is unbounded on that arc, however, note that it could actually be
485 /// bounded over the feasible flows, but this algroithm cannot handle
489 /// \see resetParams(), reset()
490 ProblemType run(int factor = 4) {
492 ProblemType pt = init();
493 if (pt != OPTIMAL) return pt;
497 /// \brief Reset all the parameters that have been given before.
499 /// This function resets all the paramaters that have been given
500 /// before using functions \ref lowerMap(), \ref upperMap(),
501 /// \ref costMap(), \ref supplyMap(), \ref stSupply().
503 /// It is useful for multiple \ref run() calls. Basically, all the given
504 /// parameters are kept for the next \ref run() call, unless
505 /// \ref resetParams() or \ref reset() is used.
506 /// If the underlying digraph was also modified after the construction
507 /// of the class or the last \ref reset() call, then the \ref reset()
508 /// function must be used, otherwise \ref resetParams() is sufficient.
512 /// CapacityScaling<ListDigraph> cs(graph);
515 /// cs.lowerMap(lower).upperMap(upper).costMap(cost)
516 /// .supplyMap(sup).run();
518 /// // Run again with modified cost map (resetParams() is not called,
519 /// // so only the cost map have to be set again)
521 /// cs.costMap(cost).run();
523 /// // Run again from scratch using resetParams()
524 /// // (the lower bounds will be set to zero on all arcs)
525 /// cs.resetParams();
526 /// cs.upperMap(capacity).costMap(cost)
527 /// .supplyMap(sup).run();
530 /// \return <tt>(*this)</tt>
532 /// \see reset(), run()
533 CapacityScaling& resetParams() {
534 for (int i = 0; i != _node_num; ++i) {
537 for (int j = 0; j != _res_arc_num; ++j) {
540 _cost[j] = _forward[j] ? 1 : -1;
546 /// \brief Reset the internal data structures and all the parameters
547 /// that have been given before.
549 /// This function resets the internal data structures and all the
550 /// paramaters that have been given before using functions \ref lowerMap(),
551 /// \ref upperMap(), \ref costMap(), \ref supplyMap(), \ref stSupply().
553 /// It is useful for multiple \ref run() calls. Basically, all the given
554 /// parameters are kept for the next \ref run() call, unless
555 /// \ref resetParams() or \ref reset() is used.
556 /// If the underlying digraph was also modified after the construction
557 /// of the class or the last \ref reset() call, then the \ref reset()
558 /// function must be used, otherwise \ref resetParams() is sufficient.
560 /// See \ref resetParams() for examples.
562 /// \return <tt>(*this)</tt>
564 /// \see resetParams(), run()
565 CapacityScaling& reset() {
567 _node_num = countNodes(_graph);
568 _arc_num = countArcs(_graph);
569 _res_arc_num = 2 * (_arc_num + _node_num);
573 _first_out.resize(_node_num + 1);
574 _forward.resize(_res_arc_num);
575 _source.resize(_res_arc_num);
576 _target.resize(_res_arc_num);
577 _reverse.resize(_res_arc_num);
579 _lower.resize(_res_arc_num);
580 _upper.resize(_res_arc_num);
581 _cost.resize(_res_arc_num);
582 _supply.resize(_node_num);
584 _res_cap.resize(_res_arc_num);
585 _pi.resize(_node_num);
586 _excess.resize(_node_num);
587 _pred.resize(_node_num);
590 int i = 0, j = 0, k = 2 * _arc_num + _node_num - 1;
591 for (NodeIt n(_graph); n != INVALID; ++n, ++i) {
595 for (NodeIt n(_graph); n != INVALID; ++n, ++i) {
597 for (OutArcIt a(_graph, n); a != INVALID; ++a, ++j) {
601 _target[j] = _node_id[_graph.runningNode(a)];
603 for (InArcIt a(_graph, n); a != INVALID; ++a, ++j) {
607 _target[j] = _node_id[_graph.runningNode(a)];
620 _first_out[_node_num] = k;
621 for (ArcIt a(_graph); a != INVALID; ++a) {
622 int fi = _arc_idf[a];
623 int bi = _arc_idb[a];
635 /// \name Query Functions
636 /// The results of the algorithm can be obtained using these
638 /// The \ref run() function must be called before using them.
642 /// \brief Return the total cost of the found flow.
644 /// This function returns the total cost of the found flow.
645 /// Its complexity is O(e).
647 /// \note The return type of the function can be specified as a
648 /// template parameter. For example,
650 /// cs.totalCost<double>();
652 /// It is useful if the total cost cannot be stored in the \c Cost
653 /// type of the algorithm, which is the default return type of the
656 /// \pre \ref run() must be called before using this function.
657 template <typename Number>
658 Number totalCost() const {
660 for (ArcIt a(_graph); a != INVALID; ++a) {
662 c += static_cast<Number>(_res_cap[i]) *
663 (-static_cast<Number>(_cost[i]));
669 Cost totalCost() const {
670 return totalCost<Cost>();
674 /// \brief Return the flow on the given arc.
676 /// This function returns the flow on the given arc.
678 /// \pre \ref run() must be called before using this function.
679 Value flow(const Arc& a) const {
680 return _res_cap[_arc_idb[a]];
683 /// \brief Copy the flow values (the primal solution) into the
686 /// This function copies the flow value on each arc into the given
687 /// map. The \c Value type of the algorithm must be convertible to
688 /// the \c Value type of the map.
690 /// \pre \ref run() must be called before using this function.
691 template <typename FlowMap>
692 void flowMap(FlowMap &map) const {
693 for (ArcIt a(_graph); a != INVALID; ++a) {
694 map.set(a, _res_cap[_arc_idb[a]]);
698 /// \brief Return the potential (dual value) of the given node.
700 /// This function returns the potential (dual value) of the
703 /// \pre \ref run() must be called before using this function.
704 Cost potential(const Node& n) const {
705 return _pi[_node_id[n]];
708 /// \brief Copy the potential values (the dual solution) into the
711 /// This function copies the potential (dual value) of each node
712 /// into the given map.
713 /// The \c Cost type of the algorithm must be convertible to the
714 /// \c Value type of the map.
716 /// \pre \ref run() must be called before using this function.
717 template <typename PotentialMap>
718 void potentialMap(PotentialMap &map) const {
719 for (NodeIt n(_graph); n != INVALID; ++n) {
720 map.set(n, _pi[_node_id[n]]);
728 // Initialize the algorithm
730 if (_node_num <= 1) return INFEASIBLE;
732 // Check the sum of supply values
734 for (int i = 0; i != _root; ++i) {
735 _sum_supply += _supply[i];
737 if (_sum_supply > 0) return INFEASIBLE;
739 // Check lower and upper bounds
740 LEMON_DEBUG(checkBoundMaps(),
741 "Upper bounds must be greater or equal to the lower bounds");
744 // Initialize vectors
745 for (int i = 0; i != _root; ++i) {
747 _excess[i] = _supply[i];
750 // Remove non-zero lower bounds
751 const Value MAX = std::numeric_limits<Value>::max();
754 for (int i = 0; i != _root; ++i) {
755 last_out = _first_out[i+1];
756 for (int j = _first_out[i]; j != last_out; ++j) {
760 _res_cap[j] = _upper[j] < MAX ? _upper[j] - c : INF;
762 _res_cap[j] = _upper[j] < MAX + c ? _upper[j] - c : INF;
765 _excess[_target[j]] += c;
772 for (int j = 0; j != _res_arc_num; ++j) {
773 _res_cap[j] = _forward[j] ? _upper[j] : 0;
777 // Handle negative costs
778 for (int i = 0; i != _root; ++i) {
779 last_out = _first_out[i+1] - 1;
780 for (int j = _first_out[i]; j != last_out; ++j) {
781 Value rc = _res_cap[j];
782 if (_cost[j] < 0 && rc > 0) {
783 if (rc >= MAX) return UNBOUNDED;
785 _excess[_target[j]] += rc;
787 _res_cap[_reverse[j]] += rc;
792 // Handle GEQ supply type
793 if (_sum_supply < 0) {
795 _excess[_root] = -_sum_supply;
796 for (int a = _first_out[_root]; a != _res_arc_num; ++a) {
797 int ra = _reverse[a];
798 _res_cap[a] = -_sum_supply + 1;
806 for (int a = _first_out[_root]; a != _res_arc_num; ++a) {
807 int ra = _reverse[a];
815 // Initialize delta value
818 Value max_sup = 0, max_dem = 0, max_cap = 0;
819 for (int i = 0; i != _root; ++i) {
820 Value ex = _excess[i];
821 if ( ex > max_sup) max_sup = ex;
822 if (-ex > max_dem) max_dem = -ex;
823 int last_out = _first_out[i+1] - 1;
824 for (int j = _first_out[i]; j != last_out; ++j) {
825 if (_res_cap[j] > max_cap) max_cap = _res_cap[j];
828 max_sup = std::min(std::min(max_sup, max_dem), max_cap);
829 for (_delta = 1; 2 * _delta <= max_sup; _delta *= 2) ;
838 // Check if the upper bound is greater than or equal to the lower bound
839 // on each forward arc.
840 bool checkBoundMaps() {
841 for (int j = 0; j != _res_arc_num; ++j) {
842 if (_forward[j] && _upper[j] < _lower[j]) return false;
847 ProblemType start() {
848 // Execute the algorithm
851 pt = startWithScaling();
853 pt = startWithoutScaling();
855 // Handle non-zero lower bounds
857 int limit = _first_out[_root];
858 for (int j = 0; j != limit; ++j) {
859 if (_forward[j]) _res_cap[_reverse[j]] += _lower[j];
863 // Shift potentials if necessary
864 Cost pr = _pi[_root];
865 if (_sum_supply < 0 || pr > 0) {
866 for (int i = 0; i != _node_num; ++i) {
874 // Execute the capacity scaling algorithm
875 ProblemType startWithScaling() {
876 // Perform capacity scaling phases
878 ResidualDijkstra _dijkstra(*this);
880 // Saturate all arcs not satisfying the optimality condition
882 for (int u = 0; u != _node_num; ++u) {
883 last_out = _sum_supply < 0 ?
884 _first_out[u+1] : _first_out[u+1] - 1;
885 for (int a = _first_out[u]; a != last_out; ++a) {
887 Cost c = _cost[a] + _pi[u] - _pi[v];
888 Value rc = _res_cap[a];
889 if (c < 0 && rc >= _delta) {
893 _res_cap[_reverse[a]] += rc;
898 // Find excess nodes and deficit nodes
899 _excess_nodes.clear();
900 _deficit_nodes.clear();
901 for (int u = 0; u != _node_num; ++u) {
902 Value ex = _excess[u];
903 if (ex >= _delta) _excess_nodes.push_back(u);
904 if (ex <= -_delta) _deficit_nodes.push_back(u);
906 int next_node = 0, next_def_node = 0;
908 // Find augmenting shortest paths
909 while (next_node < int(_excess_nodes.size())) {
910 // Check deficit nodes
912 bool delta_deficit = false;
913 for ( ; next_def_node < int(_deficit_nodes.size());
915 if (_excess[_deficit_nodes[next_def_node]] <= -_delta) {
916 delta_deficit = true;
920 if (!delta_deficit) break;
923 // Run Dijkstra in the residual network
924 s = _excess_nodes[next_node];
925 if ((t = _dijkstra.run(s, _delta)) == -1) {
933 // Augment along a shortest path from s to t
934 Value d = std::min(_excess[s], -_excess[t]);
938 while ((a = _pred[u]) != -1) {
939 if (_res_cap[a] < d) d = _res_cap[a];
944 while ((a = _pred[u]) != -1) {
946 _res_cap[_reverse[a]] += d;
952 if (_excess[s] < _delta) ++next_node;
955 if (_delta == 1) break;
956 _delta = _delta <= _factor ? 1 : _delta / _factor;
962 // Execute the successive shortest path algorithm
963 ProblemType startWithoutScaling() {
965 _excess_nodes.clear();
966 for (int i = 0; i != _node_num; ++i) {
967 if (_excess[i] > 0) _excess_nodes.push_back(i);
969 if (_excess_nodes.size() == 0) return OPTIMAL;
972 // Find shortest paths
974 ResidualDijkstra _dijkstra(*this);
975 while ( _excess[_excess_nodes[next_node]] > 0 ||
976 ++next_node < int(_excess_nodes.size()) )
978 // Run Dijkstra in the residual network
979 s = _excess_nodes[next_node];
980 if ((t = _dijkstra.run(s)) == -1) return INFEASIBLE;
982 // Augment along a shortest path from s to t
983 Value d = std::min(_excess[s], -_excess[t]);
987 while ((a = _pred[u]) != -1) {
988 if (_res_cap[a] < d) d = _res_cap[a];
993 while ((a = _pred[u]) != -1) {
995 _res_cap[_reverse[a]] += d;
1005 }; //class CapacityScaling
1011 #endif //LEMON_CAPACITY_SCALING_H