alpar@414: /* -*- mode: C++; indent-tabs-mode: nil; -*-
alpar@414:  *
alpar@414:  * This file is a part of LEMON, a generic C++ optimization library.
alpar@414:  *
alpar@414:  * Copyright (C) 2003-2008
alpar@414:  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
alpar@414:  * (Egervary Research Group on Combinatorial Optimization, EGRES).
alpar@414:  *
alpar@414:  * Permission to use, modify and distribute this software is granted
alpar@414:  * provided that this copyright notice appears in all copies. For
alpar@414:  * precise terms see the accompanying LICENSE file.
alpar@414:  *
alpar@414:  * This software is provided "AS IS" with no warranty of any kind,
alpar@414:  * express or implied, and with no claim as to its suitability for any
alpar@414:  * purpose.
alpar@414:  *
alpar@414:  */
alpar@414: 
alpar@414: #ifndef LEMON_CIRCULATION_H
alpar@414: #define LEMON_CIRCULATION_H
alpar@414: 
alpar@414: #include <iostream>
alpar@414: #include <queue>
alpar@414: #include <lemon/tolerance.h>
alpar@414: #include <lemon/elevator.h>
alpar@414: 
alpar@414: ///\ingroup max_flow
alpar@414: ///\file
alpar@414: ///\brief Push-prelabel algorithm for finding a feasible circulation.
alpar@414: ///
alpar@414: namespace lemon {
alpar@414: 
alpar@414:   /// \brief Default traits class of Circulation class.
alpar@414:   ///
alpar@414:   /// Default traits class of Circulation class.
alpar@414:   /// \param _Graph Digraph type.
alpar@414:   /// \param _CapacityMap Type of capacity map.
alpar@414:   template <typename _Graph, typename _LCapMap,
alpar@414:             typename _UCapMap, typename _DeltaMap>
alpar@414:   struct CirculationDefaultTraits {
alpar@414: 
alpar@414:     /// \brief The digraph type the algorithm runs on.
alpar@414:     typedef _Graph Digraph;
alpar@414: 
alpar@414:     /// \brief The type of the map that stores the circulation lower
alpar@414:     /// bound.
alpar@414:     ///
alpar@414:     /// The type of the map that stores the circulation lower bound.
alpar@414:     /// It must meet the \ref concepts::ReadMap "ReadMap" concept.
alpar@414:     typedef _LCapMap LCapMap;
alpar@414: 
alpar@414:     /// \brief The type of the map that stores the circulation upper
alpar@414:     /// bound.
alpar@414:     ///
alpar@414:     /// The type of the map that stores the circulation upper bound.
alpar@414:     /// It must meet the \ref concepts::ReadMap "ReadMap" concept.
alpar@414:     typedef _UCapMap UCapMap;
alpar@414: 
alpar@414:     /// \brief The type of the map that stores the upper bound of
alpar@414:     /// node excess.
alpar@414:     ///
alpar@414:     /// The type of the map that stores the lower bound of node
alpar@414:     /// excess. It must meet the \ref concepts::ReadMap "ReadMap"
alpar@414:     /// concept.
alpar@414:     typedef _DeltaMap DeltaMap;
alpar@414: 
alpar@414:     /// \brief The type of the length of the arcs.
alpar@414:     typedef typename DeltaMap::Value Value;
alpar@414: 
alpar@414:     /// \brief The map type that stores the flow values.
alpar@414:     ///
alpar@414:     /// The map type that stores the flow values.
alpar@414:     /// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept.
alpar@414:     typedef typename Digraph::template ArcMap<Value> FlowMap;
alpar@414: 
alpar@414:     /// \brief Instantiates a FlowMap.
alpar@414:     ///
alpar@414:     /// This function instantiates a \ref FlowMap.
alpar@414:     /// \param digraph The digraph, to which we would like to define
alpar@414:     /// the flow map.
alpar@414:     static FlowMap* createFlowMap(const Digraph& digraph) {
alpar@414:       return new FlowMap(digraph);
alpar@414:     }
alpar@414: 
alpar@414:     /// \brief The eleavator type used by Circulation algorithm.
alpar@414:     ///
alpar@414:     /// The elevator type used by Circulation algorithm.
alpar@414:     ///
alpar@414:     /// \sa Elevator
alpar@414:     /// \sa LinkedElevator
alpar@414:     typedef lemon::Elevator<Digraph, typename Digraph::Node> Elevator;
alpar@414: 
alpar@414:     /// \brief Instantiates an Elevator.
alpar@414:     ///
alpar@414:     /// This function instantiates a \ref Elevator.
alpar@414:     /// \param digraph The digraph, to which we would like to define
alpar@414:     /// the elevator.
alpar@414:     /// \param max_level The maximum level of the elevator.
alpar@414:     static Elevator* createElevator(const Digraph& digraph, int max_level) {
alpar@414:       return new Elevator(digraph, max_level);
alpar@414:     }
alpar@414: 
alpar@414:     /// \brief The tolerance used by the algorithm
alpar@414:     ///
alpar@414:     /// The tolerance used by the algorithm to handle inexact computation.
alpar@414:     typedef lemon::Tolerance<Value> Tolerance;
alpar@414: 
alpar@414:   };
alpar@414: 
alpar@414:   ///Push-relabel algorithm for the Network Circulation Problem.
alpar@414: 
alpar@414:   /**
alpar@414:      \ingroup max_flow
alpar@414:      This class implements a push-relabel algorithm
alpar@414:      or the Network Circulation Problem.
alpar@414:      The exact formulation of this problem is the following.
alpar@414:      \f[\sum_{e\in\rho(v)}x(e)-\sum_{e\in\delta(v)}x(e)\leq
alpar@414:      -delta(v)\quad \forall v\in V \f]
alpar@414:      \f[ lo(e)\leq x(e) \leq up(e) \quad \forall e\in E \f]
alpar@414:   */
alpar@414:   template<class _Graph,
alpar@414:            class _LCapMap=typename _Graph::template ArcMap<int>,
alpar@414:            class _UCapMap=_LCapMap,
alpar@414:            class _DeltaMap=typename _Graph::template NodeMap<
alpar@414:              typename _UCapMap::Value>,
alpar@414:            class _Traits=CirculationDefaultTraits<_Graph, _LCapMap,
alpar@414:                                                   _UCapMap, _DeltaMap> >
alpar@414:   class Circulation {
alpar@414: 
alpar@414:     typedef _Traits Traits;
alpar@414:     typedef typename Traits::Digraph Digraph;
alpar@414:     TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
alpar@414: 
alpar@414:     typedef typename Traits::Value Value;
alpar@414: 
alpar@414:     typedef typename Traits::LCapMap LCapMap;
alpar@414:     typedef typename Traits::UCapMap UCapMap;
alpar@414:     typedef typename Traits::DeltaMap DeltaMap;
alpar@414:     typedef typename Traits::FlowMap FlowMap;
alpar@414:     typedef typename Traits::Elevator Elevator;
alpar@414:     typedef typename Traits::Tolerance Tolerance;
alpar@414: 
alpar@414:     typedef typename Digraph::template NodeMap<Value> ExcessMap;
alpar@414: 
alpar@414:     const Digraph &_g;
alpar@414:     int _node_num;
alpar@414: 
alpar@414:     const LCapMap *_lo;
alpar@414:     const UCapMap *_up;
alpar@414:     const DeltaMap *_delta;
alpar@414: 
alpar@414:     FlowMap *_flow;
alpar@414:     bool _local_flow;
alpar@414: 
alpar@414:     Elevator* _level;
alpar@414:     bool _local_level;
alpar@414: 
alpar@414:     ExcessMap* _excess;
alpar@414: 
alpar@414:     Tolerance _tol;
alpar@414:     int _el;
alpar@414: 
alpar@414:   public:
alpar@414: 
alpar@414:     typedef Circulation Create;
alpar@414: 
alpar@414:     ///\name Named template parameters
alpar@414: 
alpar@414:     ///@{
alpar@414: 
alpar@414:     template <typename _FlowMap>
alpar@416:     struct SetFlowMapTraits : public Traits {
alpar@414:       typedef _FlowMap FlowMap;
alpar@414:       static FlowMap *createFlowMap(const Digraph&) {
alpar@414:         LEMON_ASSERT(false, "FlowMap is not initialized");
alpar@414:         return 0; // ignore warnings
alpar@414:       }
alpar@414:     };
alpar@414: 
alpar@414:     /// \brief \ref named-templ-param "Named parameter" for setting
alpar@414:     /// FlowMap type
alpar@414:     ///
alpar@414:     /// \ref named-templ-param "Named parameter" for setting FlowMap
alpar@414:     /// type
alpar@414:     template <typename _FlowMap>
alpar@416:     struct SetFlowMap
alpar@414:       : public Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
alpar@416:                            SetFlowMapTraits<_FlowMap> > {
alpar@414:       typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
alpar@416:                           SetFlowMapTraits<_FlowMap> > Create;
alpar@414:     };
alpar@414: 
alpar@414:     template <typename _Elevator>
alpar@416:     struct SetElevatorTraits : public Traits {
alpar@414:       typedef _Elevator Elevator;
alpar@414:       static Elevator *createElevator(const Digraph&, int) {
alpar@414:         LEMON_ASSERT(false, "Elevator is not initialized");
alpar@414:         return 0; // ignore warnings
alpar@414:       }
alpar@414:     };
alpar@414: 
alpar@414:     /// \brief \ref named-templ-param "Named parameter" for setting
alpar@414:     /// Elevator type
alpar@414:     ///
alpar@414:     /// \ref named-templ-param "Named parameter" for setting Elevator
alpar@414:     /// type
alpar@414:     template <typename _Elevator>
alpar@416:     struct SetElevator
alpar@414:       : public Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
alpar@416:                            SetElevatorTraits<_Elevator> > {
alpar@414:       typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
alpar@416:                           SetElevatorTraits<_Elevator> > Create;
alpar@414:     };
alpar@414: 
alpar@414:     template <typename _Elevator>
alpar@416:     struct SetStandardElevatorTraits : public Traits {
alpar@414:       typedef _Elevator Elevator;
alpar@414:       static Elevator *createElevator(const Digraph& digraph, int max_level) {
alpar@414:         return new Elevator(digraph, max_level);
alpar@414:       }
alpar@414:     };
alpar@414: 
alpar@414:     /// \brief \ref named-templ-param "Named parameter" for setting
alpar@414:     /// Elevator type
alpar@414:     ///
alpar@414:     /// \ref named-templ-param "Named parameter" for setting Elevator
alpar@414:     /// type. The Elevator should be standard constructor interface, ie.
alpar@414:     /// the digraph and the maximum level should be passed to it.
alpar@414:     template <typename _Elevator>
alpar@416:     struct SetStandardElevator
alpar@414:       : public Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
alpar@416:                        SetStandardElevatorTraits<_Elevator> > {
alpar@414:       typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
alpar@416:                       SetStandardElevatorTraits<_Elevator> > Create;
alpar@414:     };
alpar@414: 
alpar@414:     /// @}
alpar@414: 
alpar@414:   protected:
alpar@414: 
alpar@414:     Circulation() {}
alpar@414: 
alpar@414:   public:
alpar@414: 
alpar@414:     /// The constructor of the class.
alpar@414: 
alpar@414:     /// The constructor of the class.
alpar@414:     /// \param g The digraph the algorithm runs on.
alpar@414:     /// \param lo The lower bound capacity of the arcs.
alpar@414:     /// \param up The upper bound capacity of the arcs.
alpar@414:     /// \param delta The lower bound on node excess.
alpar@414:     Circulation(const Digraph &g,const LCapMap &lo,
alpar@414:                 const UCapMap &up,const DeltaMap &delta)
alpar@414:       : _g(g), _node_num(),
alpar@414:         _lo(&lo),_up(&up),_delta(&delta),_flow(0),_local_flow(false),
alpar@414:         _level(0), _local_level(false), _excess(0), _el() {}
alpar@414: 
alpar@414:     /// Destrcutor.
alpar@414:     ~Circulation() {
alpar@414:       destroyStructures();
alpar@414:     }
alpar@414: 
alpar@414:   private:
alpar@414: 
alpar@414:     void createStructures() {
alpar@414:       _node_num = _el = countNodes(_g);
alpar@414: 
alpar@414:       if (!_flow) {
alpar@414:         _flow = Traits::createFlowMap(_g);
alpar@414:         _local_flow = true;
alpar@414:       }
alpar@414:       if (!_level) {
alpar@414:         _level = Traits::createElevator(_g, _node_num);
alpar@414:         _local_level = true;
alpar@414:       }
alpar@414:       if (!_excess) {
alpar@414:         _excess = new ExcessMap(_g);
alpar@414:       }
alpar@414:     }
alpar@414: 
alpar@414:     void destroyStructures() {
alpar@414:       if (_local_flow) {
alpar@414:         delete _flow;
alpar@414:       }
alpar@414:       if (_local_level) {
alpar@414:         delete _level;
alpar@414:       }
alpar@414:       if (_excess) {
alpar@414:         delete _excess;
alpar@414:       }
alpar@414:     }
alpar@414: 
alpar@414:   public:
alpar@414: 
alpar@414:     /// Sets the lower bound capacity map.
alpar@414: 
alpar@414:     /// Sets the lower bound capacity map.
alpar@414:     /// \return \c (*this)
alpar@414:     Circulation& lowerCapMap(const LCapMap& map) {
alpar@414:       _lo = &map;
alpar@414:       return *this;
alpar@414:     }
alpar@414: 
alpar@414:     /// Sets the upper bound capacity map.
alpar@414: 
alpar@414:     /// Sets the upper bound capacity map.
alpar@414:     /// \return \c (*this)
alpar@414:     Circulation& upperCapMap(const LCapMap& map) {
alpar@414:       _up = &map;
alpar@414:       return *this;
alpar@414:     }
alpar@414: 
alpar@414:     /// Sets the lower bound map on excess.
alpar@414: 
alpar@414:     /// Sets the lower bound map on excess.
alpar@414:     /// \return \c (*this)
alpar@414:     Circulation& deltaMap(const DeltaMap& map) {
alpar@414:       _delta = &map;
alpar@414:       return *this;
alpar@414:     }
alpar@414: 
alpar@414:     /// Sets the flow map.
alpar@414: 
alpar@414:     /// Sets the flow map.
alpar@414:     /// \return \c (*this)
alpar@414:     Circulation& flowMap(FlowMap& map) {
alpar@414:       if (_local_flow) {
alpar@414:         delete _flow;
alpar@414:         _local_flow = false;
alpar@414:       }
alpar@414:       _flow = &map;
alpar@414:       return *this;
alpar@414:     }
alpar@414: 
alpar@414:     /// Returns the flow map.
alpar@414: 
alpar@414:     /// \return The flow map.
alpar@414:     ///
alpar@414:     const FlowMap& flowMap() {
alpar@414:       return *_flow;
alpar@414:     }
alpar@414: 
alpar@414:     /// Sets the elevator.
alpar@414: 
alpar@414:     /// Sets the elevator.
alpar@414:     /// \return \c (*this)
alpar@414:     Circulation& elevator(Elevator& elevator) {
alpar@414:       if (_local_level) {
alpar@414:         delete _level;
alpar@414:         _local_level = false;
alpar@414:       }
alpar@414:       _level = &elevator;
alpar@414:       return *this;
alpar@414:     }
alpar@414: 
alpar@414:     /// Returns the elevator.
alpar@414: 
alpar@414:     /// \return The elevator.
alpar@414:     ///
alpar@414:     const Elevator& elevator() {
alpar@414:       return *_level;
alpar@414:     }
alpar@414: 
alpar@414:     /// Sets the tolerance used by algorithm.
alpar@414: 
alpar@414:     /// Sets the tolerance used by algorithm.
alpar@414:     ///
alpar@414:     Circulation& tolerance(const Tolerance& tolerance) const {
alpar@414:       _tol = tolerance;
alpar@414:       return *this;
alpar@414:     }
alpar@414: 
alpar@414:     /// Returns the tolerance used by algorithm.
alpar@414: 
alpar@414:     /// Returns the tolerance used by algorithm.
alpar@414:     ///
alpar@414:     const Tolerance& tolerance() const {
alpar@414:       return tolerance;
alpar@414:     }
alpar@414: 
alpar@414:     /// \name Execution control
alpar@414:     /// The simplest way to execute the algorithm is to use one of the
alpar@414:     /// member functions called \c run().
alpar@414:     /// \n
alpar@414:     /// If you need more control on initial solution or execution then
alpar@414:     /// you have to call one \ref init() function and then the start()
alpar@414:     /// function.
alpar@414: 
alpar@414:     ///@{
alpar@414: 
alpar@414:     /// Initializes the internal data structures.
alpar@414: 
alpar@414:     /// Initializes the internal data structures. This function sets
alpar@414:     /// all flow values to the lower bound.
alpar@414:     /// \return This function returns false if the initialization
alpar@414:     /// process found a barrier.
alpar@414:     void init()
alpar@414:     {
alpar@414:       createStructures();
alpar@414: 
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n) {
alpar@414:         _excess->set(n, (*_delta)[n]);
alpar@414:       }
alpar@414: 
alpar@414:       for (ArcIt e(_g);e!=INVALID;++e) {
alpar@414:         _flow->set(e, (*_lo)[e]);
alpar@414:         _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_flow)[e]);
alpar@414:         _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_flow)[e]);
alpar@414:       }
alpar@414: 
alpar@414:       // global relabeling tested, but in general case it provides
alpar@414:       // worse performance for random digraphs
alpar@414:       _level->initStart();
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n)
alpar@414:         _level->initAddItem(n);
alpar@414:       _level->initFinish();
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n)
alpar@414:         if(_tol.positive((*_excess)[n]))
alpar@414:           _level->activate(n);
alpar@414:     }
alpar@414: 
alpar@414:     /// Initializes the internal data structures.
alpar@414: 
alpar@414:     /// Initializes the internal data structures. This functions uses
alpar@414:     /// greedy approach to construct the initial solution.
alpar@414:     void greedyInit()
alpar@414:     {
alpar@414:       createStructures();
alpar@414: 
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n) {
alpar@414:         _excess->set(n, (*_delta)[n]);
alpar@414:       }
alpar@414: 
alpar@414:       for (ArcIt e(_g);e!=INVALID;++e) {
alpar@414:         if (!_tol.positive((*_excess)[_g.target(e)] + (*_up)[e])) {
alpar@414:           _flow->set(e, (*_up)[e]);
alpar@414:           _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_up)[e]);
alpar@414:           _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_up)[e]);
alpar@414:         } else if (_tol.positive((*_excess)[_g.target(e)] + (*_lo)[e])) {
alpar@414:           _flow->set(e, (*_lo)[e]);
alpar@414:           _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_lo)[e]);
alpar@414:           _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_lo)[e]);
alpar@414:         } else {
alpar@414:           Value fc = -(*_excess)[_g.target(e)];
alpar@414:           _flow->set(e, fc);
alpar@414:           _excess->set(_g.target(e), 0);
alpar@414:           _excess->set(_g.source(e), (*_excess)[_g.source(e)] - fc);
alpar@414:         }
alpar@414:       }
alpar@414: 
alpar@414:       _level->initStart();
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n)
alpar@414:         _level->initAddItem(n);
alpar@414:       _level->initFinish();
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n)
alpar@414:         if(_tol.positive((*_excess)[n]))
alpar@414:           _level->activate(n);
alpar@414:     }
alpar@414: 
alpar@414:     ///Starts the algorithm
alpar@414: 
alpar@414:     ///This function starts the algorithm.
alpar@414:     ///\return This function returns true if it found a feasible circulation.
alpar@414:     ///
alpar@414:     ///\sa barrier()
alpar@414:     bool start()
alpar@414:     {
alpar@414: 
alpar@414:       Node act;
alpar@414:       Node bact=INVALID;
alpar@414:       Node last_activated=INVALID;
alpar@414:       while((act=_level->highestActive())!=INVALID) {
alpar@414:         int actlevel=(*_level)[act];
alpar@414:         int mlevel=_node_num;
alpar@414:         Value exc=(*_excess)[act];
alpar@414: 
alpar@414:         for(OutArcIt e(_g,act);e!=INVALID; ++e) {
alpar@414:           Node v = _g.target(e);
alpar@414:           Value fc=(*_up)[e]-(*_flow)[e];
alpar@414:           if(!_tol.positive(fc)) continue;
alpar@414:           if((*_level)[v]<actlevel) {
alpar@414:             if(!_tol.less(fc, exc)) {
alpar@414:               _flow->set(e, (*_flow)[e] + exc);
alpar@414:               _excess->set(v, (*_excess)[v] + exc);
alpar@414:               if(!_level->active(v) && _tol.positive((*_excess)[v]))
alpar@414:                 _level->activate(v);
alpar@414:               _excess->set(act,0);
alpar@414:               _level->deactivate(act);
alpar@414:               goto next_l;
alpar@414:             }
alpar@414:             else {
alpar@414:               _flow->set(e, (*_up)[e]);
alpar@414:               _excess->set(v, (*_excess)[v] + fc);
alpar@414:               if(!_level->active(v) && _tol.positive((*_excess)[v]))
alpar@414:                 _level->activate(v);
alpar@414:               exc-=fc;
alpar@414:             }
alpar@414:           }
alpar@414:           else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
alpar@414:         }
alpar@414:         for(InArcIt e(_g,act);e!=INVALID; ++e) {
alpar@414:           Node v = _g.source(e);
alpar@414:           Value fc=(*_flow)[e]-(*_lo)[e];
alpar@414:           if(!_tol.positive(fc)) continue;
alpar@414:           if((*_level)[v]<actlevel) {
alpar@414:             if(!_tol.less(fc, exc)) {
alpar@414:               _flow->set(e, (*_flow)[e] - exc);
alpar@414:               _excess->set(v, (*_excess)[v] + exc);
alpar@414:               if(!_level->active(v) && _tol.positive((*_excess)[v]))
alpar@414:                 _level->activate(v);
alpar@414:               _excess->set(act,0);
alpar@414:               _level->deactivate(act);
alpar@414:               goto next_l;
alpar@414:             }
alpar@414:             else {
alpar@414:               _flow->set(e, (*_lo)[e]);
alpar@414:               _excess->set(v, (*_excess)[v] + fc);
alpar@414:               if(!_level->active(v) && _tol.positive((*_excess)[v]))
alpar@414:                 _level->activate(v);
alpar@414:               exc-=fc;
alpar@414:             }
alpar@414:           }
alpar@414:           else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
alpar@414:         }
alpar@414: 
alpar@414:         _excess->set(act, exc);
alpar@414:         if(!_tol.positive(exc)) _level->deactivate(act);
alpar@414:         else if(mlevel==_node_num) {
alpar@414:           _level->liftHighestActiveToTop();
alpar@414:           _el = _node_num;
alpar@414:           return false;
alpar@414:         }
alpar@414:         else {
alpar@414:           _level->liftHighestActive(mlevel+1);
alpar@414:           if(_level->onLevel(actlevel)==0) {
alpar@414:             _el = actlevel;
alpar@414:             return false;
alpar@414:           }
alpar@414:         }
alpar@414:       next_l:
alpar@414:         ;
alpar@414:       }
alpar@414:       return true;
alpar@414:     }
alpar@414: 
alpar@414:     /// Runs the circulation algorithm.
alpar@414: 
alpar@414:     /// Runs the circulation algorithm.
alpar@414:     /// \note fc.run() is just a shortcut of the following code.
alpar@414:     /// \code
alpar@414:     ///   fc.greedyInit();
alpar@414:     ///   return fc.start();
alpar@414:     /// \endcode
alpar@414:     bool run() {
alpar@414:       greedyInit();
alpar@414:       return start();
alpar@414:     }
alpar@414: 
alpar@414:     /// @}
alpar@414: 
alpar@414:     /// \name Query Functions
alpar@414:     /// The result of the %Circulation algorithm can be obtained using
alpar@414:     /// these functions.
alpar@414:     /// \n
alpar@414:     /// Before the use of these functions,
alpar@414:     /// either run() or start() must be called.
alpar@414: 
alpar@414:     ///@{
alpar@414: 
alpar@414:     /**
alpar@414:        \brief Returns a barrier
alpar@414:        
alpar@414:        Barrier is a set \e B of nodes for which
alpar@414:        \f[ \sum_{v\in B}-delta(v)<
alpar@414:        \sum_{e\in\rho(B)}lo(e)-\sum_{e\in\delta(B)}up(e) \f]
alpar@414:        holds. The existence of a set with this property prooves that a feasible
alpar@414:        flow cannot exists.
alpar@414:        \sa checkBarrier()
alpar@414:        \sa run()
alpar@414:     */
alpar@414:     template<class GT>
alpar@414:     void barrierMap(GT &bar)
alpar@414:     {
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n)
alpar@414:         bar.set(n, (*_level)[n] >= _el);
alpar@414:     }
alpar@414: 
alpar@414:     ///Returns true if the node is in the barrier
alpar@414: 
alpar@414:     ///Returns true if the node is in the barrier
alpar@414:     ///\sa barrierMap()
alpar@414:     bool barrier(const Node& node)
alpar@414:     {
alpar@414:       return (*_level)[node] >= _el;
alpar@414:     }
alpar@414: 
alpar@414:     /// \brief Returns the flow on the arc.
alpar@414:     ///
alpar@414:     /// Sets the \c flowMap to the flow on the arcs. This method can
alpar@414:     /// be called after the second phase of algorithm.
alpar@414:     Value flow(const Arc& arc) const {
alpar@414:       return (*_flow)[arc];
alpar@414:     }
alpar@414: 
alpar@414:     /// @}
alpar@414: 
alpar@414:     /// \name Checker Functions
alpar@414:     /// The feasibility  of the results can be checked using
alpar@414:     /// these functions.
alpar@414:     /// \n
alpar@414:     /// Before the use of these functions,
alpar@414:     /// either run() or start() must be called.
alpar@414: 
alpar@414:     ///@{
alpar@414: 
alpar@414:     ///Check if the  \c flow is a feasible circulation
alpar@414:     bool checkFlow() {
alpar@414:       for(ArcIt e(_g);e!=INVALID;++e)
alpar@414:         if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false;
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n)
alpar@414:         {
alpar@414:           Value dif=-(*_delta)[n];
alpar@414:           for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e];
alpar@414:           for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e];
alpar@414:           if(_tol.negative(dif)) return false;
alpar@414:         }
alpar@414:       return true;
alpar@414:     }
alpar@414: 
alpar@414:     ///Check whether or not the last execution provides a barrier
alpar@414: 
alpar@414:     ///Check whether or not the last execution provides a barrier
alpar@414:     ///\sa barrier()
alpar@414:     bool checkBarrier()
alpar@414:     {
alpar@414:       Value delta=0;
alpar@414:       for(NodeIt n(_g);n!=INVALID;++n)
alpar@414:         if(barrier(n))
alpar@414:           delta-=(*_delta)[n];
alpar@414:       for(ArcIt e(_g);e!=INVALID;++e)
alpar@414:         {
alpar@414:           Node s=_g.source(e);
alpar@414:           Node t=_g.target(e);
alpar@414:           if(barrier(s)&&!barrier(t)) delta+=(*_up)[e];
alpar@414:           else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e];
alpar@414:         }
alpar@414:       return _tol.negative(delta);
alpar@414:     }
alpar@414: 
alpar@414:     /// @}
alpar@414: 
alpar@414:   };
alpar@414: 
alpar@414: }
alpar@414: 
alpar@414: #endif