COIN-OR::LEMON - Graph Library

source: lemon/lemon/circulation.h @ 833:e20173729589

Last change on this file since 833:e20173729589 was 833:e20173729589, checked in by Peter Kovacs <kpeter@…>, 10 years ago

Small doc fixes in several files (#331)

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 the algorithm.
461    ///
462    /// Sets the tolerance object used by the algorithm.
463    /// \return <tt>(*this)</tt>
464    Circulation& tolerance(const Tolerance& tolerance) {
465      _tol = tolerance;
466      return *this;
467    }
468
469    /// \brief Returns a const reference to the tolerance.
470    ///
471    /// Returns a const reference to the tolerance object used by
472    /// the algorithm.
473    const Tolerance& tolerance() const {
474      return _tol;
475    }
476
477    /// \name Execution Control
478    /// The simplest way to execute the algorithm is to call \ref run().\n
479    /// If you need better control on the initial solution or the execution,
480    /// you have to call one of the \ref init() functions first, then
481    /// the \ref start() function.
482
483    ///@{
484
485    /// Initializes the internal data structures.
486
487    /// Initializes the internal data structures and sets all flow values
488    /// to the lower bound.
489    void init()
490    {
491      LEMON_DEBUG(checkBoundMaps(),
492        "Upper bounds must be greater or equal to the lower bounds");
493
494      createStructures();
495
496      for(NodeIt n(_g);n!=INVALID;++n) {
497        (*_excess)[n] = (*_supply)[n];
498      }
499
500      for (ArcIt e(_g);e!=INVALID;++e) {
501        _flow->set(e, (*_lo)[e]);
502        (*_excess)[_g.target(e)] += (*_flow)[e];
503        (*_excess)[_g.source(e)] -= (*_flow)[e];
504      }
505
506      // global relabeling tested, but in general case it provides
507      // worse performance for random digraphs
508      _level->initStart();
509      for(NodeIt n(_g);n!=INVALID;++n)
510        _level->initAddItem(n);
511      _level->initFinish();
512      for(NodeIt n(_g);n!=INVALID;++n)
513        if(_tol.positive((*_excess)[n]))
514          _level->activate(n);
515    }
516
517    /// Initializes the internal data structures using a greedy approach.
518
519    /// Initializes the internal data structures using a greedy approach
520    /// to construct the initial solution.
521    void greedyInit()
522    {
523      LEMON_DEBUG(checkBoundMaps(),
524        "Upper bounds must be greater or equal to the lower bounds");
525
526      createStructures();
527
528      for(NodeIt n(_g);n!=INVALID;++n) {
529        (*_excess)[n] = (*_supply)[n];
530      }
531
532      for (ArcIt e(_g);e!=INVALID;++e) {
533        if (!_tol.less(-(*_excess)[_g.target(e)], (*_up)[e])) {
534          _flow->set(e, (*_up)[e]);
535          (*_excess)[_g.target(e)] += (*_up)[e];
536          (*_excess)[_g.source(e)] -= (*_up)[e];
537        } else if (_tol.less(-(*_excess)[_g.target(e)], (*_lo)[e])) {
538          _flow->set(e, (*_lo)[e]);
539          (*_excess)[_g.target(e)] += (*_lo)[e];
540          (*_excess)[_g.source(e)] -= (*_lo)[e];
541        } else {
542          Value fc = -(*_excess)[_g.target(e)];
543          _flow->set(e, fc);
544          (*_excess)[_g.target(e)] = 0;
545          (*_excess)[_g.source(e)] -= fc;
546        }
547      }
548
549      _level->initStart();
550      for(NodeIt n(_g);n!=INVALID;++n)
551        _level->initAddItem(n);
552      _level->initFinish();
553      for(NodeIt n(_g);n!=INVALID;++n)
554        if(_tol.positive((*_excess)[n]))
555          _level->activate(n);
556    }
557
558    ///Executes the algorithm
559
560    ///This function executes the algorithm.
561    ///
562    ///\return \c true if a feasible circulation is found.
563    ///
564    ///\sa barrier()
565    ///\sa barrierMap()
566    bool start()
567    {
568
569      Node act;
570      Node bact=INVALID;
571      Node last_activated=INVALID;
572      while((act=_level->highestActive())!=INVALID) {
573        int actlevel=(*_level)[act];
574        int mlevel=_node_num;
575        Value exc=(*_excess)[act];
576
577        for(OutArcIt e(_g,act);e!=INVALID; ++e) {
578          Node v = _g.target(e);
579          Value fc=(*_up)[e]-(*_flow)[e];
580          if(!_tol.positive(fc)) continue;
581          if((*_level)[v]<actlevel) {
582            if(!_tol.less(fc, exc)) {
583              _flow->set(e, (*_flow)[e] + exc);
584              (*_excess)[v] += exc;
585              if(!_level->active(v) && _tol.positive((*_excess)[v]))
586                _level->activate(v);
587              (*_excess)[act] = 0;
588              _level->deactivate(act);
589              goto next_l;
590            }
591            else {
592              _flow->set(e, (*_up)[e]);
593              (*_excess)[v] += fc;
594              if(!_level->active(v) && _tol.positive((*_excess)[v]))
595                _level->activate(v);
596              exc-=fc;
597            }
598          }
599          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
600        }
601        for(InArcIt e(_g,act);e!=INVALID; ++e) {
602          Node v = _g.source(e);
603          Value fc=(*_flow)[e]-(*_lo)[e];
604          if(!_tol.positive(fc)) continue;
605          if((*_level)[v]<actlevel) {
606            if(!_tol.less(fc, exc)) {
607              _flow->set(e, (*_flow)[e] - exc);
608              (*_excess)[v] += exc;
609              if(!_level->active(v) && _tol.positive((*_excess)[v]))
610                _level->activate(v);
611              (*_excess)[act] = 0;
612              _level->deactivate(act);
613              goto next_l;
614            }
615            else {
616              _flow->set(e, (*_lo)[e]);
617              (*_excess)[v] += fc;
618              if(!_level->active(v) && _tol.positive((*_excess)[v]))
619                _level->activate(v);
620              exc-=fc;
621            }
622          }
623          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
624        }
625
626        (*_excess)[act] = exc;
627        if(!_tol.positive(exc)) _level->deactivate(act);
628        else if(mlevel==_node_num) {
629          _level->liftHighestActiveToTop();
630          _el = _node_num;
631          return false;
632        }
633        else {
634          _level->liftHighestActive(mlevel+1);
635          if(_level->onLevel(actlevel)==0) {
636            _el = actlevel;
637            return false;
638          }
639        }
640      next_l:
641        ;
642      }
643      return true;
644    }
645
646    /// Runs the algorithm.
647
648    /// This function runs the algorithm.
649    ///
650    /// \return \c true if a feasible circulation is found.
651    ///
652    /// \note Apart from the return value, c.run() is just a shortcut of
653    /// the following code.
654    /// \code
655    ///   c.greedyInit();
656    ///   c.start();
657    /// \endcode
658    bool run() {
659      greedyInit();
660      return start();
661    }
662
663    /// @}
664
665    /// \name Query Functions
666    /// The results of the circulation algorithm can be obtained using
667    /// these functions.\n
668    /// Either \ref run() or \ref start() should be called before
669    /// using them.
670
671    ///@{
672
673    /// \brief Returns the flow value on the given arc.
674    ///
675    /// Returns the flow value on the given arc.
676    ///
677    /// \pre Either \ref run() or \ref init() must be called before
678    /// using this function.
679    Value flow(const Arc& arc) const {
680      return (*_flow)[arc];
681    }
682
683    /// \brief Returns a const reference to the flow map.
684    ///
685    /// Returns a const reference to the arc map storing the found flow.
686    ///
687    /// \pre Either \ref run() or \ref init() must be called before
688    /// using this function.
689    const FlowMap& flowMap() const {
690      return *_flow;
691    }
692
693    /**
694       \brief Returns \c true if the given node is in a barrier.
695
696       Barrier is a set \e B of nodes for which
697
698       \f[ \sum_{uv\in A: u\in B} upper(uv) -
699           \sum_{uv\in A: v\in B} lower(uv) < \sum_{v\in B} sup(v) \f]
700
701       holds. The existence of a set with this property prooves that a
702       feasible circualtion cannot exist.
703
704       This function returns \c true if the given node is in the found
705       barrier. If a feasible circulation is found, the function
706       gives back \c false for every node.
707
708       \pre Either \ref run() or \ref init() must be called before
709       using this function.
710
711       \sa barrierMap()
712       \sa checkBarrier()
713    */
714    bool barrier(const Node& node) const
715    {
716      return (*_level)[node] >= _el;
717    }
718
719    /// \brief Gives back a barrier.
720    ///
721    /// This function sets \c bar to the characteristic vector of the
722    /// found barrier. \c bar should be a \ref concepts::WriteMap "writable"
723    /// node map with \c bool (or convertible) value type.
724    ///
725    /// If a feasible circulation is found, the function gives back an
726    /// empty set, so \c bar[v] will be \c false for all nodes \c v.
727    ///
728    /// \note This function calls \ref barrier() for each node,
729    /// so it runs in O(n) time.
730    ///
731    /// \pre Either \ref run() or \ref init() must be called before
732    /// using this function.
733    ///
734    /// \sa barrier()
735    /// \sa checkBarrier()
736    template<class BarrierMap>
737    void barrierMap(BarrierMap &bar) const
738    {
739      for(NodeIt n(_g);n!=INVALID;++n)
740        bar.set(n, (*_level)[n] >= _el);
741    }
742
743    /// @}
744
745    /// \name Checker Functions
746    /// The feasibility of the results can be checked using
747    /// these functions.\n
748    /// Either \ref run() or \ref start() should be called before
749    /// using them.
750
751    ///@{
752
753    ///Check if the found flow is a feasible circulation
754
755    ///Check if the found flow is a feasible circulation,
756    ///
757    bool checkFlow() const {
758      for(ArcIt e(_g);e!=INVALID;++e)
759        if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false;
760      for(NodeIt n(_g);n!=INVALID;++n)
761        {
762          Value dif=-(*_supply)[n];
763          for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e];
764          for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e];
765          if(_tol.negative(dif)) return false;
766        }
767      return true;
768    }
769
770    ///Check whether or not the last execution provides a barrier
771
772    ///Check whether or not the last execution provides a barrier.
773    ///\sa barrier()
774    ///\sa barrierMap()
775    bool checkBarrier() const
776    {
777      Value delta=0;
778      Value inf_cap = std::numeric_limits<Value>::has_infinity ?
779        std::numeric_limits<Value>::infinity() :
780        std::numeric_limits<Value>::max();
781      for(NodeIt n(_g);n!=INVALID;++n)
782        if(barrier(n))
783          delta-=(*_supply)[n];
784      for(ArcIt e(_g);e!=INVALID;++e)
785        {
786          Node s=_g.source(e);
787          Node t=_g.target(e);
788          if(barrier(s)&&!barrier(t)) {
789            if (_tol.less(inf_cap - (*_up)[e], delta)) return false;
790            delta+=(*_up)[e];
791          }
792          else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e];
793        }
794      return _tol.negative(delta);
795    }
796
797    /// @}
798
799  };
800
801}
802
803#endif
Note: See TracBrowser for help on using the repository browser.