/// \brief The capacity scaling algorithm for finding a minimum cost flow.

 /// \brief Capacity scaling algorithm for finding a minimum cost flow.

 #include <vector>

   /// \brief Implementation of the capacity scaling version of the

   /// \brief Implementation of the capacity scaling algorithm for

   /// successive shortest path algorithm for finding a minimum cost

   /// finding a minimum cost flow.

   /// flow.

   /// \param Graph The directed graph type the algorithm runs on.

   /// \tparam Graph The directed graph type the algorithm runs on.

   /// \param LowerMap The type of the lower bound map.

   /// \tparam LowerMap The type of the lower bound map.

   /// \param CapacityMap The type of the capacity (upper bound) map.

   /// \tparam CapacityMap The type of the capacity (upper bound) map.

   /// \param CostMap The type of the cost (length) map.

   /// \tparam CostMap The type of the cost (length) map.

   /// \param SupplyMap The type of the supply map.

   /// \tparam SupplyMap The type of the supply map.

   /// - Edge capacities and costs should be non-negative integers.

   /// - Edge capacities and costs should be \e non-negative \e integers.

   ///   However \c CostMap::Value should be signed type.

   /// - Supply values should be \e signed \e integers.

   /// - Supply values should be signed integers.

   /// - \c LowerMap::Value must be convertible to \c CapacityMap::Value.

   /// - \c LowerMap::Value must be convertible to

   /// - \c CapacityMap::Value and \c SupplyMap::Value must be

   ///   \c CapacityMap::Value and \c CapacityMap::Value must be

   ///   convertible to each other.

   ///   convertible to \c SupplyMap::Value.

   /// - All value types must be convertible to \c CostMap::Value, which

   ///   must be signed type.

              typename CapacityMap = LowerMap,

              typename CapacityMap = typename Graph::template EdgeMap<int>,

              typename SupplyMap = typename Graph::template NodeMap

              typename SupplyMap = typename Graph::template NodeMap<int> >

                                   <typename CapacityMap::Value> >

   protected:

   private:

     /// for finding shortest paths in the residual network of the graph

     /// for finding shortest paths in the residual network.

     /// with respect to the reduced edge costs and modifying the

     ///

     /// node potentials according to the distance of the nodes.

     /// \ref ResidualDijkstra is a special implementation of the

     /// \ref Dijkstra algorithm for finding shortest paths in the

     /// residual network of the graph with respect to the reduced edge

     /// costs and modifying the node potentials according to the

     /// distance of the nodes.

     protected:

     private:

       /// The directed graph the algorithm runs on.

       // The directed graph the algorithm runs on

       const Graph &graph;

       const Graph &_graph;

       /// The flow map.

       // The main maps

       const FlowMap &flow;

       const FlowMap &_flow;

       /// The residual capacity map.

       const CapacityEdgeMap &_res_cap;

       const CapacityEdgeMap &res_cap;

       const CostMap &_cost;

       /// The cost map.

       const SupplyNodeMap &_excess;

       const CostMap &cost;

       PotentialMap &_potential;

       /// The excess map.

       const SupplyNodeMap &excess;

       // The distance map

       CostNodeMap _dist;

       /// The potential map.

       // The pred edge map

       PotentialMap &potential;

       PredMap &_pred;

       // The processed (i.e. permanently labeled) nodes

       /// The distance map.

       std::vector<Node> _proc_nodes;

       CostNodeMap dist;

       /// The map of predecessors edges.

       PredMap &pred;

       /// The processed (i.e. permanently labeled) nodes.

       std::vector<Node> proc_nodes;

       ResidualDijkstra( const Graph &_graph,

       ResidualDijkstra( const Graph &graph,

                         const FlowMap &_flow,

                         const FlowMap &flow,

                         const CapacityEdgeMap &_res_cap,

                         const CapacityEdgeMap &res_cap,

                         const CostMap &_cost,

                         const CostMap &cost,

                         const SupplyMap &_excess,

                         const SupplyMap &excess,

                         PotentialMap &_potential,

                         PotentialMap &potential,

                         PredMap &_pred ) :

                         PredMap &pred ) :

         graph(_graph), flow(_flow), res_cap(_res_cap), cost(_cost),

         _graph(graph), _flow(flow), _res_cap(res_cap), _cost(cost),

         excess(_excess), potential(_potential), dist(_graph),

         _excess(excess), _potential(potential), _dist(graph),

         pred(_pred)

         _pred(pred)

         HeapCrossRef heap_cross_ref(graph, Heap::PRE_HEAP);

         HeapCrossRef heap_cross_ref(_graph, Heap::PRE_HEAP);

         pred[s] = INVALID;

         _pred[s] = INVALID;

         proc_nodes.clear();

         _proc_nodes.clear();

         while (!heap.empty() && excess[heap.top()] > -delta) {

         while (!heap.empty() && _excess[heap.top()] > -delta) {

           Cost d = heap.prio() - potential[u], nd;

           Cost d = heap.prio() + _potential[u], nd;

           dist[u] = heap.prio();

           _dist[u] = heap.prio();

           proc_nodes.push_back(u);

           _proc_nodes.push_back(u);

           for (OutEdgeIt e(graph, u); e != INVALID; ++e) {

           for (OutEdgeIt e(_graph, u); e != INVALID; ++e) {

             if (res_cap[e] >= delta) {

             if (_res_cap[e] >= delta) {

               v = graph.target(e);

               v = _graph.target(e);

                 heap.push(v, d + cost[e] + potential[v]);

                 heap.push(v, d + _cost[e] - _potential[v]);

                 pred[v] = e;

                 _pred[v] = e;

                 nd = d + cost[e] + potential[v];

                 nd = d + _cost[e] - _potential[v];

                   pred[v] = e;

                   _pred[v] = e;

           for (InEdgeIt e(graph, u); e != INVALID; ++e) {

           for (InEdgeIt e(_graph, u); e != INVALID; ++e) {

             if (flow[e] >= delta) {

             if (_flow[e] >= delta) {

               v = graph.source(e);

               v = _graph.source(e);

                 heap.push(v, d - cost[e] + potential[v]);

                 heap.push(v, d - _cost[e] - _potential[v]);

                 pred[v] = e;

                 _pred[v] = e;

                 nd = d - cost[e] + potential[v];

                 nd = d - _cost[e] - _potential[v];

                   pred[v] = e;

                   _pred[v] = e;

         Cost dt = heap.prio();

         Cost t_dist = heap.prio();

         for (int i = 0; i < proc_nodes.size(); ++i)

         for (int i = 0; i < int(_proc_nodes.size()); ++i)

           potential[proc_nodes[i]] -= dist[proc_nodes[i]] - dt;

           _potential[_proc_nodes[i]] += _dist[_proc_nodes[i]] - t_dist;

   protected:

   private:

     /// The directed graph the algorithm runs on.

     // The directed graph the algorithm runs on

     const Graph &graph;

     const Graph &_graph;

     /// The original lower bound map.

     // The original lower bound map

     const LowerMap *lower;

     const LowerMap *_lower;

     /// The modified capacity map.

     // The modified capacity map

     CapacityEdgeMap capacity;

     CapacityEdgeMap _capacity;

     /// The cost map.

     // The original cost map

     const CostMap &cost;

     const CostMap &_cost;

     /// The modified supply map.

     // The modified supply map

     SupplyNodeMap supply;

     SupplyNodeMap _supply;

     bool valid_supply;

     bool _valid_supply;

     /// The edge map of the current flow.

     // Edge map of the current flow

     FlowMap flow;

     FlowMap _flow;

     /// The potential node map.

     // Node map of the current potentials

     PotentialMap potential;

     PotentialMap _potential;

     /// The residual capacity map.

     // The residual capacity map

     CapacityEdgeMap res_cap;

     CapacityEdgeMap _res_cap;

     /// The excess map.

     // The excess map

     SupplyNodeMap excess;

     SupplyNodeMap _excess;

     /// The excess nodes (i.e. the nodes with positive excess).

     // The excess nodes (i.e. nodes with positive excess)

     std::vector<Node> excess_nodes;

     std::vector<Node> _excess_nodes;

     /// The deficit nodes (i.e. the nodes with negative excess).

     // The deficit nodes (i.e. nodes with negative excess)

     std::vector<Node> deficit_nodes;

     std::vector<Node> _deficit_nodes;

     /// The scaling status (enabled or disabled).

     // The delta parameter used for capacity scaling

     ScalingEnum scaling;

     Capacity _delta;

     /// The \c delta parameter used for capacity scaling.

     // The maximum number of phases

     Capacity delta;

     int _phase_num;

     /// The maximum number of phases.

     int phase_num;

     // The pred edge map

     PredMap _pred;

     /// \brief Implementation of the \ref Dijkstra algorithm for

     // Implementation of the Dijkstra algorithm for finding augmenting

     /// finding augmenting shortest paths in the residual network.

     // shortest paths in the residual network

     ResidualDijkstra dijkstra;

     ResidualDijkstra _dijkstra;

     /// The map of predecessors edges.

     PredMap pred;

     /// \param _graph The directed graph the algorithm runs on.

     /// \param graph The directed graph the algorithm runs on.

     /// \param _lower The lower bounds of the edges.

     /// \param lower The lower bounds of the edges.

     /// \param _capacity The capacities (upper bounds) of the edges.

     /// \param capacity The capacities (upper bounds) of the edges.

     /// \param _cost The cost (length) values of the edges.

     /// \param cost The cost (length) values of the edges.

     /// \param _supply The supply values of the nodes (signed).

     /// \param supply The supply values of the nodes (signed).

     CapacityScaling( const Graph &_graph,

     CapacityScaling( const Graph &graph,

                      const LowerMap &_lower,

                      const LowerMap &lower,

                      const CapacityMap &_capacity,

                      const CapacityMap &capacity,

                      const CostMap &_cost,

                      const CostMap &cost,

                      const SupplyMap &_supply ) :

                      const SupplyMap &supply ) :

       graph(_graph), lower(&_lower), capacity(_graph), cost(_cost),

       _graph(graph), _lower(&lower), _capacity(graph), _cost(cost),

       supply(_graph), flow(_graph, 0), potential(_graph, 0),

       _supply(graph), _flow(graph, 0), _potential(graph, 0),

       res_cap(_graph), excess(_graph), pred(_graph),

       _res_cap(graph), _excess(graph), _pred(graph),

       dijkstra(graph, flow, res_cap, cost, excess, potential, pred)

       _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred)

       capacity = subMap(_capacity, _lower);

       _capacity = subMap(capacity, lower);

       res_cap = capacity;

       _res_cap = _capacity;

       for (NodeIt n(graph); n != INVALID; ++n) {

       for (NodeIt n(_graph); n != INVALID; ++n) {

         Supply s = _supply[n];

         Supply s = supply[n];

         for (InEdgeIt e(graph, n); e != INVALID; ++e)

         for (InEdgeIt e(_graph, n); e != INVALID; ++e)

           s += _lower[e];

           s += lower[e];

         for (OutEdgeIt e(graph, n); e != INVALID; ++e)

         for (OutEdgeIt e(_graph, n); e != INVALID; ++e)

           s -= _lower[e];

           s -= lower[e];

         supply[n] = s;

         _supply[n] = s;

       valid_supply = sum == 0;

       _valid_supply = sum == 0;

     /// \param _graph The directed graph the algorithm runs on.

     /// \param graph The directed graph the algorithm runs on.

     /// \param _capacity The capacities (upper bounds) of the edges.

     /// \param capacity The capacities (upper bounds) of the edges.

     /// \param _cost The cost (length) values of the edges.

     /// \param cost The cost (length) values of the edges.

     /// \param _supply The supply values of the nodes (signed).

     /// \param supply The supply values of the nodes (signed).

     CapacityScaling( const Graph &_graph,

     CapacityScaling( const Graph &graph,

                      const CapacityMap &_capacity,

                      const CapacityMap &capacity,

                      const CostMap &_cost,

                      const CostMap &cost,

                      const SupplyMap &_supply ) :

                      const SupplyMap &supply ) :

       graph(_graph), lower(NULL), capacity(_capacity), cost(_cost),

       _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost),

       supply(_supply), flow(_graph, 0), potential(_graph, 0),

       _supply(supply), _flow(graph, 0), _potential(graph, 0),

       res_cap(_capacity), excess(_graph), pred(_graph),

       _res_cap(capacity), _excess(graph), _pred(graph),

       dijkstra(graph, flow, res_cap, cost, excess, potential)

       _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred)

       for (NodeIt n(graph); n != INVALID; ++n) sum += supply[n];

       for (NodeIt n(_graph); n != INVALID; ++n) sum += _supply[n];

       valid_supply = sum == 0;

       _valid_supply = sum == 0;

     /// \param _graph The directed graph the algorithm runs on.

     /// \param graph The directed graph the algorithm runs on.

     /// \param _lower The lower bounds of the edges.

     /// \param lower The lower bounds of the edges.

     /// \param _capacity The capacities (upper bounds) of the edges.

     /// \param capacity The capacities (upper bounds) of the edges.

     /// \param _cost The cost (length) values of the edges.

     /// \param cost The cost (length) values of the edges.

     /// \param _s The source node.

     /// \param s The source node.

     /// \param _t The target node.

     /// \param t The target node.

     /// \param _flow_value The required amount of flow from node \c _s

     /// \param flow_value The required amount of flow from node \c s

     /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t).

     /// to node \c t (i.e. the supply of \c s and the demand of \c t).

     CapacityScaling( const Graph &_graph,

     CapacityScaling( const Graph &graph,

                      const LowerMap &_lower,

                      const LowerMap &lower,

                      const CapacityMap &_capacity,

                      const CapacityMap &capacity,

                      const CostMap &_cost,

                      const CostMap &cost,

                      Node _s, Node _t,

                      Node s, Node t,

                      Supply _flow_value ) :

                      Supply flow_value ) :

       graph(_graph), lower(&_lower), capacity(_graph), cost(_cost),

       _graph(graph), _lower(&lower), _capacity(graph), _cost(cost),

       supply(_graph), flow(_graph, 0), potential(_graph, 0),

       _supply(graph), _flow(graph, 0), _potential(graph, 0),

       res_cap(_graph), excess(_graph), pred(_graph),

       _res_cap(graph), _excess(graph), _pred(graph),

       dijkstra(graph, flow, res_cap, cost, excess, potential)

       _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred)

       capacity = subMap(_capacity, _lower);

       _capacity = subMap(capacity, lower);

       res_cap = capacity;

       _res_cap = _capacity;

       for (NodeIt n(graph); n != INVALID; ++n) {

       for (NodeIt n(_graph); n != INVALID; ++n) {

         Supply s = 0;

         Supply sum = 0;

         if (n == _s) s =  _flow_value;

         if (n == s) sum =  flow_value;

         if (n == _t) s = -_flow_value;

         if (n == t) sum = -flow_value;

         for (InEdgeIt e(graph, n); e != INVALID; ++e)

         for (InEdgeIt e(_graph, n); e != INVALID; ++e)

           s += _lower[e];

           sum += lower[e];

         for (OutEdgeIt e(graph, n); e != INVALID; ++e)

         for (OutEdgeIt e(_graph, n); e != INVALID; ++e)

           s -= _lower[e];

           sum -= lower[e];

         supply[n] = s;

         _supply[n] = sum;

       }

       }

       valid_supply = true;

       _valid_supply = true;

     /// \param _graph The directed graph the algorithm runs on.

     /// \param graph The directed graph the algorithm runs on.

     /// \param _capacity The capacities (upper bounds) of the edges.

     /// \param capacity The capacities (upper bounds) of the edges.

     /// \param _cost The cost (length) values of the edges.

     /// \param cost The cost (length) values of the edges.

     /// \param _s The source node.

     /// \param s The source node.

     /// \param _t The target node.

     /// \param t The target node.

     /// \param _flow_value The required amount of flow from node \c _s

     /// \param flow_value The required amount of flow from node \c s

     /// to node \c _t (i.e. the supply of \c _s and the demand of \c _t).

     /// to node \c t (i.e. the supply of \c s and the demand of \c t).

     CapacityScaling( const Graph &_graph,

     CapacityScaling( const Graph &graph,

                      const CapacityMap &_capacity,

                      const CapacityMap &capacity,

                      const CostMap &_cost,

                      const CostMap &cost,

                      Node _s, Node _t,

                      Node s, Node t,

                      Supply _flow_value ) :

                      Supply flow_value ) :

       graph(_graph), lower(NULL), capacity(_capacity), cost(_cost),

       _graph(graph), _lower(NULL), _capacity(capacity), _cost(cost),

       supply(_graph, 0), flow(_graph, 0), potential(_graph, 0),

       _supply(graph, 0), _flow(graph, 0), _potential(graph, 0),

       res_cap(_capacity), excess(_graph), pred(_graph),

       _res_cap(capacity), _excess(graph), _pred(graph),

       dijkstra(graph, flow, res_cap, cost, excess, potential)

       _dijkstra(_graph, _flow, _res_cap, _cost, _excess, _potential, _pred)

       supply[_s] =  _flow_value;

       _supply[s] =  flow_value;

       supply[_t] = -_flow_value;

       _supply[t] = -flow_value;

       valid_supply = true;

       _valid_supply = true;

     /// \param scaling_mode The scaling mode. In case of WITH_SCALING

     /// \param scaling Enable or disable capacity scaling.

     /// capacity scaling is enabled in the algorithm (this is the

     /// default) otherwise it is disabled.

     /// is small, the algorithm could be slightly faster without

     /// is rather small, the algorithm could be slightly faster without

     bool run(int scaling_mode = WITH_SCALING) {

     bool run(bool scaling = true) {

       return init(scaling_mode) && start();

       return init(scaling) && start();

     }

     }

     /// \brief Returns a const reference to the flow map.

     /// \brief Returns a const reference to the edge map storing the

     ///

     /// found flow.

     /// Returns a const reference to the flow map.

     ///

     /// Returns a const reference to the edge map storing the found flow.

       return flow;

       return _flow;

     }

     }

     /// \brief Returns a const reference to the potential map (the dual

     /// \brief Returns a const reference to the node map storing the

     /// solution).

     /// found potentials (the dual solution).

     ///

     ///

     /// Returns a const reference to the potential map (the dual

     /// Returns a const reference to the node map storing the found

     /// solution).

     /// potentials (the dual solution).

       return potential;

       return _potential;

       for (EdgeIt e(graph); e != INVALID; ++e)

       for (EdgeIt e(_graph); e != INVALID; ++e)

         c += flow[e] * cost[e];

         c += _flow[e] * _cost[e];

   protected:

   private:

     bool init(int scaling_mode) {

     bool init(bool scaling) {

       if (!valid_supply) return false;

       if (!_valid_supply) return false;

       excess = supply;

       _excess = _supply;

       if (scaling_mode == WITH_SCALING) {

       if (scaling) {

         for (NodeIt n(graph); n != INVALID; ++n) {

         for (NodeIt n(_graph); n != INVALID; ++n) {

           if ( supply[n] > max_sup) max_sup =  supply[n];

           if ( _supply[n] > max_sup) max_sup =  _supply[n];

           if (-supply[n] > max_dem) max_dem = -supply[n];

           if (-_supply[n] > max_dem) max_dem = -_supply[n];

         phase_num = 0;

         _phase_num = 0;

         for (delta = 1; 2 * delta <= max_sup; delta *= 2)

         for (_delta = 1; 2 * _delta <= max_sup; _delta *= 2)

           ++phase_num;

           ++_phase_num;

         delta = 1;

         _delta = 1;

       if (delta > 1)

       if (_delta > 1)

     /// \brief Executes the capacity scaling version of the successive

     /// Executes the capacity scaling algorithm.

     /// shortest path algorithm.

         for (EdgeIt e(graph); e != INVALID; ++e) {

         for (EdgeIt e(_graph); e != INVALID; ++e) {

           Node u = graph.source(e), v = graph.target(e);

           Node u = _graph.source(e), v = _graph.target(e);

           Cost c = cost[e] - potential[u] + potential[v];

           Cost c = _cost[e] + _potential[u] - _potential[v];

           if (c < 0 && res_cap[e] >= delta) {

           if (c < 0 && _res_cap[e] >= _delta) {

             excess[u] -= res_cap[e];

             _excess[u] -= _res_cap[e];

             excess[v] += res_cap[e];

             _excess[v] += _res_cap[e];

             flow[e] = capacity[e];

             _flow[e] = _capacity[e];

             res_cap[e] = 0;

             _res_cap[e] = 0;

           }

           }

           else if (c > 0 && flow[e] >= delta) {

           else if (c > 0 && _flow[e] >= _delta) {

             excess[u] += flow[e];

             _excess[u] += _flow[e];

             excess[v] -= flow[e];

             _excess[v] -= _flow[e];

             flow[e] = 0;

             _flow[e] = 0;

             res_cap[e] = capacity[e];

             _res_cap[e] = _capacity[e];

         excess_nodes.clear();

         _excess_nodes.clear();

         deficit_nodes.clear();

         _deficit_nodes.clear();

         for (NodeIt n(graph); n != INVALID; ++n) {

         for (NodeIt n(_graph); n != INVALID; ++n) {

           if (excess[n] >=  delta) excess_nodes.push_back(n);

           if (_excess[n] >=  _delta) _excess_nodes.push_back(n);

           if (excess[n] <= -delta) deficit_nodes.push_back(n);

           if (_excess[n] <= -_delta) _deficit_nodes.push_back(n);

         while (next_node < excess_nodes.size()) {

         while (next_node < int(_excess_nodes.size())) {

           if (delta > 1) {

           if (_delta > 1) {

             for (int i = 0; i < deficit_nodes.size(); ++i) {

             for (int i = 0; i < int(_deficit_nodes.size()); ++i) {

               if (excess[deficit_nodes[i]] <= -delta) {

               if (_excess[_deficit_nodes[i]] <= -_delta) {

           s = excess_nodes[next_node];

           s = _excess_nodes[next_node];

           if ((t = dijkstra.run(s, delta)) == INVALID) {

           if ((t = _dijkstra.run(s, _delta)) == INVALID) {

             if (delta > 1) {

             if (_delta > 1) {

           Capacity d = excess[s] < -excess[t] ? excess[s] : -excess[t];

           Capacity d = _excess[s] < -_excess[t] ? _excess[s] : -_excess[t];

           if (d > delta) {

           if (d > _delta) {

             while ((e = pred[u]) != INVALID) {

             while ((e = _pred[u]) != INVALID) {

               if (u == graph.target(e)) {

               if (u == _graph.target(e)) {

                 rc = res_cap[e];

                 rc = _res_cap[e];

                 u = graph.source(e);

                 u = _graph.source(e);

                 rc = flow[e];

                 rc = _flow[e];

                 u = graph.target(e);

                 u = _graph.target(e);

           while ((e = pred[u]) != INVALID) {

           while ((e = _pred[u]) != INVALID) {

             if (u == graph.target(e)) {

             if (u == _graph.target(e)) {

               flow[e] += d;

               _flow[e] += d;

               res_cap[e] -= d;

               _res_cap[e] -= d;

               u = graph.source(e);

               u = _graph.source(e);

               flow[e] -= d;

               _flow[e] -= d;

               res_cap[e] += d;

               _res_cap[e] += d;

               u = graph.target(e);

               u = _graph.target(e);

           excess[s] -= d;

           _excess[s] -= d;

           excess[t] += d;

           _excess[t] += d;

           if (excess[s] < delta) ++next_node;

           if (_excess[s] < _delta) ++next_node;

         if (delta == 1) break;

         if (_delta == 1) break;

         if (++phase_cnt > phase_num / 4) factor = 2;

         if (++phase_cnt > _phase_num / 4) factor = 2;

         delta = delta <= factor ? 1 : delta / factor;

         _delta = _delta <= factor ? 1 : _delta / factor;

       if (lower) {

       if (_lower) {

         for (EdgeIt e(graph); e != INVALID; ++e)

         for (EdgeIt e(_graph); e != INVALID; ++e)

           flow[e] += (*lower)[e];

           _flow[e] += (*_lower)[e];

     /// \brief Executes the successive shortest path algorithm without

     /// Executes the successive shortest path algorithm.

     /// capacity scaling.

       for (NodeIt n(graph); n != INVALID; ++n)

       for (NodeIt n(_graph); n != INVALID; ++n)

         if (excess[n] > 0) excess_nodes.push_back(n);

         if (_excess[n] > 0) _excess_nodes.push_back(n);

       if (excess_nodes.size() == 0) return true;

       if (_excess_nodes.size() == 0) return true;

       while ( excess[excess_nodes[next_node]] > 0 ||

       while ( _excess[_excess_nodes[next_node]] > 0 ||

               ++next_node < excess_nodes.size() )

               ++next_node < int(_excess_nodes.size()) )

         s = excess_nodes[next_node];

         s = _excess_nodes[next_node];

         if ((t = dijkstra.run(s, 1)) == INVALID)

         if ((t = _dijkstra.run(s, 1)) == INVALID)

         Capacity d = excess[s] < -excess[t] ? excess[s] : -excess[t];

         Capacity d = _excess[s] < -_excess[t] ? _excess[s] : -_excess[t];

         while ((e = pred[u]) != INVALID) {

         while ((e = _pred[u]) != INVALID) {

           if (u == graph.target(e)) {

           if (u == _graph.target(e)) {

             rc = res_cap[e];

             rc = _res_cap[e];

             u = graph.source(e);

             u = _graph.source(e);

             rc = flow[e];

             rc = _flow[e];

             u = graph.target(e);

             u = _graph.target(e);

         while ((e = pred[u]) != INVALID) {

         while ((e = _pred[u]) != INVALID) {

           if (u == graph.target(e)) {

           if (u == _graph.target(e)) {

             flow[e] += d;

             _flow[e] += d;

             res_cap[e] -= d;

             _res_cap[e] -= d;

             u = graph.source(e);

             u = _graph.source(e);

             flow[e] -= d;

             _flow[e] -= d;

             res_cap[e] += d;

             _res_cap[e] += d;

             u = graph.target(e);

             u = _graph.target(e);

         excess[s] -= d;

         _excess[s] -= d;

         excess[t] += d;

         _excess[t] += d;

       if (lower) {

       if (_lower) {

         for (EdgeIt e(graph); e != INVALID; ++e)

         for (EdgeIt e(_graph); e != INVALID; ++e)

           flow[e] += (*lower)[e];

           _flow[e] += (*_lower)[e];