alpar@399: /* -*- mode: C++; indent-tabs-mode: nil; -*- alpar@399: * alpar@399: * This file is a part of LEMON, a generic C++ optimization library. alpar@399: * alpar@399: * Copyright (C) 2003-2008 alpar@399: * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport alpar@399: * (Egervary Research Group on Combinatorial Optimization, EGRES). alpar@399: * alpar@399: * Permission to use, modify and distribute this software is granted alpar@399: * provided that this copyright notice appears in all copies. For alpar@399: * precise terms see the accompanying LICENSE file. alpar@399: * alpar@399: * This software is provided "AS IS" with no warranty of any kind, alpar@399: * express or implied, and with no claim as to its suitability for any alpar@399: * purpose. alpar@399: * alpar@399: */ alpar@399: alpar@399: #ifndef LEMON_CIRCULATION_H alpar@399: #define LEMON_CIRCULATION_H alpar@399: alpar@399: #include alpar@399: #include alpar@399: alpar@399: ///\ingroup max_flow alpar@399: ///\file kpeter@402: ///\brief Push-relabel algorithm for finding a feasible circulation. alpar@399: /// alpar@399: namespace lemon { alpar@399: alpar@399: /// \brief Default traits class of Circulation class. alpar@399: /// alpar@399: /// Default traits class of Circulation class. kpeter@402: /// \tparam _Diraph Digraph type. kpeter@402: /// \tparam _LCapMap Lower bound capacity map type. kpeter@402: /// \tparam _UCapMap Upper bound capacity map type. kpeter@402: /// \tparam _DeltaMap Delta map type. kpeter@402: template alpar@399: struct CirculationDefaultTraits { alpar@399: kpeter@402: /// \brief The type of the digraph the algorithm runs on. kpeter@402: typedef _Diraph Digraph; alpar@399: alpar@399: /// \brief The type of the map that stores the circulation lower alpar@399: /// bound. alpar@399: /// alpar@399: /// The type of the map that stores the circulation lower bound. alpar@399: /// It must meet the \ref concepts::ReadMap "ReadMap" concept. alpar@399: typedef _LCapMap LCapMap; alpar@399: alpar@399: /// \brief The type of the map that stores the circulation upper alpar@399: /// bound. alpar@399: /// alpar@399: /// The type of the map that stores the circulation upper bound. alpar@399: /// It must meet the \ref concepts::ReadMap "ReadMap" concept. alpar@399: typedef _UCapMap UCapMap; alpar@399: kpeter@402: /// \brief The type of the map that stores the lower bound for kpeter@402: /// the supply of the nodes. alpar@399: /// kpeter@402: /// The type of the map that stores the lower bound for the supply kpeter@402: /// of the nodes. It must meet the \ref concepts::ReadMap "ReadMap" alpar@399: /// concept. alpar@399: typedef _DeltaMap DeltaMap; alpar@399: kpeter@402: /// \brief The type of the flow values. alpar@399: typedef typename DeltaMap::Value Value; alpar@399: kpeter@402: /// \brief The type of the map that stores the flow values. alpar@399: /// kpeter@402: /// The type of the map that stores the flow values. alpar@399: /// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept. alpar@399: typedef typename Digraph::template ArcMap FlowMap; alpar@399: alpar@399: /// \brief Instantiates a FlowMap. alpar@399: /// alpar@399: /// This function instantiates a \ref FlowMap. alpar@399: /// \param digraph The digraph, to which we would like to define alpar@399: /// the flow map. alpar@399: static FlowMap* createFlowMap(const Digraph& digraph) { alpar@399: return new FlowMap(digraph); alpar@399: } alpar@399: kpeter@402: /// \brief The elevator type used by the algorithm. alpar@399: /// kpeter@402: /// The elevator type used by the algorithm. alpar@399: /// alpar@399: /// \sa Elevator alpar@399: /// \sa LinkedElevator alpar@399: typedef lemon::Elevator Elevator; alpar@399: alpar@399: /// \brief Instantiates an Elevator. alpar@399: /// kpeter@402: /// This function instantiates an \ref Elevator. alpar@399: /// \param digraph The digraph, to which we would like to define alpar@399: /// the elevator. alpar@399: /// \param max_level The maximum level of the elevator. alpar@399: static Elevator* createElevator(const Digraph& digraph, int max_level) { alpar@399: return new Elevator(digraph, max_level); alpar@399: } alpar@399: alpar@399: /// \brief The tolerance used by the algorithm alpar@399: /// alpar@399: /// The tolerance used by the algorithm to handle inexact computation. alpar@399: typedef lemon::Tolerance Tolerance; alpar@399: alpar@399: }; alpar@399: kpeter@402: /** kpeter@402: \brief Push-relabel algorithm for the network circulation problem. alpar@399: alpar@399: \ingroup max_flow kpeter@402: This class implements a push-relabel algorithm for the network kpeter@402: circulation problem. kpeter@402: It is to find a feasible circulation when lower and upper bounds kpeter@402: are given for the flow values on the arcs and lower bounds kpeter@402: are given for the supply values of the nodes. kpeter@402: alpar@399: The exact formulation of this problem is the following. kpeter@402: Let \f$G=(V,A)\f$ be a digraph, kpeter@402: \f$lower, upper: A\rightarrow\mathbf{R}^+_0\f$, kpeter@402: \f$delta: V\rightarrow\mathbf{R}\f$. Find a feasible circulation kpeter@402: \f$f: A\rightarrow\mathbf{R}^+_0\f$ so that kpeter@402: \f[ \sum_{a\in\delta_{out}(v)} f(a) - \sum_{a\in\delta_{in}(v)} f(a) kpeter@402: \geq delta(v) \quad \forall v\in V, \f] kpeter@402: \f[ lower(a)\leq f(a) \leq upper(a) \quad \forall a\in A. \f] kpeter@402: \note \f$delta(v)\f$ specifies a lower bound for the supply of node kpeter@402: \f$v\f$. It can be either positive or negative, however note that kpeter@402: \f$\sum_{v\in V}delta(v)\f$ should be zero or negative in order to kpeter@402: have a feasible solution. kpeter@402: kpeter@402: \note A special case of this problem is when kpeter@402: \f$\sum_{v\in V}delta(v) = 0\f$. Then the supply of each node \f$v\f$ kpeter@402: will be \e equal \e to \f$delta(v)\f$, if a circulation can be found. kpeter@402: Thus a feasible solution for the kpeter@402: \ref min_cost_flow "minimum cost flow" problem can be calculated kpeter@402: in this way. kpeter@402: kpeter@402: \tparam _Digraph The type of the digraph the algorithm runs on. kpeter@402: \tparam _LCapMap The type of the lower bound capacity map. The default kpeter@402: map type is \ref concepts::Digraph::ArcMap "_Digraph::ArcMap". kpeter@402: \tparam _UCapMap The type of the upper bound capacity map. The default kpeter@402: map type is \c _LCapMap. kpeter@402: \tparam _DeltaMap The type of the map that stores the lower bound kpeter@402: for the supply of the nodes. The default map type is kpeter@402: \c _Digraph::ArcMap<_UCapMap::Value>. alpar@399: */ kpeter@402: #ifdef DOXYGEN kpeter@402: template< typename _Digraph, kpeter@402: typename _LCapMap, kpeter@402: typename _UCapMap, kpeter@402: typename _DeltaMap, kpeter@402: typename _Traits > kpeter@402: #else kpeter@402: template< typename _Digraph, kpeter@402: typename _LCapMap = typename _Digraph::template ArcMap, kpeter@402: typename _UCapMap = _LCapMap, kpeter@402: typename _DeltaMap = typename _Digraph:: kpeter@402: template NodeMap, kpeter@402: typename _Traits=CirculationDefaultTraits<_Digraph, _LCapMap, kpeter@402: _UCapMap, _DeltaMap> > kpeter@402: #endif alpar@399: class Circulation { kpeter@402: public: alpar@399: kpeter@402: ///The \ref CirculationDefaultTraits "traits class" of the algorithm. alpar@399: typedef _Traits Traits; kpeter@402: ///The type of the digraph the algorithm runs on. alpar@399: typedef typename Traits::Digraph Digraph; kpeter@402: ///The type of the flow values. alpar@399: typedef typename Traits::Value Value; alpar@399: kpeter@402: /// The type of the lower bound capacity map. alpar@399: typedef typename Traits::LCapMap LCapMap; kpeter@402: /// The type of the upper bound capacity map. alpar@399: typedef typename Traits::UCapMap UCapMap; kpeter@402: /// \brief The type of the map that stores the lower bound for kpeter@402: /// the supply of the nodes. alpar@399: typedef typename Traits::DeltaMap DeltaMap; kpeter@402: ///The type of the flow map. alpar@399: typedef typename Traits::FlowMap FlowMap; kpeter@402: kpeter@402: ///The type of the elevator. alpar@399: typedef typename Traits::Elevator Elevator; kpeter@402: ///The type of the tolerance. alpar@399: typedef typename Traits::Tolerance Tolerance; alpar@399: kpeter@402: private: kpeter@402: kpeter@402: TEMPLATE_DIGRAPH_TYPEDEFS(Digraph); alpar@399: alpar@399: const Digraph &_g; alpar@399: int _node_num; alpar@399: alpar@399: const LCapMap *_lo; alpar@399: const UCapMap *_up; alpar@399: const DeltaMap *_delta; alpar@399: alpar@399: FlowMap *_flow; alpar@399: bool _local_flow; alpar@399: alpar@399: Elevator* _level; alpar@399: bool _local_level; alpar@399: kpeter@402: typedef typename Digraph::template NodeMap ExcessMap; alpar@399: ExcessMap* _excess; alpar@399: alpar@399: Tolerance _tol; alpar@399: int _el; alpar@399: alpar@399: public: alpar@399: alpar@399: typedef Circulation Create; alpar@399: kpeter@402: ///\name Named Template Parameters alpar@399: alpar@399: ///@{ alpar@399: alpar@399: template alpar@401: struct SetFlowMapTraits : public Traits { alpar@399: typedef _FlowMap FlowMap; alpar@399: static FlowMap *createFlowMap(const Digraph&) { alpar@399: LEMON_ASSERT(false, "FlowMap is not initialized"); alpar@399: return 0; // ignore warnings alpar@399: } alpar@399: }; alpar@399: alpar@399: /// \brief \ref named-templ-param "Named parameter" for setting alpar@399: /// FlowMap type alpar@399: /// alpar@399: /// \ref named-templ-param "Named parameter" for setting FlowMap kpeter@402: /// type. alpar@399: template alpar@401: struct SetFlowMap alpar@399: : public Circulation > { alpar@399: typedef Circulation > Create; alpar@399: }; alpar@399: alpar@399: template alpar@401: struct SetElevatorTraits : public Traits { alpar@399: typedef _Elevator Elevator; alpar@399: static Elevator *createElevator(const Digraph&, int) { alpar@399: LEMON_ASSERT(false, "Elevator is not initialized"); alpar@399: return 0; // ignore warnings alpar@399: } alpar@399: }; alpar@399: alpar@399: /// \brief \ref named-templ-param "Named parameter" for setting alpar@399: /// Elevator type alpar@399: /// alpar@399: /// \ref named-templ-param "Named parameter" for setting Elevator kpeter@402: /// type. If this named parameter is used, then an external kpeter@402: /// elevator object must be passed to the algorithm using the kpeter@402: /// \ref elevator(Elevator&) "elevator()" function before calling kpeter@402: /// \ref run() or \ref init(). kpeter@402: /// \sa SetStandardElevator alpar@399: template alpar@401: struct SetElevator alpar@399: : public Circulation > { alpar@399: typedef Circulation > Create; alpar@399: }; alpar@399: alpar@399: template alpar@401: struct SetStandardElevatorTraits : public Traits { alpar@399: typedef _Elevator Elevator; alpar@399: static Elevator *createElevator(const Digraph& digraph, int max_level) { alpar@399: return new Elevator(digraph, max_level); alpar@399: } alpar@399: }; alpar@399: alpar@399: /// \brief \ref named-templ-param "Named parameter" for setting kpeter@402: /// Elevator type with automatic allocation alpar@399: /// alpar@399: /// \ref named-templ-param "Named parameter" for setting Elevator kpeter@402: /// type with automatic allocation. kpeter@402: /// The Elevator should have standard constructor interface to be kpeter@402: /// able to automatically created by the algorithm (i.e. the kpeter@402: /// digraph and the maximum level should be passed to it). kpeter@402: /// However an external elevator object could also be passed to the kpeter@402: /// algorithm with the \ref elevator(Elevator&) "elevator()" function kpeter@402: /// before calling \ref run() or \ref init(). kpeter@402: /// \sa SetElevator alpar@399: template alpar@401: struct SetStandardElevator alpar@399: : public Circulation > { alpar@399: typedef Circulation > Create; alpar@399: }; alpar@399: alpar@399: /// @} alpar@399: alpar@399: protected: alpar@399: alpar@399: Circulation() {} alpar@399: alpar@399: public: alpar@399: alpar@399: /// The constructor of the class. alpar@399: alpar@399: /// The constructor of the class. alpar@399: /// \param g The digraph the algorithm runs on. alpar@399: /// \param lo The lower bound capacity of the arcs. alpar@399: /// \param up The upper bound capacity of the arcs. kpeter@402: /// \param delta The lower bound for the supply of the nodes. alpar@399: Circulation(const Digraph &g,const LCapMap &lo, alpar@399: const UCapMap &up,const DeltaMap &delta) alpar@399: : _g(g), _node_num(), alpar@399: _lo(&lo),_up(&up),_delta(&delta),_flow(0),_local_flow(false), alpar@399: _level(0), _local_level(false), _excess(0), _el() {} alpar@399: kpeter@402: /// Destructor. alpar@399: ~Circulation() { alpar@399: destroyStructures(); alpar@399: } alpar@399: kpeter@402: alpar@399: private: alpar@399: alpar@399: void createStructures() { alpar@399: _node_num = _el = countNodes(_g); alpar@399: alpar@399: if (!_flow) { alpar@399: _flow = Traits::createFlowMap(_g); alpar@399: _local_flow = true; alpar@399: } alpar@399: if (!_level) { alpar@399: _level = Traits::createElevator(_g, _node_num); alpar@399: _local_level = true; alpar@399: } alpar@399: if (!_excess) { alpar@399: _excess = new ExcessMap(_g); alpar@399: } alpar@399: } alpar@399: alpar@399: void destroyStructures() { alpar@399: if (_local_flow) { alpar@399: delete _flow; alpar@399: } alpar@399: if (_local_level) { alpar@399: delete _level; alpar@399: } alpar@399: if (_excess) { alpar@399: delete _excess; alpar@399: } alpar@399: } alpar@399: alpar@399: public: alpar@399: alpar@399: /// Sets the lower bound capacity map. alpar@399: alpar@399: /// Sets the lower bound capacity map. kpeter@402: /// \return (*this) alpar@399: Circulation& lowerCapMap(const LCapMap& map) { alpar@399: _lo = ↦ alpar@399: return *this; alpar@399: } alpar@399: alpar@399: /// Sets the upper bound capacity map. alpar@399: alpar@399: /// Sets the upper bound capacity map. kpeter@402: /// \return (*this) alpar@399: Circulation& upperCapMap(const LCapMap& map) { alpar@399: _up = ↦ alpar@399: return *this; alpar@399: } alpar@399: kpeter@402: /// Sets the lower bound map for the supply of the nodes. alpar@399: kpeter@402: /// Sets the lower bound map for the supply of the nodes. kpeter@402: /// \return (*this) alpar@399: Circulation& deltaMap(const DeltaMap& map) { alpar@399: _delta = ↦ alpar@399: return *this; alpar@399: } alpar@399: kpeter@402: /// \brief Sets the flow map. kpeter@402: /// alpar@399: /// Sets the flow map. kpeter@402: /// If you don't use this function before calling \ref run() or kpeter@402: /// \ref init(), an instance will be allocated automatically. kpeter@402: /// The destructor deallocates this automatically allocated map, kpeter@402: /// of course. kpeter@402: /// \return (*this) alpar@399: Circulation& flowMap(FlowMap& map) { alpar@399: if (_local_flow) { alpar@399: delete _flow; alpar@399: _local_flow = false; alpar@399: } alpar@399: _flow = ↦ alpar@399: return *this; alpar@399: } alpar@399: kpeter@402: /// \brief Sets the elevator used by algorithm. alpar@399: /// kpeter@402: /// Sets the elevator used by algorithm. kpeter@402: /// If you don't use this function before calling \ref run() or kpeter@402: /// \ref init(), an instance will be allocated automatically. kpeter@402: /// The destructor deallocates this automatically allocated elevator, kpeter@402: /// of course. kpeter@402: /// \return (*this) alpar@399: Circulation& elevator(Elevator& elevator) { alpar@399: if (_local_level) { alpar@399: delete _level; alpar@399: _local_level = false; alpar@399: } alpar@399: _level = &elevator; alpar@399: return *this; alpar@399: } alpar@399: kpeter@402: /// \brief Returns a const reference to the elevator. alpar@399: /// kpeter@402: /// Returns a const reference to the elevator. kpeter@402: /// kpeter@402: /// \pre Either \ref run() or \ref init() must be called before kpeter@402: /// using this function. alpar@399: const Elevator& elevator() { alpar@399: return *_level; alpar@399: } alpar@399: kpeter@402: /// \brief Sets the tolerance used by algorithm. kpeter@402: /// alpar@399: /// Sets the tolerance used by algorithm. alpar@399: Circulation& tolerance(const Tolerance& tolerance) const { alpar@399: _tol = tolerance; alpar@399: return *this; alpar@399: } alpar@399: kpeter@402: /// \brief Returns a const reference to the tolerance. alpar@399: /// kpeter@402: /// Returns a const reference to the tolerance. alpar@399: const Tolerance& tolerance() const { alpar@399: return tolerance; alpar@399: } alpar@399: kpeter@402: /// \name Execution Control kpeter@402: /// The simplest way to execute the algorithm is to call \ref run().\n kpeter@402: /// If you need more control on the initial solution or the execution, kpeter@402: /// first you have to call one of the \ref init() functions, then kpeter@402: /// the \ref start() function. alpar@399: alpar@399: ///@{ alpar@399: alpar@399: /// Initializes the internal data structures. alpar@399: kpeter@402: /// Initializes the internal data structures and sets all flow values kpeter@402: /// to the lower bound. alpar@399: void init() alpar@399: { alpar@399: createStructures(); alpar@399: alpar@399: for(NodeIt n(_g);n!=INVALID;++n) { alpar@399: _excess->set(n, (*_delta)[n]); alpar@399: } alpar@399: alpar@399: for (ArcIt e(_g);e!=INVALID;++e) { alpar@399: _flow->set(e, (*_lo)[e]); alpar@399: _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_flow)[e]); alpar@399: _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_flow)[e]); alpar@399: } alpar@399: alpar@399: // global relabeling tested, but in general case it provides alpar@399: // worse performance for random digraphs alpar@399: _level->initStart(); alpar@399: for(NodeIt n(_g);n!=INVALID;++n) alpar@399: _level->initAddItem(n); alpar@399: _level->initFinish(); alpar@399: for(NodeIt n(_g);n!=INVALID;++n) alpar@399: if(_tol.positive((*_excess)[n])) alpar@399: _level->activate(n); alpar@399: } alpar@399: kpeter@402: /// Initializes the internal data structures using a greedy approach. alpar@399: kpeter@402: /// Initializes the internal data structures using a greedy approach kpeter@402: /// to construct the initial solution. alpar@399: void greedyInit() alpar@399: { alpar@399: createStructures(); alpar@399: alpar@399: for(NodeIt n(_g);n!=INVALID;++n) { alpar@399: _excess->set(n, (*_delta)[n]); alpar@399: } alpar@399: alpar@399: for (ArcIt e(_g);e!=INVALID;++e) { alpar@399: if (!_tol.positive((*_excess)[_g.target(e)] + (*_up)[e])) { alpar@399: _flow->set(e, (*_up)[e]); alpar@399: _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_up)[e]); alpar@399: _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_up)[e]); alpar@399: } else if (_tol.positive((*_excess)[_g.target(e)] + (*_lo)[e])) { alpar@399: _flow->set(e, (*_lo)[e]); alpar@399: _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_lo)[e]); alpar@399: _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_lo)[e]); alpar@399: } else { alpar@399: Value fc = -(*_excess)[_g.target(e)]; alpar@399: _flow->set(e, fc); alpar@399: _excess->set(_g.target(e), 0); alpar@399: _excess->set(_g.source(e), (*_excess)[_g.source(e)] - fc); alpar@399: } alpar@399: } alpar@399: alpar@399: _level->initStart(); alpar@399: for(NodeIt n(_g);n!=INVALID;++n) alpar@399: _level->initAddItem(n); alpar@399: _level->initFinish(); alpar@399: for(NodeIt n(_g);n!=INVALID;++n) alpar@399: if(_tol.positive((*_excess)[n])) alpar@399: _level->activate(n); alpar@399: } alpar@399: kpeter@402: ///Executes the algorithm alpar@399: kpeter@402: ///This function executes the algorithm. kpeter@402: /// kpeter@402: ///\return \c true if a feasible circulation is found. alpar@399: /// alpar@399: ///\sa barrier() kpeter@402: ///\sa barrierMap() alpar@399: bool start() alpar@399: { alpar@399: alpar@399: Node act; alpar@399: Node bact=INVALID; alpar@399: Node last_activated=INVALID; alpar@399: while((act=_level->highestActive())!=INVALID) { alpar@399: int actlevel=(*_level)[act]; alpar@399: int mlevel=_node_num; alpar@399: Value exc=(*_excess)[act]; alpar@399: alpar@399: for(OutArcIt e(_g,act);e!=INVALID; ++e) { alpar@399: Node v = _g.target(e); alpar@399: Value fc=(*_up)[e]-(*_flow)[e]; alpar@399: if(!_tol.positive(fc)) continue; alpar@399: if((*_level)[v]set(e, (*_flow)[e] + exc); alpar@399: _excess->set(v, (*_excess)[v] + exc); alpar@399: if(!_level->active(v) && _tol.positive((*_excess)[v])) alpar@399: _level->activate(v); alpar@399: _excess->set(act,0); alpar@399: _level->deactivate(act); alpar@399: goto next_l; alpar@399: } alpar@399: else { alpar@399: _flow->set(e, (*_up)[e]); alpar@399: _excess->set(v, (*_excess)[v] + fc); alpar@399: if(!_level->active(v) && _tol.positive((*_excess)[v])) alpar@399: _level->activate(v); alpar@399: exc-=fc; alpar@399: } alpar@399: } alpar@399: else if((*_level)[v]set(e, (*_flow)[e] - exc); alpar@399: _excess->set(v, (*_excess)[v] + exc); alpar@399: if(!_level->active(v) && _tol.positive((*_excess)[v])) alpar@399: _level->activate(v); alpar@399: _excess->set(act,0); alpar@399: _level->deactivate(act); alpar@399: goto next_l; alpar@399: } alpar@399: else { alpar@399: _flow->set(e, (*_lo)[e]); alpar@399: _excess->set(v, (*_excess)[v] + fc); alpar@399: if(!_level->active(v) && _tol.positive((*_excess)[v])) alpar@399: _level->activate(v); alpar@399: exc-=fc; alpar@399: } alpar@399: } alpar@399: else if((*_level)[v]set(act, exc); alpar@399: if(!_tol.positive(exc)) _level->deactivate(act); alpar@399: else if(mlevel==_node_num) { alpar@399: _level->liftHighestActiveToTop(); alpar@399: _el = _node_num; alpar@399: return false; alpar@399: } alpar@399: else { alpar@399: _level->liftHighestActive(mlevel+1); alpar@399: if(_level->onLevel(actlevel)==0) { alpar@399: _el = actlevel; alpar@399: return false; alpar@399: } alpar@399: } alpar@399: next_l: alpar@399: ; alpar@399: } alpar@399: return true; alpar@399: } alpar@399: kpeter@402: /// Runs the algorithm. alpar@399: kpeter@402: /// This function runs the algorithm. kpeter@402: /// kpeter@402: /// \return \c true if a feasible circulation is found. kpeter@402: /// kpeter@402: /// \note Apart from the return value, c.run() is just a shortcut of kpeter@402: /// the following code. alpar@399: /// \code kpeter@402: /// c.greedyInit(); kpeter@402: /// c.start(); alpar@399: /// \endcode alpar@399: bool run() { alpar@399: greedyInit(); alpar@399: return start(); alpar@399: } alpar@399: alpar@399: /// @} alpar@399: alpar@399: /// \name Query Functions kpeter@402: /// The results of the circulation algorithm can be obtained using kpeter@402: /// these functions.\n kpeter@402: /// Either \ref run() or \ref start() should be called before kpeter@402: /// using them. alpar@399: alpar@399: ///@{ alpar@399: kpeter@402: /// \brief Returns the flow on the given arc. kpeter@402: /// kpeter@402: /// Returns the flow on the given arc. kpeter@402: /// kpeter@402: /// \pre Either \ref run() or \ref init() must be called before kpeter@402: /// using this function. kpeter@402: Value flow(const Arc& arc) const { kpeter@402: return (*_flow)[arc]; kpeter@402: } kpeter@402: kpeter@402: /// \brief Returns a const reference to the flow map. kpeter@402: /// kpeter@402: /// Returns a const reference to the arc map storing the found flow. kpeter@402: /// kpeter@402: /// \pre Either \ref run() or \ref init() must be called before kpeter@402: /// using this function. kpeter@402: const FlowMap& flowMap() { kpeter@402: return *_flow; kpeter@402: } kpeter@402: alpar@399: /** kpeter@402: \brief Returns \c true if the given node is in a barrier. kpeter@402: alpar@399: Barrier is a set \e B of nodes for which kpeter@402: kpeter@402: \f[ \sum_{a\in\delta_{out}(B)} upper(a) - kpeter@402: \sum_{a\in\delta_{in}(B)} lower(a) < \sum_{v\in B}delta(v) \f] kpeter@402: kpeter@402: holds. The existence of a set with this property prooves that a kpeter@402: feasible circualtion cannot exist. kpeter@402: kpeter@402: This function returns \c true if the given node is in the found kpeter@402: barrier. If a feasible circulation is found, the function kpeter@402: gives back \c false for every node. kpeter@402: kpeter@402: \pre Either \ref run() or \ref init() must be called before kpeter@402: using this function. kpeter@402: kpeter@402: \sa barrierMap() alpar@399: \sa checkBarrier() alpar@399: */ kpeter@402: bool barrier(const Node& node) kpeter@402: { kpeter@402: return (*_level)[node] >= _el; kpeter@402: } kpeter@402: kpeter@402: /// \brief Gives back a barrier. kpeter@402: /// kpeter@402: /// This function sets \c bar to the characteristic vector of the kpeter@402: /// found barrier. \c bar should be a \ref concepts::WriteMap "writable" kpeter@402: /// node map with \c bool (or convertible) value type. kpeter@402: /// kpeter@402: /// If a feasible circulation is found, the function gives back an kpeter@402: /// empty set, so \c bar[v] will be \c false for all nodes \c v. kpeter@402: /// kpeter@402: /// \note This function calls \ref barrier() for each node, kpeter@402: /// so it runs in \f$O(n)\f$ time. kpeter@402: /// kpeter@402: /// \pre Either \ref run() or \ref init() must be called before kpeter@402: /// using this function. kpeter@402: /// kpeter@402: /// \sa barrier() kpeter@402: /// \sa checkBarrier() kpeter@402: template kpeter@402: void barrierMap(BarrierMap &bar) alpar@399: { alpar@399: for(NodeIt n(_g);n!=INVALID;++n) alpar@399: bar.set(n, (*_level)[n] >= _el); alpar@399: } alpar@399: alpar@399: /// @} alpar@399: alpar@399: /// \name Checker Functions kpeter@402: /// The feasibility of the results can be checked using kpeter@402: /// these functions.\n kpeter@402: /// Either \ref run() or \ref start() should be called before kpeter@402: /// using them. alpar@399: alpar@399: ///@{ alpar@399: kpeter@402: ///Check if the found flow is a feasible circulation kpeter@402: kpeter@402: ///Check if the found flow is a feasible circulation, kpeter@402: /// alpar@399: bool checkFlow() { alpar@399: for(ArcIt e(_g);e!=INVALID;++e) alpar@399: if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false; alpar@399: for(NodeIt n(_g);n!=INVALID;++n) alpar@399: { alpar@399: Value dif=-(*_delta)[n]; alpar@399: for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e]; alpar@399: for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e]; alpar@399: if(_tol.negative(dif)) return false; alpar@399: } alpar@399: return true; alpar@399: } alpar@399: alpar@399: ///Check whether or not the last execution provides a barrier alpar@399: kpeter@402: ///Check whether or not the last execution provides a barrier. alpar@399: ///\sa barrier() kpeter@402: ///\sa barrierMap() alpar@399: bool checkBarrier() alpar@399: { alpar@399: Value delta=0; alpar@399: for(NodeIt n(_g);n!=INVALID;++n) alpar@399: if(barrier(n)) alpar@399: delta-=(*_delta)[n]; alpar@399: for(ArcIt e(_g);e!=INVALID;++e) alpar@399: { alpar@399: Node s=_g.source(e); alpar@399: Node t=_g.target(e); alpar@399: if(barrier(s)&&!barrier(t)) delta+=(*_up)[e]; alpar@399: else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e]; alpar@399: } alpar@399: return _tol.negative(delta); alpar@399: } alpar@399: alpar@399: /// @} alpar@399: alpar@399: }; alpar@399: alpar@399: } alpar@399: alpar@399: #endif