RhodeCode

 *

 * 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_CIRCULATION_H

#define LEMON_CIRCULATION_H

#include <lemon/tolerance.h>

#include <lemon/elevator.h>

///\ingroup max_flow

///\file

///\brief Push-relabel algorithm for finding a feasible circulation.

///

namespace lemon {

  /// \brief Default traits class of Circulation class.

  ///

  /// Default traits class of Circulation class.

  ///

  /// \tparam GR Type of the digraph the algorithm runs on.

  /**

     \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

     at the nodes.

     The exact formulation of this problem is the following.

     holds for all \f$uv\in A\f$, 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.

     solution of the following problem.

     \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.

     If you need the opposite inequalities in the supply/demand constraints

     (i.e. the total demand is less than the total supply and all the demands

     have to be satisfied while there could be supplies that are not used),

     then you could easily transform the problem to the above form by reversing

     the direction of the arcs and taking the negative of the supply values

     (e.g. using \ref ReverseDigraph and \ref NegMap adaptors).

     Note that this algorithm also provides a feasible solution for the

     \ref min_cost_flow "minimum cost flow problem".

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

     \tparam LM The type of the lower bound map. The default

     map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".

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

     The default map type is \c LM.

     \tparam SM The type of the supply map. The default map type is

     \ref concepts::Digraph::NodeMap "GR::NodeMap<UM::Value>".

  */

#ifdef DOXYGEN

      : _g(graph), _lo(&lower), _up(&upper), _supply(&supply),

        _flow(NULL), _local_flow(false), _level(NULL), _local_level(false),

        _excess(NULL) {}

    /// Destructor.

    ~Circulation() {

      destroyStructures();

    }

  private:

    void createStructures() {

      _node_num = _el = countNodes(_g);

      if (!_flow) {

        _flow = Traits::createFlowMap(_g);

        _local_flow = true;

      }

      if (!_level) {

        _level = Traits::createElevator(_g, _node_num);

        _local_level = true;

      }

      if (!_excess) {

    /// Sets the lower bound map.

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

    Circulation& lowerMap(const LowerMap& map) {

      _lo = ↦

      return *this;

    }

    /// Sets the upper bound (capacity) map.

    /// Sets the upper bound (capacity) map.

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

      _up = ↦

      return *this;

    }

    /// Sets the supply map.

    /// Sets the supply map.

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

    Circulation& supplyMap(const SupplyMap& map) {

      _supply = ↦

      return *this;

    }

    /// If you need more control on the initial solution or the execution,

    /// first you have to call one of the \ref init() functions, then

    /// the \ref start() function.

    ///@{

    /// Initializes the internal data structures.

    /// Initializes the internal data structures and sets all flow values

    /// to the lower bound.

    void init()

    {

      createStructures();

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

        (*_excess)[n] = (*_supply)[n];

      }

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

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

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

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

      }

      _level->initFinish();

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

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

          _level->activate(n);

    }

    /// Initializes the internal data structures using a greedy approach.

    /// Initializes the internal data structures using a greedy approach

    /// to construct the initial solution.

    void greedyInit()

    {

      createStructures();

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

        (*_excess)[n] = (*_supply)[n];

      }

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

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

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

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

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

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

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

        } else {

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

          _flow->set(e, fc);

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

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

        }

      }

      _level->initStart();

        }

      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

    {

      Flow delta=0;

      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);

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

        }

      return _tol.negative(delta);

    }

    /// @}

  };

}

#endif