COIN-OR::LEMON - Graph Library

source: lemon/lemon/circulation.h @ 760:4ac30454f1c1

Last change on this file since 760:4ac30454f1c1 was 760:4ac30454f1c1, checked in by Peter Kovacs <kpeter@…>, 15 years ago

Small doc improvements

File size: 24.6 KB
Line 
1/* -*- mode: C++; indent-tabs-mode: nil; -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library.
4 *
5 * Copyright (C) 2003-2009
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 *
9 * Permission to use, modify and distribute this software is granted
10 * provided that this copyright notice appears in all copies. For
11 * precise terms see the accompanying LICENSE file.
12 *
13 * This software is provided "AS IS" with no warranty of any kind,
14 * express or implied, and with no claim as to its suitability for any
15 * purpose.
16 *
17 */
18
19#ifndef LEMON_CIRCULATION_H
20#define LEMON_CIRCULATION_H
21
22#include <lemon/tolerance.h>
23#include <lemon/elevator.h>
24#include <limits>
25
26///\ingroup max_flow
27///\file
28///\brief Push-relabel algorithm for finding a feasible circulation.
29///
30namespace lemon {
31
32  /// \brief Default traits class of Circulation class.
33  ///
34  /// Default traits class of Circulation class.
35  ///
36  /// \tparam GR Type of the digraph the algorithm runs on.
37  /// \tparam LM The type of the lower bound map.
38  /// \tparam UM The type of the upper bound (capacity) map.
39  /// \tparam SM The type of the supply map.
40  template <typename GR, typename LM,
41            typename UM, typename SM>
42  struct CirculationDefaultTraits {
43
44    /// \brief The type of the digraph the algorithm runs on.
45    typedef GR Digraph;
46
47    /// \brief The type of the lower bound map.
48    ///
49    /// The type of the map that stores the lower bounds on the arcs.
50    /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
51    typedef LM LowerMap;
52
53    /// \brief The type of the upper bound (capacity) map.
54    ///
55    /// The type of the map that stores the upper bounds (capacities)
56    /// on the arcs.
57    /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
58    typedef UM UpperMap;
59
60    /// \brief The type of supply map.
61    ///
62    /// The type of the map that stores the signed supply values of the
63    /// nodes.
64    /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
65    typedef SM SupplyMap;
66
67    /// \brief The type of the flow and supply values.
68    typedef typename SupplyMap::Value Value;
69
70    /// \brief The type of the map that stores the flow values.
71    ///
72    /// The type of the map that stores the flow values.
73    /// It must conform to the \ref concepts::ReadWriteMap "ReadWriteMap"
74    /// concept.
75#ifdef DOXYGEN
76    typedef GR::ArcMap<Value> FlowMap;
77#else
78    typedef typename Digraph::template ArcMap<Value> FlowMap;
79#endif
80
81    /// \brief Instantiates a FlowMap.
82    ///
83    /// This function instantiates a \ref FlowMap.
84    /// \param digraph The digraph for which we would like to define
85    /// the flow map.
86    static FlowMap* createFlowMap(const Digraph& digraph) {
87      return new FlowMap(digraph);
88    }
89
90    /// \brief The elevator type used by the algorithm.
91    ///
92    /// The elevator type used by the algorithm.
93    ///
94    /// \sa Elevator, LinkedElevator
95#ifdef DOXYGEN
96    typedef lemon::Elevator<GR, GR::Node> Elevator;
97#else
98    typedef lemon::Elevator<Digraph, typename Digraph::Node> Elevator;
99#endif
100
101    /// \brief Instantiates an Elevator.
102    ///
103    /// This function instantiates an \ref Elevator.
104    /// \param digraph The digraph for which we would like to define
105    /// the elevator.
106    /// \param max_level The maximum level of the elevator.
107    static Elevator* createElevator(const Digraph& digraph, int max_level) {
108      return new Elevator(digraph, max_level);
109    }
110
111    /// \brief The tolerance used by the algorithm
112    ///
113    /// The tolerance used by the algorithm to handle inexact computation.
114    typedef lemon::Tolerance<Value> Tolerance;
115
116  };
117
118  /**
119     \brief Push-relabel algorithm for the network circulation problem.
120
121     \ingroup max_flow
122     This class implements a push-relabel algorithm for the \e network
123     \e circulation problem.
124     It is to find a feasible circulation when lower and upper bounds
125     are given for the flow values on the arcs and lower bounds are
126     given for the difference between the outgoing and incoming flow
127     at the nodes.
128
129     The exact formulation of this problem is the following.
130     Let \f$G=(V,A)\f$ be a digraph, \f$lower: A\rightarrow\mathbf{R}\f$
131     \f$upper: A\rightarrow\mathbf{R}\cup\{\infty\}\f$ denote the lower and
132     upper bounds on the arcs, for which \f$lower(uv) \leq upper(uv)\f$
133     holds for all \f$uv\in A\f$, and \f$sup: V\rightarrow\mathbf{R}\f$
134     denotes the signed supply values of the nodes.
135     If \f$sup(u)>0\f$, then \f$u\f$ is a supply node with \f$sup(u)\f$
136     supply, if \f$sup(u)<0\f$, then \f$u\f$ is a demand node with
137     \f$-sup(u)\f$ demand.
138     A feasible circulation is an \f$f: A\rightarrow\mathbf{R}\f$
139     solution of the following problem.
140
141     \f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu)
142     \geq sup(u) \quad \forall u\in V, \f]
143     \f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A. \f]
144     
145     The sum of the supply values, i.e. \f$\sum_{u\in V} sup(u)\f$ must be
146     zero or negative in order to have a feasible solution (since the sum
147     of the expressions on the left-hand side of the inequalities is zero).
148     It means that the total demand must be greater or equal to the total
149     supply and all the supplies have to be carried out from the supply nodes,
150     but there could be demands that are not satisfied.
151     If \f$\sum_{u\in V} sup(u)\f$ is zero, then all the supply/demand
152     constraints have to be satisfied with equality, i.e. all demands
153     have to be satisfied and all supplies have to be used.
154     
155     If you need the opposite inequalities in the supply/demand constraints
156     (i.e. the total demand is less than the total supply and all the demands
157     have to be satisfied while there could be supplies that are not used),
158     then you could easily transform the problem to the above form by reversing
159     the direction of the arcs and taking the negative of the supply values
160     (e.g. using \ref ReverseDigraph and \ref NegMap adaptors).
161
162     This algorithm either calculates a feasible circulation, or provides
163     a \ref barrier() "barrier", which prooves that a feasible soultion
164     cannot exist.
165
166     Note that this algorithm also provides a feasible solution for the
167     \ref min_cost_flow "minimum cost flow problem".
168
169     \tparam GR The type of the digraph the algorithm runs on.
170     \tparam LM The type of the lower bound map. The default
171     map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".
172     \tparam UM The type of the upper bound (capacity) map.
173     The default map type is \c LM.
174     \tparam SM The type of the supply map. The default map type is
175     \ref concepts::Digraph::NodeMap "GR::NodeMap<UM::Value>".
176  */
177#ifdef DOXYGEN
178template< typename GR,
179          typename LM,
180          typename UM,
181          typename SM,
182          typename TR >
183#else
184template< typename GR,
185          typename LM = typename GR::template ArcMap<int>,
186          typename UM = LM,
187          typename SM = typename GR::template NodeMap<typename UM::Value>,
188          typename TR = CirculationDefaultTraits<GR, LM, UM, SM> >
189#endif
190  class Circulation {
191  public:
192
193    ///The \ref CirculationDefaultTraits "traits class" of the algorithm.
194    typedef TR Traits;
195    ///The type of the digraph the algorithm runs on.
196    typedef typename Traits::Digraph Digraph;
197    ///The type of the flow and supply values.
198    typedef typename Traits::Value Value;
199
200    ///The type of the lower bound map.
201    typedef typename Traits::LowerMap LowerMap;
202    ///The type of the upper bound (capacity) map.
203    typedef typename Traits::UpperMap UpperMap;
204    ///The type of the supply map.
205    typedef typename Traits::SupplyMap SupplyMap;
206    ///The type of the flow map.
207    typedef typename Traits::FlowMap FlowMap;
208
209    ///The type of the elevator.
210    typedef typename Traits::Elevator Elevator;
211    ///The type of the tolerance.
212    typedef typename Traits::Tolerance Tolerance;
213
214  private:
215
216    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
217
218    const Digraph &_g;
219    int _node_num;
220
221    const LowerMap *_lo;
222    const UpperMap *_up;
223    const SupplyMap *_supply;
224
225    FlowMap *_flow;
226    bool _local_flow;
227
228    Elevator* _level;
229    bool _local_level;
230
231    typedef typename Digraph::template NodeMap<Value> ExcessMap;
232    ExcessMap* _excess;
233
234    Tolerance _tol;
235    int _el;
236
237  public:
238
239    typedef Circulation Create;
240
241    ///\name Named Template Parameters
242
243    ///@{
244
245    template <typename T>
246    struct SetFlowMapTraits : public Traits {
247      typedef T FlowMap;
248      static FlowMap *createFlowMap(const Digraph&) {
249        LEMON_ASSERT(false, "FlowMap is not initialized");
250        return 0; // ignore warnings
251      }
252    };
253
254    /// \brief \ref named-templ-param "Named parameter" for setting
255    /// FlowMap type
256    ///
257    /// \ref named-templ-param "Named parameter" for setting FlowMap
258    /// type.
259    template <typename T>
260    struct SetFlowMap
261      : public Circulation<Digraph, LowerMap, UpperMap, SupplyMap,
262                           SetFlowMapTraits<T> > {
263      typedef Circulation<Digraph, LowerMap, UpperMap, SupplyMap,
264                          SetFlowMapTraits<T> > Create;
265    };
266
267    template <typename T>
268    struct SetElevatorTraits : public Traits {
269      typedef T Elevator;
270      static Elevator *createElevator(const Digraph&, int) {
271        LEMON_ASSERT(false, "Elevator is not initialized");
272        return 0; // ignore warnings
273      }
274    };
275
276    /// \brief \ref named-templ-param "Named parameter" for setting
277    /// Elevator type
278    ///
279    /// \ref named-templ-param "Named parameter" for setting Elevator
280    /// type. If this named parameter is used, then an external
281    /// elevator object must be passed to the algorithm using the
282    /// \ref elevator(Elevator&) "elevator()" function before calling
283    /// \ref run() or \ref init().
284    /// \sa SetStandardElevator
285    template <typename T>
286    struct SetElevator
287      : public Circulation<Digraph, LowerMap, UpperMap, SupplyMap,
288                           SetElevatorTraits<T> > {
289      typedef Circulation<Digraph, LowerMap, UpperMap, SupplyMap,
290                          SetElevatorTraits<T> > Create;
291    };
292
293    template <typename T>
294    struct SetStandardElevatorTraits : public Traits {
295      typedef T Elevator;
296      static Elevator *createElevator(const Digraph& digraph, int max_level) {
297        return new Elevator(digraph, max_level);
298      }
299    };
300
301    /// \brief \ref named-templ-param "Named parameter" for setting
302    /// Elevator type with automatic allocation
303    ///
304    /// \ref named-templ-param "Named parameter" for setting Elevator
305    /// type with automatic allocation.
306    /// The Elevator should have standard constructor interface to be
307    /// able to automatically created by the algorithm (i.e. the
308    /// digraph and the maximum level should be passed to it).
309    /// However an external elevator object could also be passed to the
310    /// algorithm with the \ref elevator(Elevator&) "elevator()" function
311    /// before calling \ref run() or \ref init().
312    /// \sa SetElevator
313    template <typename T>
314    struct SetStandardElevator
315      : public Circulation<Digraph, LowerMap, UpperMap, SupplyMap,
316                       SetStandardElevatorTraits<T> > {
317      typedef Circulation<Digraph, LowerMap, UpperMap, SupplyMap,
318                      SetStandardElevatorTraits<T> > Create;
319    };
320
321    /// @}
322
323  protected:
324
325    Circulation() {}
326
327  public:
328
329    /// Constructor.
330
331    /// The constructor of the class.
332    ///
333    /// \param graph The digraph the algorithm runs on.
334    /// \param lower The lower bounds for the flow values on the arcs.
335    /// \param upper The upper bounds (capacities) for the flow values
336    /// on the arcs.
337    /// \param supply The signed supply values of the nodes.
338    Circulation(const Digraph &graph, const LowerMap &lower,
339                const UpperMap &upper, const SupplyMap &supply)
340      : _g(graph), _lo(&lower), _up(&upper), _supply(&supply),
341        _flow(NULL), _local_flow(false), _level(NULL), _local_level(false),
342        _excess(NULL) {}
343
344    /// Destructor.
345    ~Circulation() {
346      destroyStructures();
347    }
348
349
350  private:
351
352    bool checkBoundMaps() {
353      for (ArcIt e(_g);e!=INVALID;++e) {
354        if (_tol.less((*_up)[e], (*_lo)[e])) return false;
355      }
356      return true;
357    }
358
359    void createStructures() {
360      _node_num = _el = countNodes(_g);
361
362      if (!_flow) {
363        _flow = Traits::createFlowMap(_g);
364        _local_flow = true;
365      }
366      if (!_level) {
367        _level = Traits::createElevator(_g, _node_num);
368        _local_level = true;
369      }
370      if (!_excess) {
371        _excess = new ExcessMap(_g);
372      }
373    }
374
375    void destroyStructures() {
376      if (_local_flow) {
377        delete _flow;
378      }
379      if (_local_level) {
380        delete _level;
381      }
382      if (_excess) {
383        delete _excess;
384      }
385    }
386
387  public:
388
389    /// Sets the lower bound map.
390
391    /// Sets the lower bound map.
392    /// \return <tt>(*this)</tt>
393    Circulation& lowerMap(const LowerMap& map) {
394      _lo = &map;
395      return *this;
396    }
397
398    /// Sets the upper bound (capacity) map.
399
400    /// Sets the upper bound (capacity) map.
401    /// \return <tt>(*this)</tt>
402    Circulation& upperMap(const UpperMap& map) {
403      _up = &map;
404      return *this;
405    }
406
407    /// Sets the supply map.
408
409    /// Sets the supply map.
410    /// \return <tt>(*this)</tt>
411    Circulation& supplyMap(const SupplyMap& map) {
412      _supply = &map;
413      return *this;
414    }
415
416    /// \brief Sets the flow map.
417    ///
418    /// Sets the flow map.
419    /// If you don't use this function before calling \ref run() or
420    /// \ref init(), an instance will be allocated automatically.
421    /// The destructor deallocates this automatically allocated map,
422    /// of course.
423    /// \return <tt>(*this)</tt>
424    Circulation& flowMap(FlowMap& map) {
425      if (_local_flow) {
426        delete _flow;
427        _local_flow = false;
428      }
429      _flow = &map;
430      return *this;
431    }
432
433    /// \brief Sets the elevator used by algorithm.
434    ///
435    /// Sets the elevator used by algorithm.
436    /// If you don't use this function before calling \ref run() or
437    /// \ref init(), an instance will be allocated automatically.
438    /// The destructor deallocates this automatically allocated elevator,
439    /// of course.
440    /// \return <tt>(*this)</tt>
441    Circulation& elevator(Elevator& elevator) {
442      if (_local_level) {
443        delete _level;
444        _local_level = false;
445      }
446      _level = &elevator;
447      return *this;
448    }
449
450    /// \brief Returns a const reference to the elevator.
451    ///
452    /// Returns a const reference to the elevator.
453    ///
454    /// \pre Either \ref run() or \ref init() must be called before
455    /// using this function.
456    const Elevator& elevator() const {
457      return *_level;
458    }
459
460    /// \brief Sets the tolerance used by algorithm.
461    ///
462    /// Sets the tolerance used by algorithm.
463    Circulation& tolerance(const Tolerance& tolerance) const {
464      _tol = tolerance;
465      return *this;
466    }
467
468    /// \brief Returns a const reference to the tolerance.
469    ///
470    /// Returns a const reference to the tolerance.
471    const Tolerance& tolerance() const {
472      return tolerance;
473    }
474
475    /// \name Execution Control
476    /// The simplest way to execute the algorithm is to call \ref run().\n
477    /// If you need better control on the initial solution or the execution,
478    /// you have to call one of the \ref init() functions first, then
479    /// the \ref start() function.
480
481    ///@{
482
483    /// Initializes the internal data structures.
484
485    /// Initializes the internal data structures and sets all flow values
486    /// to the lower bound.
487    void init()
488    {
489      LEMON_DEBUG(checkBoundMaps(),
490        "Upper bounds must be greater or equal to the lower bounds");
491
492      createStructures();
493
494      for(NodeIt n(_g);n!=INVALID;++n) {
495        (*_excess)[n] = (*_supply)[n];
496      }
497
498      for (ArcIt e(_g);e!=INVALID;++e) {
499        _flow->set(e, (*_lo)[e]);
500        (*_excess)[_g.target(e)] += (*_flow)[e];
501        (*_excess)[_g.source(e)] -= (*_flow)[e];
502      }
503
504      // global relabeling tested, but in general case it provides
505      // worse performance for random digraphs
506      _level->initStart();
507      for(NodeIt n(_g);n!=INVALID;++n)
508        _level->initAddItem(n);
509      _level->initFinish();
510      for(NodeIt n(_g);n!=INVALID;++n)
511        if(_tol.positive((*_excess)[n]))
512          _level->activate(n);
513    }
514
515    /// Initializes the internal data structures using a greedy approach.
516
517    /// Initializes the internal data structures using a greedy approach
518    /// to construct the initial solution.
519    void greedyInit()
520    {
521      LEMON_DEBUG(checkBoundMaps(),
522        "Upper bounds must be greater or equal to the lower bounds");
523
524      createStructures();
525
526      for(NodeIt n(_g);n!=INVALID;++n) {
527        (*_excess)[n] = (*_supply)[n];
528      }
529
530      for (ArcIt e(_g);e!=INVALID;++e) {
531        if (!_tol.less(-(*_excess)[_g.target(e)], (*_up)[e])) {
532          _flow->set(e, (*_up)[e]);
533          (*_excess)[_g.target(e)] += (*_up)[e];
534          (*_excess)[_g.source(e)] -= (*_up)[e];
535        } else if (_tol.less(-(*_excess)[_g.target(e)], (*_lo)[e])) {
536          _flow->set(e, (*_lo)[e]);
537          (*_excess)[_g.target(e)] += (*_lo)[e];
538          (*_excess)[_g.source(e)] -= (*_lo)[e];
539        } else {
540          Value fc = -(*_excess)[_g.target(e)];
541          _flow->set(e, fc);
542          (*_excess)[_g.target(e)] = 0;
543          (*_excess)[_g.source(e)] -= fc;
544        }
545      }
546
547      _level->initStart();
548      for(NodeIt n(_g);n!=INVALID;++n)
549        _level->initAddItem(n);
550      _level->initFinish();
551      for(NodeIt n(_g);n!=INVALID;++n)
552        if(_tol.positive((*_excess)[n]))
553          _level->activate(n);
554    }
555
556    ///Executes the algorithm
557
558    ///This function executes the algorithm.
559    ///
560    ///\return \c true if a feasible circulation is found.
561    ///
562    ///\sa barrier()
563    ///\sa barrierMap()
564    bool start()
565    {
566
567      Node act;
568      Node bact=INVALID;
569      Node last_activated=INVALID;
570      while((act=_level->highestActive())!=INVALID) {
571        int actlevel=(*_level)[act];
572        int mlevel=_node_num;
573        Value exc=(*_excess)[act];
574
575        for(OutArcIt e(_g,act);e!=INVALID; ++e) {
576          Node v = _g.target(e);
577          Value fc=(*_up)[e]-(*_flow)[e];
578          if(!_tol.positive(fc)) continue;
579          if((*_level)[v]<actlevel) {
580            if(!_tol.less(fc, exc)) {
581              _flow->set(e, (*_flow)[e] + exc);
582              (*_excess)[v] += exc;
583              if(!_level->active(v) && _tol.positive((*_excess)[v]))
584                _level->activate(v);
585              (*_excess)[act] = 0;
586              _level->deactivate(act);
587              goto next_l;
588            }
589            else {
590              _flow->set(e, (*_up)[e]);
591              (*_excess)[v] += fc;
592              if(!_level->active(v) && _tol.positive((*_excess)[v]))
593                _level->activate(v);
594              exc-=fc;
595            }
596          }
597          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
598        }
599        for(InArcIt e(_g,act);e!=INVALID; ++e) {
600          Node v = _g.source(e);
601          Value fc=(*_flow)[e]-(*_lo)[e];
602          if(!_tol.positive(fc)) continue;
603          if((*_level)[v]<actlevel) {
604            if(!_tol.less(fc, exc)) {
605              _flow->set(e, (*_flow)[e] - exc);
606              (*_excess)[v] += exc;
607              if(!_level->active(v) && _tol.positive((*_excess)[v]))
608                _level->activate(v);
609              (*_excess)[act] = 0;
610              _level->deactivate(act);
611              goto next_l;
612            }
613            else {
614              _flow->set(e, (*_lo)[e]);
615              (*_excess)[v] += fc;
616              if(!_level->active(v) && _tol.positive((*_excess)[v]))
617                _level->activate(v);
618              exc-=fc;
619            }
620          }
621          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
622        }
623
624        (*_excess)[act] = exc;
625        if(!_tol.positive(exc)) _level->deactivate(act);
626        else if(mlevel==_node_num) {
627          _level->liftHighestActiveToTop();
628          _el = _node_num;
629          return false;
630        }
631        else {
632          _level->liftHighestActive(mlevel+1);
633          if(_level->onLevel(actlevel)==0) {
634            _el = actlevel;
635            return false;
636          }
637        }
638      next_l:
639        ;
640      }
641      return true;
642    }
643
644    /// Runs the algorithm.
645
646    /// This function runs the algorithm.
647    ///
648    /// \return \c true if a feasible circulation is found.
649    ///
650    /// \note Apart from the return value, c.run() is just a shortcut of
651    /// the following code.
652    /// \code
653    ///   c.greedyInit();
654    ///   c.start();
655    /// \endcode
656    bool run() {
657      greedyInit();
658      return start();
659    }
660
661    /// @}
662
663    /// \name Query Functions
664    /// The results of the circulation algorithm can be obtained using
665    /// these functions.\n
666    /// Either \ref run() or \ref start() should be called before
667    /// using them.
668
669    ///@{
670
671    /// \brief Returns the flow value on the given arc.
672    ///
673    /// Returns the flow value on the given arc.
674    ///
675    /// \pre Either \ref run() or \ref init() must be called before
676    /// using this function.
677    Value flow(const Arc& arc) const {
678      return (*_flow)[arc];
679    }
680
681    /// \brief Returns a const reference to the flow map.
682    ///
683    /// Returns a const reference to the arc map storing the found flow.
684    ///
685    /// \pre Either \ref run() or \ref init() must be called before
686    /// using this function.
687    const FlowMap& flowMap() const {
688      return *_flow;
689    }
690
691    /**
692       \brief Returns \c true if the given node is in a barrier.
693
694       Barrier is a set \e B of nodes for which
695
696       \f[ \sum_{uv\in A: u\in B} upper(uv) -
697           \sum_{uv\in A: v\in B} lower(uv) < \sum_{v\in B} sup(v) \f]
698
699       holds. The existence of a set with this property prooves that a
700       feasible circualtion cannot exist.
701
702       This function returns \c true if the given node is in the found
703       barrier. If a feasible circulation is found, the function
704       gives back \c false for every node.
705
706       \pre Either \ref run() or \ref init() must be called before
707       using this function.
708
709       \sa barrierMap()
710       \sa checkBarrier()
711    */
712    bool barrier(const Node& node) const
713    {
714      return (*_level)[node] >= _el;
715    }
716
717    /// \brief Gives back a barrier.
718    ///
719    /// This function sets \c bar to the characteristic vector of the
720    /// found barrier. \c bar should be a \ref concepts::WriteMap "writable"
721    /// node map with \c bool (or convertible) value type.
722    ///
723    /// If a feasible circulation is found, the function gives back an
724    /// empty set, so \c bar[v] will be \c false for all nodes \c v.
725    ///
726    /// \note This function calls \ref barrier() for each node,
727    /// so it runs in O(n) time.
728    ///
729    /// \pre Either \ref run() or \ref init() must be called before
730    /// using this function.
731    ///
732    /// \sa barrier()
733    /// \sa checkBarrier()
734    template<class BarrierMap>
735    void barrierMap(BarrierMap &bar) const
736    {
737      for(NodeIt n(_g);n!=INVALID;++n)
738        bar.set(n, (*_level)[n] >= _el);
739    }
740
741    /// @}
742
743    /// \name Checker Functions
744    /// The feasibility of the results can be checked using
745    /// these functions.\n
746    /// Either \ref run() or \ref start() should be called before
747    /// using them.
748
749    ///@{
750
751    ///Check if the found flow is a feasible circulation
752
753    ///Check if the found flow is a feasible circulation,
754    ///
755    bool checkFlow() const {
756      for(ArcIt e(_g);e!=INVALID;++e)
757        if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false;
758      for(NodeIt n(_g);n!=INVALID;++n)
759        {
760          Value dif=-(*_supply)[n];
761          for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e];
762          for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e];
763          if(_tol.negative(dif)) return false;
764        }
765      return true;
766    }
767
768    ///Check whether or not the last execution provides a barrier
769
770    ///Check whether or not the last execution provides a barrier.
771    ///\sa barrier()
772    ///\sa barrierMap()
773    bool checkBarrier() const
774    {
775      Value delta=0;
776      Value inf_cap = std::numeric_limits<Value>::has_infinity ?
777        std::numeric_limits<Value>::infinity() :
778        std::numeric_limits<Value>::max();
779      for(NodeIt n(_g);n!=INVALID;++n)
780        if(barrier(n))
781          delta-=(*_supply)[n];
782      for(ArcIt e(_g);e!=INVALID;++e)
783        {
784          Node s=_g.source(e);
785          Node t=_g.target(e);
786          if(barrier(s)&&!barrier(t)) {
787            if (_tol.less(inf_cap - (*_up)[e], delta)) return false;
788            delta+=(*_up)[e];
789          }
790          else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e];
791        }
792      return _tol.negative(delta);
793    }
794
795    /// @}
796
797  };
798
799}
800
801#endif
Note: See TracBrowser for help on using the repository browser.