RhodeCode

/* -*- mode: C++; indent-tabs-mode: nil; -*-

 *

 * This file is a part of LEMON, a generic C++ optimization library.

 *

 * Copyright (C) 2003-2009

 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport

 * (Egervary Research Group on Combinatorial Optimization, EGRES).

 *

 * Permission to use, modify and distribute this software is granted

 * provided that this copyright notice appears in all copies. For

 * precise terms see the accompanying LICENSE file.

 *

 * This software is provided "AS IS" with no warranty of any kind,

 * express or implied, and with no claim as to its suitability for any

 * purpose.

 *

 */

namespace lemon {

/**

\page min_cost_flow Minimum Cost Flow Problem

\section mcf_def Definition (GEQ form)

The \e minimum \e cost \e flow \e problem is to find a feasible flow of

minimum total cost from a set of supply nodes to a set of demand nodes

in a network with capacity constraints (lower and upper bounds)

and arc costs.

Formally, let \f$G=(V,A)\f$ be a digraph, \f$lower: A\rightarrow\mathbf{R}\f$,

\f$upper: A\rightarrow\mathbf{R}\cup\{+\infty\}\f$ denote the lower and

upper bounds for the flow values on the arcs, for which

\f$lower(uv) \leq upper(uv)\f$ must hold for all \f$uv\in A\f$,

\f$cost: A\rightarrow\mathbf{R}\f$ denotes the cost per unit flow

on the arcs and \f$sup: V\rightarrow\mathbf{R}\f$ denotes the

signed supply values of the nodes.

If \f$sup(u)>0\f$, then \f$u\f$ is a supply node with \f$sup(u)\f$

supply, if \f$sup(u)<0\f$, then \f$u\f$ is a demand node with

\f$-sup(u)\f$ demand.

A minimum cost flow is an \f$f: A\rightarrow\mathbf{R}\f$ solution

of the following optimization problem.

\f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]

\f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \geq

    sup(u) \quad \forall u\in V \f]

\f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]

The sum of the supply values, i.e. \f$\sum_{u\in V} sup(u)\f$ must be

zero or negative in order to have a feasible solution (since the sum

of the expressions on the left-hand side of the inequalities is zero).

It means that the total demand must be greater or equal to the total

supply and all the supplies have to be carried out from the supply nodes,

but there could be demands that are not satisfied.

If \f$\sum_{u\in V} sup(u)\f$ is zero, then all the supply/demand

constraints have to be satisfied with equality, i.e. all demands

have to be satisfied and all supplies have to be used.

\section mcf_algs Algorithms

LEMON contains several algorithms for solving this problem, for more

information see \ref min_cost_flow_algs "Minimum Cost Flow Algorithms".

A feasible solution for this problem can be found using \ref Circulation.

\section mcf_dual Dual Solution

The dual solution of the minimum cost flow problem is represented by

node potentials \f$\pi: V\rightarrow\mathbf{R}\f$.

An \f$f: A\rightarrow\mathbf{R}\f$ primal feasible solution is optimal

if and only if for some \f$\pi: V\rightarrow\mathbf{R}\f$ node potentials

the following \e complementary \e slackness optimality conditions hold.

 - For all \f$uv\in A\f$ arcs:

   - if \f$cost^\pi(uv)>0\f$, then \f$f(uv)=lower(uv)\f$;

   - if \f$lower(uv)<f(uv)<upper(uv)\f$, then \f$cost^\pi(uv)=0\f$;

   - if \f$cost^\pi(uv)<0\f$, then \f$f(uv)=upper(uv)\f$.

 - For all \f$u\in V\f$ nodes:

   - if \f$\sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \neq sup(u)\f$,

     then \f$\pi(u)=0\f$.

Here \f$cost^\pi(uv)\f$ denotes the \e reduced \e cost of the arc

\f$uv\in A\f$ with respect to the potential function \f$\pi\f$, i.e.

\f[ cost^\pi(uv) = cost(uv) + \pi(u) - \pi(v).\f]

All algorithms provide dual solution (node potentials), as well,

if an optimal flow is found.

\section mcf_eq Equality Form

The above \ref mcf_def "definition" is actually more general than the

usual formulation of the minimum cost flow problem, in which strict

equalities are required in the supply/demand contraints.

\f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]

\f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) =

    sup(u) \quad \forall u\in V \f]

\f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]

However if the sum of the supply values is zero, then these two problems

are equivalent.

The \ref min_cost_flow_algs "algorithms" in LEMON support the general

form, so if you need the equality form, you have to ensure this additional

contraint manually.

\section mcf_leq Opposite Inequalites (LEQ Form)

Another possible definition of the minimum cost flow problem is

when there are <em>"less or equal"</em> (LEQ) supply/demand constraints,

instead of the <em>"greater or equal"</em> (GEQ) constraints.

\f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]

\f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \leq

    sup(u) \quad \forall u\in V \f]

\f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]

It means that the total demand must be less or equal to the 

total supply (i.e. \f$\sum_{u\in V} sup(u)\f$ must be zero or

positive) and all the demands have to be satisfied, but there

could be supplies that are not carried out from the supply

nodes.

The equality form is also a special case of this form, of course.

You could easily transform this case to the \ref mcf_def "GEQ form"

of the problem by reversing the direction of the arcs and taking the

negative of the supply values (e.g. using \ref ReverseDigraph and

\ref NegMap adaptors).

However \ref NetworkSimplex algorithm also supports this form directly

for the sake of convenience.

Note that the optimality conditions for this supply constraint type are

slightly differ from the conditions that are discussed for the GEQ form,

namely the potentials have to be non-negative instead of non-positive.

An \f$f: A\rightarrow\mathbf{R}\f$ feasible solution of this problem

is optimal if and only if for some \f$\pi: V\rightarrow\mathbf{R}\f$

node potentials the following conditions hold.

 - For all \f$uv\in A\f$ arcs:

   - if \f$cost^\pi(uv)>0\f$, then \f$f(uv)=lower(uv)\f$;

   - if \f$lower(uv)<f(uv)<upper(uv)\f$, then \f$cost^\pi(uv)=0\f$;

   - if \f$cost^\pi(uv)<0\f$, then \f$f(uv)=upper(uv)\f$.

 - For all \f$u\in V\f$ nodes:

   - if \f$\sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \neq sup(u)\f$,

     then \f$\pi(u)=0\f$.

*/

}

    typedef BellmanFordDefaultOperationTraits<Value> OperationTraits;

    /// \brief The type of the map that stores the last arcs of the 

    /// shortest paths.

    /// 

    /// The type of the map that stores the last

    /// arcs of the shortest paths.

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

    typedef typename GR::template NodeMap<typename GR::Arc> PredMap;

    /// \brief Instantiates a \c PredMap.

    /// 

    /// This function instantiates a \ref PredMap. 

    /// \param g is the digraph to which we would like to define the

    /// \ref PredMap.

    static PredMap *createPredMap(const GR& g) {

      return new PredMap(g);

    }

    /// \brief The type of the map that stores the distances of the nodes.

    ///

    /// The type of the map that stores the distances of the nodes.

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

    typedef typename GR::template NodeMap<typename LEN::Value> DistMap;

    /// \brief Instantiates a \c DistMap.

    ///

    /// This function instantiates a \ref DistMap. 

    /// \param g is the digraph to which we would like to define the 

    /// \ref DistMap.

    static DistMap *createDistMap(const GR& g) {

      return new DistMap(g);

    }

  };

  /// \brief %BellmanFord algorithm class.

  ///

  /// \ingroup shortest_path

  /// This class provides an efficient implementation of the Bellman-Ford 

  /// algorithm. The maximum time complexity of the algorithm is

  /// <tt>O(ne)</tt>.

  ///

  /// The Bellman-Ford algorithm solves the single-source shortest path

  /// problem when the arcs can have negative lengths, but the digraph

  /// should not contain directed cycles with negative total length.

  /// If all arc costs are non-negative, consider to use the Dijkstra

  /// algorithm instead, since it is more efficient.

  ///

  /// The arc lengths are passed to the algorithm using a

  /// \ref concepts::ReadMap "ReadMap", so it is easy to change it to any 

  /// kind of length. The type of the length values is determined by the

  /// \ref concepts::ReadMap::Value "Value" type of the length map.

  ///

  /// There is also a \ref bellmanFord() "function-type interface" for the

  /// Bellman-Ford algorithm, which is convenient in the simplier cases and

  /// it can be used easier.

  ///

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

  /// The default type is \ref ListDigraph.

  /// \tparam LEN A \ref concepts::ReadMap "readable" arc map that specifies

  /// the lengths of the arcs. The default map type is

  /// \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".

