RhodeCode

    /// The type of the map that stores the upper bounds (capacities)

    /// on the arcs.

    /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.

    typedef UM UpperMap;

    /// \brief The type of supply map.

    ///

    /// The type of the map that stores the signed supply values of the 

    /// nodes. 

    /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.

    typedef SM SupplyMap;

    /// \brief The type of the map that stores the flow values.

    ///

    /// The type of the map that stores the flow values.

    /// It must conform to the \ref concepts::ReadWriteMap "ReadWriteMap"

    /// concept.

    /// \brief Instantiates a FlowMap.

    ///

    /// This function instantiates a \ref FlowMap.

    /// \param digraph The digraph for which we would like to define

    /// the flow map.

    static FlowMap* createFlowMap(const Digraph& digraph) {

      return new FlowMap(digraph);

    }

    /// \brief The elevator type used by the algorithm.

    ///

    ///

    /// This function instantiates an \ref Elevator.

    /// \param digraph The digraph for which we would like to define

    /// the elevator.

    /// \param max_level The maximum level of the elevator.

    static Elevator* createElevator(const Digraph& digraph, int max_level) {

      return new Elevator(digraph, max_level);

    }

    /// \brief The tolerance used by the algorithm

    ///

    /// The tolerance used by the algorithm to handle inexact computation.

  };

  /**

     \brief Push-relabel algorithm for the network circulation problem.

     \ingroup max_flow

     This class implements a push-relabel algorithm for the \e network

     \e circulation problem.

     It is to find a feasible circulation when lower and upper bounds

     are given for the flow values on the arcs and lower bounds are

     given for the difference between the outgoing and incoming flow

          typename LM = typename GR::template ArcMap<int>,

          typename UM = LM,

          typename SM = typename GR::template NodeMap<typename UM::Value>,

          typename TR = CirculationDefaultTraits<GR, LM, UM, SM> >

#endif

  class Circulation {

  public:

    ///The \ref CirculationDefaultTraits "traits class" of the algorithm.

    typedef TR Traits;

    ///The type of the digraph the algorithm runs on.

    typedef typename Traits::Digraph Digraph;

    ///The type of the lower bound map.

    typedef typename Traits::LowerMap LowerMap;

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

    typedef typename Traits::UpperMap UpperMap;

    ///The type of the supply map.

    typedef typename Traits::SupplyMap SupplyMap;

    ///The type of the flow map.

    typedef typename Traits::FlowMap FlowMap;

    ///The type of the elevator.

    typedef typename Traits::Elevator Elevator;

    int _node_num;

    const LowerMap *_lo;

    const UpperMap *_up;

    const SupplyMap *_supply;

    FlowMap *_flow;

    bool _local_flow;

    Elevator* _level;

    bool _local_level;

    ExcessMap* _excess;

    Tolerance _tol;

    int _el;

  public:

    typedef Circulation Create;

    ///\name Named Template Parameters

    ///@{

      }

      for (ArcIt e(_g);e!=INVALID;++e) {

        if (!_tol.less(-(*_excess)[_g.target(e)], (*_up)[e])) {

          _flow->set(e, (*_up)[e]);

          (*_excess)[_g.target(e)] += (*_up)[e];

          (*_excess)[_g.source(e)] -= (*_up)[e];

        } else if (_tol.less(-(*_excess)[_g.target(e)], (*_lo)[e])) {

          _flow->set(e, (*_lo)[e]);

          (*_excess)[_g.target(e)] += (*_lo)[e];

          (*_excess)[_g.source(e)] -= (*_lo)[e];

        } else {

          _flow->set(e, fc);

          (*_excess)[_g.target(e)] = 0;

          (*_excess)[_g.source(e)] -= fc;

        }

      }

      _level->initStart();

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

        _level->initAddItem(n);

      _level->initFinish();

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

        if(_tol.positive((*_excess)[n]))

    ///

    ///\sa barrier()

    ///\sa barrierMap()

    bool start()

    {

      Node act;

      Node bact=INVALID;

      Node last_activated=INVALID;

      while((act=_level->highestActive())!=INVALID) {

        int actlevel=(*_level)[act];

        int mlevel=_node_num;

        for(OutArcIt e(_g,act);e!=INVALID; ++e) {

          Node v = _g.target(e);

          if(!_tol.positive(fc)) continue;

          if((*_level)[v]<actlevel) {

            if(!_tol.less(fc, exc)) {

              _flow->set(e, (*_flow)[e] + exc);

              (*_excess)[v] += exc;

              if(!_level->active(v) && _tol.positive((*_excess)[v]))

                _level->activate(v);

              (*_excess)[act] = 0;

              _level->deactivate(act);

              goto next_l;

            }

            else {

              _flow->set(e, (*_up)[e]);

              (*_excess)[v] += fc;

              if(!_level->active(v) && _tol.positive((*_excess)[v]))

                _level->activate(v);

              exc-=fc;

            }

          }

          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];

        }

        for(InArcIt e(_g,act);e!=INVALID; ++e) {

          Node v = _g.source(e);

          if(!_tol.positive(fc)) continue;

          if((*_level)[v]<actlevel) {

            if(!_tol.less(fc, exc)) {

              _flow->set(e, (*_flow)[e] - exc);

              (*_excess)[v] += exc;

              if(!_level->active(v) && _tol.positive((*_excess)[v]))

                _level->activate(v);

              (*_excess)[act] = 0;

              _level->deactivate(act);

              goto next_l;

            }

            else {

    }

    /// @}

    /// \name Query Functions

    /// The results of the circulation algorithm can be obtained using

    /// these functions.\n

    /// Either \ref run() or \ref start() should be called before

    /// using them.

    ///@{

    ///

    ///

    /// \pre Either \ref run() or \ref init() must be called before

    /// using this function.

      return (*_flow)[arc];

    }

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

    ///

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

    ///

    /// \pre Either \ref run() or \ref init() must be called before

    /// using this function.

    const FlowMap& flowMap() const {

      return *_flow;

    }

    ///@{

    ///Check if the found flow is a feasible circulation

    ///Check if the found flow is a feasible circulation,

    ///

    bool checkFlow() const {

      for(ArcIt e(_g);e!=INVALID;++e)

        if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false;

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

        {

          for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e];

          for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e];

          if(_tol.negative(dif)) return false;

        }

      return true;

    }

    ///Check whether or not the last execution provides a barrier

    ///Check whether or not the last execution provides a barrier.

    ///\sa barrier()

    ///\sa barrierMap()

    bool checkBarrier() const

    {

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

        if(barrier(n))

          delta-=(*_supply)[n];

      for(ArcIt e(_g);e!=INVALID;++e)

        {

          Node s=_g.source(e);

          Node t=_g.target(e);

          if(barrier(s)&&!barrier(t)) {

            if (_tol.less(inf_cap - (*_up)[e], delta)) return false;

            delta+=(*_up)[e];

          }

          else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e];

  ///

  /// In general this class is the fastest implementation available

  /// in LEMON for the minimum cost flow problem.

  /// Moreover it supports both directions of the supply/demand inequality

  /// constraints. For more information see \ref SupplyType.

  ///

  /// Most of the parameters of the problem (except for the digraph)

  /// can be given using separate functions, and the algorithm can be

  /// executed using the \ref run() function. If some parameters are not

  /// specified, then default values will be used.

  ///

  /// \tparam GR The digraph type the algorithm runs on.

  /// and supply values in the algorithm. By default it is \c int.

  /// \tparam C The value type used for costs and potentials in the

  ///

  /// \warning Both value types must be signed and all input data must

  /// be integer.

  ///

  /// \note %NetworkSimplex provides five different pivot rule

  /// implementations, from which the most efficient one is used

  /// by default. For more information see \ref PivotRule.

  class NetworkSimplex

  {

  public:

    /// The flow type of the algorithm

    /// The cost type of the algorithm

    typedef C Cost;

#ifdef DOXYGEN

    /// The type of the flow map

    /// The type of the potential map

    typedef GR::NodeMap<Cost> PotentialMap;

#else

    /// The type of the flow map

    /// The type of the potential map

    typedef typename GR::template NodeMap<Cost> PotentialMap;

#endif

  public:

    /// \brief Problem type constants for the \c run() function.

    ///

    /// Enum type containing the problem type constants that can be

    /// returned by the \ref run() function of the algorithm.

    enum ProblemType {

      /// The problem has no feasible solution (flow).

      /// The Altering Candidate List pivot rule.

      /// It is a modified version of the Candidate List method.

      /// It keeps only the several best eligible arcs from the former

      /// candidate list and extends this list in every iteration.

      ALTERING_LIST

    };

  private:

    TEMPLATE_DIGRAPH_TYPEDEFS(GR);

    typedef typename GR::template ArcMap<Cost> CostArcMap;

    typedef std::vector<Arc> ArcVector;

    typedef std::vector<Node> NodeVector;

    typedef std::vector<int> IntVector;

    typedef std::vector<bool> BoolVector;

    typedef std::vector<Cost> CostVector;

    // State constants for arcs

    enum ArcStateEnum {

      STATE_UPPER = -1,

      STATE_TREE  =  0,

      STATE_LOWER =  1

    };

  private:

    // Data related to the underlying digraph

    const GR &_graph;

    int _node_num;

    int _arc_num;

    // Parameters of the problem

    CostArcMap *_pcost;

    bool _pstsup;

    Node _psource, _ptarget;

    SupplyType _stype;

    // Result maps

    FlowMap *_flow_map;

    PotentialMap *_potential_map;

    bool _local_flow;

    bool _local_potential;

    // Data structures for storing the digraph

    IntNodeMap _node_id;

    ArcVector _arc_ref;

    IntVector _source;

    IntVector _target;

    IntVector _rev_thread;

    IntVector _succ_num;

    IntVector _last_succ;

    IntVector _dirty_revs;

    BoolVector _forward;

    IntVector _state;

    int _root;

    // Temporary data used in the current pivot iteration

    int in_arc, join, u_in, v_in, u_out, v_out;

    int first, second, right, last;

    int stem, par_stem, new_stem;

  public:

    /// \brief Constant for infinite upper bounds (capacities).

    ///

    /// Constant for infinite upper bounds (capacities).

  private:

    // Implementation of the First Eligible pivot rule

    class FirstEligiblePivotRule

    {

    private:

      // References to the NetworkSimplex class

      const IntVector  &_source;

      const IntVector  &_target;

      const CostVector &_cost;

    /// \brief Constructor.

    ///

    /// The constructor of the class.

    ///

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

    NetworkSimplex(const GR& graph) :

      _graph(graph),

      _plower(NULL), _pupper(NULL), _pcost(NULL),

      _psupply(NULL), _pstsup(false), _stype(GEQ),

      _flow_map(NULL), _potential_map(NULL),

      _local_flow(false), _local_potential(false),

      _node_id(graph),

    {

      // Check the value types

        "The flow type of NetworkSimplex must be signed");

      LEMON_ASSERT(std::numeric_limits<Cost>::is_signed,

        "The cost type of NetworkSimplex must be signed");

    }

    /// Destructor.

    ~NetworkSimplex() {

      if (_local_flow) delete _flow_map;

      if (_local_potential) delete _potential_map;

    }

    /// \name Parameters

    /// The parameters of the algorithm can be specified using these

    /// functions.

    /// @{

    /// \brief Set the lower bounds on the arcs.

    ///

    /// This function sets the lower bounds on the arcs.

    /// If it is not used before calling \ref run(), the lower bounds

    /// will be set to zero on all arcs.

    ///

    /// \param map An arc map storing the lower bounds.

    /// of the algorithm.

    ///

    /// \return <tt>(*this)</tt>

    template <typename LowerMap>

    NetworkSimplex& lowerMap(const LowerMap& map) {

      delete _plower;

      for (ArcIt a(_graph); a != INVALID; ++a) {

        (*_plower)[a] = map[a];

      }

      return *this;

    }

    /// \brief Set the upper bounds (capacities) on the arcs.

    ///

    /// This function sets the upper bounds (capacities) on the arcs.

    /// If it is not used before calling \ref run(), the upper bounds

    /// will be set to \ref INF on all arcs (i.e. the flow value will be

    /// unbounded from above on each arc).

    ///

    /// \param map An arc map storing the upper bounds.

    /// of the algorithm.

    ///

    /// \return <tt>(*this)</tt>

    template<typename UpperMap>

    NetworkSimplex& upperMap(const UpperMap& map) {

      delete _pupper;

      for (ArcIt a(_graph); a != INVALID; ++a) {

        (*_pupper)[a] = map[a];

      }

      return *this;

    }

    /// \brief Set the costs of the arcs.

    ///

    /// This function sets the costs of the arcs.

    /// If it is not used before calling \ref run(), the costs

    /// will be set to \c 1 on all arcs.

    ///

      }

      return *this;

    }

    /// \brief Set the supply values of the nodes.

    ///

    /// This function sets the supply values of the nodes.

    /// If neither this function nor \ref stSupply() is used before

    /// calling \ref run(), the supply of each node will be set to zero.

    /// (It makes sense only if non-zero lower bounds are given.)

    ///

    /// \param map A node map storing the supply values.

    /// of the algorithm.

    ///

    /// \return <tt>(*this)</tt>

    template<typename SupplyMap>

    NetworkSimplex& supplyMap(const SupplyMap& map) {

      delete _psupply;

      _pstsup = false;

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

        (*_psupply)[n] = map[n];

      }

      return *this;

    }

    /// \brief Set single source and target nodes and a supply value.

    ///

    /// This function sets a single source node and a single target node

    /// and the required flow value.

    /// If neither this function nor \ref supplyMap() is used before

    /// calling \ref run(), the supply of each node will be set to zero.

    /// (It makes sense only if non-zero lower bounds are given.)

    ///

    /// Using this function has the same effect as using \ref supplyMap()

    /// with such a map in which \c k is assigned to \c s, \c -k is

    /// assigned to \c t and all other nodes have zero supply value.

    ///

    /// \param s The source node.

    /// \param t The target node.

    /// \param k 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).

    ///

    /// \return <tt>(*this)</tt>

      delete _psupply;

      _psupply = NULL;

      _pstsup = true;

      _psource = s;

      _ptarget = t;

      _pstflow = k;

      return *this;

    }

    /// \brief Set the type of the supply constraints.

    ///

    /// This function sets the type of the supply/demand constraints.

#ifndef DOXYGEN

    Cost totalCost() const {

      return totalCost<Cost>();

    }

#endif

    /// \brief Return the flow on the given arc.

    ///

    /// This function returns the flow on the given arc.

    ///

    /// \pre \ref run() must be called before using this function.

      return (*_flow_map)[a];

    }

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

    ///

    /// This function returns a const reference to an arc map storing

    /// the found flow.

    ///

    /// \pre \ref run() must be called before using this function.

    const FlowMap& flowMap() const {

      return *_flow_map;

    }

        if (_pcost) {

          for (int i = 0; i != _arc_num; ++i)

            _cost[i] = (*_pcost)[_arc_ref[i]];

        } else {

          for (int i = 0; i != _arc_num; ++i)

            _cost[i] = 1;

        }

      }

      // Remove non-zero lower bounds

      if (_plower) {

        for (int i = 0; i != _arc_num; ++i) {

          if (c > 0) {

            if (_cap[i] < INF) _cap[i] -= c;

            _supply[_source[i]] -= c;

            _supply[_target[i]] += c;

          }

          else if (c < 0) {

            if (_cap[i] < INF + c) {

              _cap[i] -= c;

            } else {

              _cap[i] = INF;

            }

            _supply[_source[i]] -= c;

    bool findLeavingArc() {

      // Initialize first and second nodes according to the direction

      // of the cycle

      if (_state[in_arc] == STATE_LOWER) {

        first  = _source[in_arc];

        second = _target[in_arc];

      } else {

        first  = _target[in_arc];

        second = _source[in_arc];

      }

      delta = _cap[in_arc];

      int result = 0;

      int e;

      // Search the cycle along the path form the first node to the root

      for (int u = first; u != join; u = _parent[u]) {

        e = _pred[u];

        d = _forward[u] ?

          _flow[e] : (_cap[e] == INF ? INF : _cap[e] - _flow[e]);

        if (d < delta) {

          delta = d;

          u_out = u;

          result = 1;

        }

        v_in = second;

      } else {

        u_in = second;

        v_in = first;

      }

      return result != 0;

    }

    // Change _flow and _state vectors

    void changeFlow(bool change) {

      // Augment along the cycle

      if (delta > 0) {

        _flow[in_arc] += val;

        for (int u = _source[in_arc]; u != join; u = _parent[u]) {

          _flow[_pred[u]] += _forward[u] ? -val : val;

        }

        for (int u = _target[in_arc]; u != join; u = _parent[u]) {

          _flow[_pred[u]] += _forward[u] ? val : -val;

        }

      }

      // Update the state of the entering and leaving arcs

      if (change) {

        _state[in_arc] = STATE_TREE;

        _state[_pred[u_out]] =

  struct PreflowDefaultTraits {

    /// \brief The type of the digraph the algorithm runs on.

    typedef GR Digraph;

    /// \brief The type of the map that stores the arc capacities.

    ///

    /// The type of the map that stores the arc capacities.

    /// It must meet the \ref concepts::ReadMap "ReadMap" concept.

    typedef CAP CapacityMap;

    /// \brief The type of the flow values.

    /// \brief The type of the map that stores the flow values.

    ///

    /// The type of the map that stores the flow values.

    /// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept.

    /// \brief Instantiates a FlowMap.

    ///

    /// This function instantiates a \ref FlowMap.

    /// \param digraph The digraph for which we would like to define

    /// the flow map.

    static FlowMap* createFlowMap(const Digraph& digraph) {

      return new FlowMap(digraph);

    }

    /// \brief The elevator type used by Preflow algorithm.

    ///

    ///

    /// This function instantiates an \ref Elevator.

    /// \param digraph The digraph for which we would like to define

    /// the elevator.

    /// \param max_level The maximum level of the elevator.

    static Elevator* createElevator(const Digraph& digraph, int max_level) {

      return new Elevator(digraph, max_level);

    }

    /// \brief The tolerance used by the algorithm

    ///

    /// The tolerance used by the algorithm to handle inexact computation.

  };

  /// \ingroup max_flow

  ///

  /// \brief %Preflow algorithm class.

  ///

  /// This class provides an implementation of Goldberg-Tarjan's \e preflow

  /// \e push-relabel algorithm producing a \ref max_flow

  /// "flow of maximum value" in a digraph.

  /// The preflow algorithms are the fastest known maximum

            typename TR = PreflowDefaultTraits<GR, CAP> >

#endif

  class Preflow {

  public:

    ///The \ref PreflowDefaultTraits "traits class" of the algorithm.

    typedef TR Traits;

    ///The type of the digraph the algorithm runs on.

    typedef typename Traits::Digraph Digraph;

    ///The type of the capacity map.

    typedef typename Traits::CapacityMap CapacityMap;

    ///The type of the flow values.

    ///The type of the flow map.

    typedef typename Traits::FlowMap FlowMap;

    ///The type of the elevator.

    typedef typename Traits::Elevator Elevator;

    ///The type of the tolerance.

    typedef typename Traits::Tolerance Tolerance;

  private:

    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);

    const CapacityMap* _capacity;

    int _node_num;

    Node _source, _target;

    FlowMap* _flow;

    bool _local_flow;

    Elevator* _level;

    bool _local_level;

    ExcessMap* _excess;

    Tolerance _tolerance;

    bool _phase;

    void createStructures() {

      _node_num = countNodes(_graph);

      if (!_flow) {

        _flow = Traits::createFlowMap(_graph);

    /// source node the incoming flow should greater or equal to the

    /// outgoing flow.

    /// \return \c false if the given \c flowMap is not a preflow.

    template <typename FlowMap>

    bool init(const FlowMap& flowMap) {

      createStructures();

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

        _flow->set(e, flowMap[e]);

      }

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

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

          excess += (*_flow)[e];

        }

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

          excess -= (*_flow)[e];

        }

        if (excess < 0 && n != _source) return false;

        (*_excess)[n] = excess;

      }

      typename Digraph::template NodeMap<bool> reached(_graph, false);

            if (!reached[v] && _tolerance.positive((*_flow)[e])) {

              reached[v] = true;

              _level->initAddItem(v);

              nqueue.push_back(v);

            }

          }

        }

        queue.swap(nqueue);

      }

      _level->initFinish();

      for (OutArcIt e(_graph, _source); e != INVALID; ++e) {

        if (_tolerance.positive(rem)) {

          Node u = _graph.target(e);

          if ((*_level)[u] == _level->maxLevel()) continue;

          _flow->set(e, (*_capacity)[e]);

          (*_excess)[u] += rem;

          if (u != _target && !_level->active(u)) {

            _level->activate(u);

          }

        }

      }

      for (InArcIt e(_graph, _source); e != INVALID; ++e) {

        if (_tolerance.positive(rem)) {

          Node v = _graph.source(e);

          if ((*_level)[v] == _level->maxLevel()) continue;

          _flow->set(e, 0);

          (*_excess)[v] += rem;

          if (v != _target && !_level->active(v)) {

            _level->activate(v);

          }

        }

      }

      return true;

    }

    /// minCut() returns a minimum cut.

    /// \pre One of the \ref init() functions must be called before

    /// using this function.

    void startFirstPhase() {

      _phase = true;

      Node n = _level->highestActive();

      int level = _level->highestActiveLevel();

      while (n != INVALID) {

        int num = _node_num;

        while (num > 0 && n != INVALID) {

          int new_level = _level->maxLevel();

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

            if (!_tolerance.positive(rem)) continue;

            Node v = _graph.target(e);

            if ((*_level)[v] < level) {

              if (!_level->active(v) && v != _target) {

                _level->activate(v);

              }

              if (!_tolerance.less(rem, excess)) {

                _flow->set(e, (*_flow)[e] + excess);

                (*_excess)[v] += excess;

                excess = 0;

                goto no_more_push_1;

              } else {

                excess -= rem;

                (*_excess)[v] += rem;

                _flow->set(e, (*_capacity)[e]);

              }

            } else if (new_level > (*_level)[v]) {

              new_level = (*_level)[v];

            }

          }

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

            if (!_tolerance.positive(rem)) continue;

            Node v = _graph.source(e);

            if ((*_level)[v] < level) {

              if (!_level->active(v) && v != _target) {

                _level->activate(v);

              }

              if (!_tolerance.less(rem, excess)) {

                _flow->set(e, (*_flow)[e] - excess);

                (*_excess)[v] += excess;

                excess = 0;

                goto no_more_push_1;

              } else {

            }

          } else {

            _level->deactivate(n);

          }

          n = _level->highestActive();

          level = _level->highestActiveLevel();

          --num;

        }

        num = _node_num * 20;

        while (num > 0 && n != INVALID) {

          int new_level = _level->maxLevel();

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

            if (!_tolerance.positive(rem)) continue;

            Node v = _graph.target(e);

            if ((*_level)[v] < level) {

              if (!_level->active(v) && v != _target) {

                _level->activate(v);

              }

              if (!_tolerance.less(rem, excess)) {

                _flow->set(e, (*_flow)[e] + excess);

                (*_excess)[v] += excess;

                excess = 0;

                goto no_more_push_2;

              } else {

                excess -= rem;

                (*_excess)[v] += rem;

                _flow->set(e, (*_capacity)[e]);

              }

            } else if (new_level > (*_level)[v]) {

              new_level = (*_level)[v];

            }

          }

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

            if (!_tolerance.positive(rem)) continue;

            Node v = _graph.source(e);

            if ((*_level)[v] < level) {

              if (!_level->active(v) && v != _target) {

                _level->activate(v);

              }

              if (!_tolerance.less(rem, excess)) {

                _flow->set(e, (*_flow)[e] - excess);

                (*_excess)[v] += excess;

                excess = 0;

                goto no_more_push_2;

              } else {

      _level->initFinish();

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

        if (!reached[n]) {

          _level->dirtyTopButOne(n);

        } else if ((*_excess)[n] > 0 && _target != n) {

          _level->activate(n);

        }

      }

      Node n;

      while ((n = _level->highestActive()) != INVALID) {

        int level = _level->highestActiveLevel();

        int new_level = _level->maxLevel();

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

          if (!_tolerance.positive(rem)) continue;

          Node v = _graph.target(e);

          if ((*_level)[v] < level) {

            if (!_level->active(v) && v != _source) {

              _level->activate(v);

            }

            if (!_tolerance.less(rem, excess)) {

              _flow->set(e, (*_flow)[e] + excess);

              (*_excess)[v] += excess;

              excess = 0;

              goto no_more_push;

            } else {

              excess -= rem;

              (*_excess)[v] += rem;

              _flow->set(e, (*_capacity)[e]);

            }

          } else if (new_level > (*_level)[v]) {

            new_level = (*_level)[v];

          }

        }

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

          if (!_tolerance.positive(rem)) continue;

          Node v = _graph.source(e);

          if ((*_level)[v] < level) {

            if (!_level->active(v) && v != _source) {

              _level->activate(v);

            }

            if (!_tolerance.less(rem, excess)) {

              _flow->set(e, (*_flow)[e] - excess);

              (*_excess)[v] += excess;

              excess = 0;

              goto no_more_push;

            } else {

    /// before using them.

    ///@{

    /// \brief Returns the value of the maximum flow.

    ///

    /// Returns the value of the maximum flow by returning the excess

    /// of the target node. This value equals to the value of

    /// the maximum flow already after the first phase of the algorithm.

    ///

    /// \pre Either \ref run() or \ref init() must be called before

    /// using this function.

      return (*_excess)[_target];

    }

    ///