#ifdef DOXYGEN

  template <typename GR, typename LEN, typename TR>

#else

  template <typename GR=ListDigraph,

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

            typename TR=BellmanFordDefaultTraits<GR,LEN> >

#endif

  class BellmanFord {

  public:

    ///The type of the underlying digraph.

    typedef typename TR::Digraph Digraph;

    /// \brief The type of the arc lengths.

    typedef typename TR::LengthMap::Value Value;

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

    typedef typename TR::LengthMap LengthMap;

    /// \brief The type of the map that stores the last

    /// arcs of the shortest paths.

    typedef typename TR::PredMap PredMap;

    /// \brief The type of the map that stores the distances of the nodes.

    typedef typename TR::DistMap DistMap;

    /// The type of the paths.

    typedef PredMapPath<Digraph, PredMap> Path;

    ///\brief The \ref BellmanFordDefaultOperationTraits

    /// "operation traits class" of the algorithm.

    typedef typename TR::OperationTraits OperationTraits;

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

    typedef TR Traits;

  private:

    typedef typename Digraph::Node Node;

    typedef typename Digraph::NodeIt NodeIt;

    typedef typename Digraph::Arc Arc;

    typedef typename Digraph::OutArcIt OutArcIt;

    // Pointer to the underlying digraph.

    const Digraph *_gr;

    // Pointer to the length map

    const LengthMap *_length;

    // Pointer to the map of predecessors arcs.

    PredMap *_pred;

    // Indicates if _pred is locally allocated (true) or not.

    bool _local_pred;

    // Pointer to the map of distances.

    DistMap *_dist;

    // Indicates if _dist is locally allocated (true) or not.

    bool _local_dist;

    typedef typename Digraph::template NodeMap<bool> MaskMap;

    MaskMap *_mask;

    std::vector<Node> _process;

    // Creates the maps if necessary.

    void create_maps() {

      if(!_pred) {

	_local_pred = true;

	_pred = Traits::createPredMap(*_gr);

      }

      if(!_dist) {

	_local_dist = true;

	_dist = Traits::createDistMap(*_gr);

      }

      _mask = new MaskMap(*_gr, false);

    }

  public :

    typedef BellmanFord Create;

    /// \name Named Template Parameters

    ///@{

    template <class T>

    struct SetPredMapTraits : public Traits {

      typedef T PredMap;

      static PredMap *createPredMap(const Digraph&) {

        LEMON_ASSERT(false, "PredMap is not initialized");

        return 0; // ignore warnings

      }

    };

    /// \brief \ref named-templ-param "Named parameter" for setting

    /// \c PredMap type.

    ///

    /// \ref named-templ-param "Named parameter" for setting

    /// \c PredMap type.

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

    template <class T>

    struct SetPredMap 

      : public BellmanFord< Digraph, LengthMap, SetPredMapTraits<T> > {

      typedef BellmanFord< Digraph, LengthMap, SetPredMapTraits<T> > Create;

    };

    template <class T>

    struct SetDistMapTraits : public Traits {

      typedef T DistMap;

      static DistMap *createDistMap(const Digraph&) {

        LEMON_ASSERT(false, "DistMap is not initialized");

        return 0; // ignore warnings

      }

    };

    /// \brief \ref named-templ-param "Named parameter" for setting

    /// \c DistMap type.

    ///

    /// \ref named-templ-param "Named parameter" for setting

    /// \c DistMap type.

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

    template <class T>

    struct SetDistMap 

      : public BellmanFord< Digraph, LengthMap, SetDistMapTraits<T> > {

      typedef BellmanFord< Digraph, LengthMap, SetDistMapTraits<T> > Create;

    };

    template <class T>

    struct SetOperationTraitsTraits : public Traits {

      typedef T OperationTraits;

    };

    /// \brief \ref named-templ-param "Named parameter" for setting 

    /// \c OperationTraits type.

    ///

    /// \ref named-templ-param "Named parameter" for setting

    /// \c OperationTraits type.

    template <class T>

    struct SetOperationTraits

      : public BellmanFord< Digraph, LengthMap, SetOperationTraitsTraits<T> > {

      typedef BellmanFord< Digraph, LengthMap, SetOperationTraitsTraits<T> >

      Create;

    };

    ///@}

  protected:

    BellmanFord() {}

  public:      

    /// \brief Constructor.

    ///

    /// Constructor.

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

    /// \param length The length map used by the algorithm.

    BellmanFord(const Digraph& g, const LengthMap& length) :

      _gr(&g), _length(&length),

      _pred(0), _local_pred(false),

      _dist(0), _local_dist(false), _mask(0) {}

    ///Destructor.

    ~BellmanFord() {

      if(_local_pred) delete _pred;

      if(_local_dist) delete _dist;

      if(_mask) delete _mask;

    }

    /// \brief Sets the length map.

    ///

    /// Sets the length map.

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

    BellmanFord &lengthMap(const LengthMap &map) {

      _length = ↦

      return *this;

    }

    /// \brief Sets the map that stores the predecessor arcs.

    ///

    /// Sets the map that stores the predecessor arcs.

    /// If you don't use this function before calling \ref run()

    /// or \ref init(), an instance will be allocated automatically.

    /// The destructor deallocates this automatically allocated map,

    /// of course.

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

    BellmanFord &predMap(PredMap &map) {

      if(_local_pred) {

	delete _pred;

	_local_pred=false;

      }

      _pred = ↦

      return *this;

    }

    /// \brief Sets the map that stores the distances of the nodes.

    ///

    /// Sets the map that stores the distances of the nodes calculated

    /// by the algorithm.

    /// If you don't use this function before calling \ref run()

    /// or \ref init(), an instance will be allocated automatically.

    /// The destructor deallocates this automatically allocated map,

    /// of course.

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

    BellmanFord &distMap(DistMap &map) {

      if(_local_dist) {

	delete _dist;

	_local_dist=false;

      }

      _dist = ↦

      return *this;

    }

    /// \name Execution Control

    /// The simplest way to execute the Bellman-Ford algorithm is to use

    /// one of the member functions called \ref run().\n

    /// If you need better control on the execution, you have to call

    /// \ref init() first, then you can add several source nodes

    /// with \ref addSource(). Finally the actual path computation can be

    /// performed with \ref start(), \ref checkedStart() or

    /// \ref limitedStart().

    ///@{

    /// \brief Initializes the internal data structures.

    /// 

    /// Initializes the internal data structures. The optional parameter

    /// is the initial distance of each node.

    void init(const Value value = OperationTraits::infinity()) {

      create_maps();

      for (NodeIt it(*_gr); it != INVALID; ++it) {

	_pred->set(it, INVALID);

	_dist->set(it, value);

      }

      _process.clear();

      if (OperationTraits::less(value, OperationTraits::infinity())) {

	for (NodeIt it(*_gr); it != INVALID; ++it) {

	  _process.push_back(it);

	  _mask->set(it, true);

	}

      }

    }

    /// \brief Adds a new source node.

    ///

    /// This function adds a new source node. The optional second parameter

    /// is the initial distance of the node.

    void addSource(Node source, Value dst = OperationTraits::zero()) {

      _dist->set(source, dst);

      if (!(*_mask)[source]) {

	_process.push_back(source);

	_mask->set(source, true);

      }

    }

    /// \brief Executes one round from the Bellman-Ford algorithm.

    ///

    /// If the algoritm calculated the distances in the previous round

    /// exactly for the paths of at most \c k arcs, then this function

    /// will calculate the distances exactly for the paths of at most

    /// <tt>k+1</tt> arcs. Performing \c k iterations using this function

    /// calculates the shortest path distances exactly for the paths

    /// consisting of at most \c k arcs.

    ///

    /// \warning The paths with limited arc number cannot be retrieved

    /// easily with \ref path() or \ref predArc() functions. If you also

    /// need the shortest paths and not only the distances, you should

    /// store the \ref predMap() "predecessor map" after each iteration

    /// and build the path manually.

    ///

    /// \return \c true when the algorithm have not found more shorter

    /// paths.

    ///

    /// \see ActiveIt

    bool processNextRound() {

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

	_mask->set(_process[i], false);

      }

      std::vector<Node> nextProcess;

      std::vector<Value> values(_process.size());

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

	values[i] = (*_dist)[_process[i]];

      }

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

	for (OutArcIt it(*_gr, _process[i]); it != INVALID; ++it) {

	  Node target = _gr->target(it);

	  Value relaxed = OperationTraits::plus(values[i], (*_length)[it]);

	  if (OperationTraits::less(relaxed, (*_dist)[target])) {

	    _pred->set(target, it);

	    _dist->set(target, relaxed);

	    if (!(*_mask)[target]) {

	      _mask->set(target, true);

	      nextProcess.push_back(target);

	    }

	  }	  

	}

      }

      _process.swap(nextProcess);

      return _process.empty();

    }

    /// \brief Executes one weak round from the Bellman-Ford algorithm.

    ///

    /// If the algorithm calculated the distances in the previous round

    /// at least for the paths of at most \c k arcs, then this function

    /// will calculate the distances at least for the paths of at most

    /// <tt>k+1</tt> arcs.

    /// This function does not make it possible to calculate the shortest

    /// path distances exactly for paths consisting of at most \c k arcs,

    /// this is why it is called weak round.

    ///

    /// \return \c true when the algorithm have not found more shorter

    /// paths.

    ///

    /// \see ActiveIt

    bool processNextWeakRound() {

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

	_mask->set(_process[i], false);

      }

      std::vector<Node> nextProcess;

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

	for (OutArcIt it(*_gr, _process[i]); it != INVALID; ++it) {

	  Node target = _gr->target(it);

	  Value relaxed = 

	    OperationTraits::plus((*_dist)[_process[i]], (*_length)[it]);

	  if (OperationTraits::less(relaxed, (*_dist)[target])) {

	    _pred->set(target, it);

	    _dist->set(target, relaxed);

	    if (!(*_mask)[target]) {

    ///

    /// This method runs the Bellman-Ford algorithm from the root node(s)

    /// in order to compute the shortest path to each node and also checks

    /// if the digraph contains cycles with negative total length.

    ///

    /// The algorithm computes 

    /// - the shortest path tree (forest),

    /// - the distance of each node from the root(s).

    /// 

    /// \return \c false if there is a negative cycle in the digraph.

    ///

    /// \pre init() must be called and at least one root node should be

    /// added with addSource() before using this function. 

    bool checkedStart() {

      int num = countNodes(*_gr);

      for (int i = 0; i < num; ++i) {

	if (processNextWeakRound()) return true;

      }

      return _process.empty();

    }

    /// \brief Executes the algorithm with arc number limit.

    ///

    /// Executes the algorithm with arc number limit.

    ///

    /// This method runs the Bellman-Ford algorithm from the root node(s)

    /// in order to compute the shortest path distance for each node

    /// using only the paths consisting of at most \c num arcs.

    ///

    /// The algorithm computes

    /// - the limited distance of each node from the root(s),

    /// - the predecessor arc for each node.

    ///

    /// \warning The paths with limited arc number cannot be retrieved

    /// easily with \ref path() or \ref predArc() functions. If you also

    /// need the shortest paths and not only the distances, you should

    /// store the \ref predMap() "predecessor map" after each iteration

    /// and build the path manually.

    ///

    /// \pre init() must be called and at least one root node should be

    /// added with addSource() before using this function. 

    void limitedStart(int num) {

      for (int i = 0; i < num; ++i) {

	if (processNextRound()) break;

      }

    }

    /// \brief Runs the algorithm from the given root node.

    ///    

    /// This method runs the Bellman-Ford algorithm from the given root

    /// node \c s in order to compute the shortest path to each node.

    ///

    /// The algorithm computes

    /// - the shortest path tree (forest),

    /// - the distance of each node from the root(s).

    ///

    /// \note bf.run(s) is just a shortcut of the following code.

    /// \code

    ///   bf.init();

    ///   bf.addSource(s);

    ///   bf.start();

    /// \endcode

    void run(Node s) {

      init();

      addSource(s);

      start();

    }

    /// \brief Runs the algorithm from the given root node with arc

    /// number limit.

    ///    

    /// This method runs the Bellman-Ford algorithm from the given root

    /// node \c s in order to compute the shortest path distance for each

    /// node using only the paths consisting of at most \c num arcs.

    ///

    /// The algorithm computes

    /// - the limited distance of each node from the root(s),

    /// - the predecessor arc for each node.

    ///

    /// \warning The paths with limited arc number cannot be retrieved

    /// easily with \ref path() or \ref predArc() functions. If you also

    /// need the shortest paths and not only the distances, you should

    /// store the \ref predMap() "predecessor map" after each iteration

    /// and build the path manually.

    ///

    /// \note bf.run(s, num) is just a shortcut of the following code.

    /// \code

    ///   bf.init();

    ///   bf.addSource(s);

    ///   bf.limitedStart(num);

    /// \endcode

    void run(Node s, int num) {

      init();

      addSource(s);

      limitedStart(num);

    }

    ///@}

    /// \brief LEMON iterator for getting the active nodes.

    ///

    /// This class provides a common style LEMON iterator that traverses

    /// the active nodes of the Bellman-Ford algorithm after the last

    /// phase. These nodes should be checked in the next phase to

    /// find augmenting arcs outgoing from them.

    class ActiveIt {

    public:

      /// \brief Constructor.

      ///

      /// Constructor for getting the active nodes of the given BellmanFord

      /// instance. 

      ActiveIt(const BellmanFord& algorithm) : _algorithm(&algorithm)

      {

        _index = _algorithm->_process.size() - 1;

      }

      /// \brief Invalid constructor.

      ///

      /// Invalid constructor.

      ActiveIt(Invalid) : _algorithm(0), _index(-1) {}

      /// \brief Conversion to \c Node.

      ///

      /// Conversion to \c Node.

      operator Node() const { 

        return _index >= 0 ? _algorithm->_process[_index] : INVALID;

      }

      /// \brief Increment operator.

      ///

      /// Increment operator.

      ActiveIt& operator++() {

        --_index;

        return *this; 

      }

      bool operator==(const ActiveIt& it) const { 

        return static_cast<Node>(*this) == static_cast<Node>(it); 

      }

      bool operator!=(const ActiveIt& it) const { 

        return static_cast<Node>(*this) != static_cast<Node>(it); 

      }

      bool operator<(const ActiveIt& it) const { 

        return static_cast<Node>(*this) < static_cast<Node>(it); 

      }

    private:

      const BellmanFord* _algorithm;

      int _index;

    };

    /// \name Query Functions

    /// The result of the Bellman-Ford algorithm can be obtained using these

    /// functions.\n

    /// Either \ref run() or \ref init() should be called before using them.

    ///@{

    /// \brief The shortest path to the given node.

    ///    

    /// Gives back the shortest path to the given node from the root(s).

    ///

    /// \warning \c t should be reached from the root(s).

    ///

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

    /// using this function.

    Path path(Node t) const

    {

      return Path(*_gr, *_pred, t);

    }

    /// \brief The distance of the given node from the root(s).

    ///

    /// Returns the distance of the given node from the root(s).

    ///

    /// \warning If node \c v is not reached from the root(s), then

    /// the return value of this function is undefined.

    ///

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

    /// using this function.

    Value dist(Node v) const { return (*_dist)[v]; }

    /// \brief Returns the 'previous arc' of the shortest path tree for

    /// the given node.

    ///

    /// This function returns the 'previous arc' of the shortest path

    /// tree for node \c v, i.e. it returns the last arc of a

    /// shortest path from a root to \c v. It is \c INVALID if \c v

    /// is not reached from the root(s) or if \c v is a root.

    ///

    /// The shortest path tree used here is equal to the shortest path

    ///

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

    /// using this function.

    Arc predArc(Node v) const { return (*_pred)[v]; }

    /// \brief Returns the 'previous node' of the shortest path tree for

    /// the given node.

    ///

    /// This function returns the 'previous node' of the shortest path

    /// tree for node \c v, i.e. it returns the last but one node of

    /// a shortest path from a root to \c v. It is \c INVALID if \c v

    /// is not reached from the root(s) or if \c v is a root.

    ///

    /// The shortest path tree used here is equal to the shortest path

    ///

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

    /// using this function.

    Node predNode(Node v) const { 

      return (*_pred)[v] == INVALID ? INVALID : _gr->source((*_pred)[v]); 

    }

    /// \brief Returns a const reference to the node map that stores the

    /// distances of the nodes.

    ///

    /// Returns a const reference to the node map that stores the distances

    /// of the nodes calculated by the algorithm.

    ///

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

    /// using this function.

    const DistMap &distMap() const { return *_dist;}

    /// \brief Returns a const reference to the node map that stores the

    /// predecessor arcs.

    ///

    /// Returns a const reference to the node map that stores the predecessor

    /// arcs, which form the shortest path tree (forest).

    ///

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

    /// using this function.

    const PredMap &predMap() const { return *_pred; }

    /// \brief Checks if a node is reached from the root(s).

    ///

    /// Returns \c true if \c v is reached from the root(s).

    ///

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

    /// using this function.

    bool reached(Node v) const {

      return (*_dist)[v] != OperationTraits::infinity();

    }

    /// \brief Gives back a negative cycle.

    ///    

    /// This function gives back a directed cycle with negative total

    /// length if the algorithm has already found one.

    /// Otherwise it gives back an empty path.

    lemon::Path<Digraph> negativeCycle() {

      typename Digraph::template NodeMap<int> state(*_gr, -1);

      lemon::Path<Digraph> cycle;

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

        if (state[_process[i]] != -1) continue;

        for (Node v = _process[i]; (*_pred)[v] != INVALID;

             v = _gr->source((*_pred)[v])) {

          if (state[v] == i) {

            cycle.addFront((*_pred)[v]);

            for (Node u = _gr->source((*_pred)[v]); u != v;

                 u = _gr->source((*_pred)[u])) {

              cycle.addFront((*_pred)[u]);

            }

            return cycle;

          }

          else if (state[v] >= 0) {

            break;

          }

          state[v] = i;

        }

      }

      return cycle;

    }

    ///@}

  };

  /// \brief Default traits class of bellmanFord() function.

  ///

  /// Default traits class of bellmanFord() function.

  /// \tparam GR The type of the digraph.

  /// \tparam LEN The type of the length map.

  template <typename GR, typename LEN>

  struct BellmanFordWizardDefaultTraits {

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

    typedef GR Digraph;

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

    ///

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

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

    typedef LEN LengthMap;

    /// The type of the arc lengths.

    typedef typename LEN::Value Value;

    /// \brief Operation traits for Bellman-Ford algorithm.

    ///

    /// It defines the used operations and the infinity value for the

    /// given \c Value type.

    /// \see BellmanFordDefaultOperationTraits

    typedef BellmanFordDefaultOperationTraits<Value> OperationTraits;

    /// \brief The type of the map that stores the last

    /// arcs of the shortest paths.

    /// 

    /// The type of the map that stores the last arcs of the shortest paths.

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

    typedef typename GR::template NodeMap<typename GR::Arc> PredMap;

    /// \brief Instantiates a \c PredMap.

    /// 

    /// This function instantiates a \ref PredMap.

    /// \param g is the digraph to which we would like to define the

    /// \ref PredMap.

    static PredMap *createPredMap(const GR &g) {

      return new PredMap(g);

    }

    /// \brief The type of the map that stores the distances of the nodes.

    ///

    /// The type of the map that stores the distances of the nodes.

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

    typedef typename GR::template NodeMap<Value> DistMap;

    /// \brief Instantiates a \c DistMap.

    ///

    /// This function instantiates a \ref DistMap. 

    /// \param g is the digraph to which we would like to define the

    /// \ref DistMap.

    static DistMap *createDistMap(const GR &g) {

      return new DistMap(g);

    }

    ///The type of the shortest paths.

    ///The type of the shortest paths.

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

    typedef lemon::Path<Digraph> Path;

  };

  /// \brief Default traits class used by BellmanFordWizard.

  ///

  /// Default traits class used by BellmanFordWizard.

  /// \tparam GR The type of the digraph.

  /// \tparam LEN The type of the length map.

  template <typename GR, typename LEN>

  class BellmanFordWizardBase 

    : public BellmanFordWizardDefaultTraits<GR, LEN> {

    typedef BellmanFordWizardDefaultTraits<GR, LEN> Base;

  protected:

    // Type of the nodes in the digraph.

    typedef typename Base::Digraph::Node Node;

    // Pointer to the underlying digraph.

    void *_graph;

    // Pointer to the length map

    void *_length;

    // Pointer to the map of predecessors arcs.

    void *_pred;

    // Pointer to the map of distances.

    void *_dist;

    //Pointer to the shortest path to the target node.

    void *_path;

    //Pointer to the distance of the target node.

    void *_di;

    public:

    /// Constructor.

    /// This constructor does not require parameters, it initiates

    /// all of the attributes to default values \c 0.

    BellmanFordWizardBase() :

      _graph(0), _length(0), _pred(0), _dist(0), _path(0), _di(0) {}

    /// Constructor.

    /// This constructor requires two parameters,

    /// others are initiated to \c 0.

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

    /// \param len The length map.

    BellmanFordWizardBase(const GR& gr, 

			  const LEN& len) :

      _graph(reinterpret_cast<void*>(const_cast<GR*>(&gr))), 

      _length(reinterpret_cast<void*>(const_cast<LEN*>(&len))), 

      _pred(0), _dist(0), _path(0), _di(0) {}

  };

  /// \brief Auxiliary class for the function-type interface of the

  /// \ref BellmanFord "Bellman-Ford" algorithm.

  ///

  /// This auxiliary class is created to implement the

  /// \ref bellmanFord() "function-type interface" of the

  /// \ref BellmanFord "Bellman-Ford" algorithm.

  /// It does not have own \ref run() method, it uses the

  /// functions and features of the plain \ref BellmanFord.

  ///

  /// This class should only be used through the \ref bellmanFord()

/* -*- mode: C++; indent-tabs-mode: nil; -*-

 *

 * This file is a part of LEMON, a generic C++ optimization library.

 *

 * Copyright (C) 2003-2009

 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport

 * (Egervary Research Group on Combinatorial Optimization, EGRES).

 *

 * Permission to use, modify and distribute this software is granted

 * provided that this copyright notice appears in all copies. For

 * precise terms see the accompanying LICENSE file.

 *

 * This software is provided "AS IS" with no warranty of any kind,

 * express or implied, and with no claim as to its suitability for any

 * purpose.

 *

 */

#ifndef LEMON_BFS_H

#define LEMON_BFS_H

///\ingroup search

///\file

///\brief BFS algorithm.

#include <lemon/list_graph.h>

#include <lemon/bits/path_dump.h>

#include <lemon/core.h>

#include <lemon/error.h>

#include <lemon/maps.h>

#include <lemon/path.h>

namespace lemon {

  ///Default traits class of Bfs class.

  ///Default traits class of Bfs class.

  ///\tparam GR Digraph type.

  template<class GR>

  struct BfsDefaultTraits

  {

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

    typedef GR Digraph;

    ///\brief The type of the map that stores the predecessor

    ///arcs of the shortest paths.

    ///

    ///The type of the map that stores the predecessor

    ///arcs of the shortest paths.

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

    typedef typename Digraph::template NodeMap<typename Digraph::Arc> PredMap;

    ///Instantiates a \c PredMap.

    ///This function instantiates a \ref PredMap.

    ///\param g is the digraph, to which we would like to define the

    ///\ref PredMap.

    static PredMap *createPredMap(const Digraph &g)

    {

      return new PredMap(g);

    }

    ///The type of the map that indicates which nodes are processed.

    ///The type of the map that indicates which nodes are processed.

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

    typedef NullMap<typename Digraph::Node,bool> ProcessedMap;

    ///Instantiates a \c ProcessedMap.

    ///This function instantiates a \ref ProcessedMap.

    ///\param g is the digraph, to which

    ///we would like to define the \ref ProcessedMap

#ifdef DOXYGEN

    static ProcessedMap *createProcessedMap(const Digraph &g)

#else

    static ProcessedMap *createProcessedMap(const Digraph &)

#endif

    {

      return new ProcessedMap();

    }

    ///The type of the map that indicates which nodes are reached.

    ///The type of the map that indicates which nodes are reached.

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

    typedef typename Digraph::template NodeMap<bool> ReachedMap;

    ///Instantiates a \c ReachedMap.

    ///This function instantiates a \ref ReachedMap.

    ///\param g is the digraph, to which

    ///we would like to define the \ref ReachedMap.

    static ReachedMap *createReachedMap(const Digraph &g)

    {

      return new ReachedMap(g);

    }

    ///The type of the map that stores the distances of the nodes.

    ///The type of the map that stores the distances of the nodes.

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

    typedef typename Digraph::template NodeMap<int> DistMap;

    ///Instantiates a \c DistMap.

    ///This function instantiates a \ref DistMap.

    ///\param g is the digraph, to which we would like to define the

    ///\ref DistMap.

    static DistMap *createDistMap(const Digraph &g)

    {

      return new DistMap(g);

    }

  };

  ///%BFS algorithm class.

  ///\ingroup search

  ///This class provides an efficient implementation of the %BFS algorithm.

  ///

  ///There is also a \ref bfs() "function-type interface" for the BFS

  ///algorithm, which is convenient in the simplier cases and it can be

  ///used easier.

  ///

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

  ///The default type is \ref ListDigraph.

#ifdef DOXYGEN

  template <typename GR,

            typename TR>

#else

  template <typename GR=ListDigraph,

            typename TR=BfsDefaultTraits<GR> >

#endif

  class Bfs {

  public:

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

    typedef typename TR::Digraph Digraph;

    ///\brief The type of the map that stores the predecessor arcs of the

    ///shortest paths.

    typedef typename TR::PredMap PredMap;

    ///The type of the map that stores the distances of the nodes.

    typedef typename TR::DistMap DistMap;

    ///The type of the map that indicates which nodes are reached.

    typedef typename TR::ReachedMap ReachedMap;

    ///The type of the map that indicates which nodes are processed.