Port Circulation from svn -r3516 (#175)
authorAlpar Juttner <alpar@cs.elte.hu>
Fri, 21 Nov 2008 14:42:47 +0000
changeset 4145a7dbeaed70e
parent 399 fc6954b4fce4
child 415 fa341dd6ab23
Port Circulation from svn -r3516 (#175)
Namely,
- port the files
- apply the migrate script
- apply the unify script
- fix the compilation
- strip the demo input file
- break long lines
demo/Makefile.am
demo/circulation-input.lgf
demo/circulation_demo.cc
lemon/Makefile.am
lemon/circulation.h
     1.1 --- a/demo/Makefile.am	Fri Nov 21 10:49:39 2008 +0000
     1.2 +++ b/demo/Makefile.am	Fri Nov 21 14:42:47 2008 +0000
     1.3 @@ -1,16 +1,19 @@
     1.4  EXTRA_DIST += \
     1.5  	demo/CMakeLists.txt \
     1.6 +	demo/circulation-input.lgf \
     1.7  	demo/digraph.lgf
     1.8  
     1.9  if WANT_DEMO
    1.10  
    1.11  noinst_PROGRAMS += \
    1.12  	demo/arg_parser_demo \
    1.13 +	demo/circulation_demo \
    1.14  	demo/graph_to_eps_demo \
    1.15  	demo/lgf_demo
    1.16  
    1.17  endif WANT_DEMO
    1.18  
    1.19  demo_arg_parser_demo_SOURCES = demo/arg_parser_demo.cc
    1.20 +demo_circulation_demo_SOURCES = demo/circulation_demo.cc
    1.21  demo_graph_to_eps_demo_SOURCES = demo/graph_to_eps_demo.cc
    1.22  demo_lgf_demo_SOURCES = demo/lgf_demo.cc
     2.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     2.2 +++ b/demo/circulation-input.lgf	Fri Nov 21 14:42:47 2008 +0000
     2.3 @@ -0,0 +1,32 @@
     2.4 +@nodes 
     2.5 +coordinates_x   coordinates_y   delta   label   
     2.6 +-396.638        -311.798        0       0       
     2.7 +154.409         -214.714        13      1       
     2.8 +-378.119        -135.808        0       2       
     2.9 +-138.182        -58.0452        0       3       
    2.10 +55              -76.1018        0       4       
    2.11 +-167.302        169.88          0       5       
    2.12 +71.6876         38.7452         0       6       
    2.13 +-328.784        257.777         0       7       
    2.14 +354.242         67.9628         -13     8       
    2.15 +147             266             0       9       
    2.16 +@edges 
    2.17 +                label   lo_cap  up_cap  
    2.18 +0       1       0       0       20      
    2.19 +0       2       1       0       0       
    2.20 +1       1       2       0       3       
    2.21 +1       2       3       0       8       
    2.22 +1       3       4       0       8       
    2.23 +2       5       5       0       5       
    2.24 +3       2       6       0       5       
    2.25 +3       5       7       0       5       
    2.26 +3       6       8       0       5       
    2.27 +4       3       9       0       3       
    2.28 +5       7       10      0       3       
    2.29 +5       6       11      0       10      
    2.30 +5       8       12      0       10      
    2.31 +6       8       13      0       8       
    2.32 +8       9       14      0       20      
    2.33 +8       1       15      0       5       
    2.34 +9       5       16      0       5       
    2.35 +@end
     3.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     3.2 +++ b/demo/circulation_demo.cc	Fri Nov 21 14:42:47 2008 +0000
     3.3 @@ -0,0 +1,107 @@
     3.4 +/* -*- mode: C++; indent-tabs-mode: nil; -*-
     3.5 + *
     3.6 + * This file is a part of LEMON, a generic C++ optimization library.
     3.7 + *
     3.8 + * Copyright (C) 2003-2008
     3.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    3.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    3.11 + *
    3.12 + * Permission to use, modify and distribute this software is granted
    3.13 + * provided that this copyright notice appears in all copies. For
    3.14 + * precise terms see the accompanying LICENSE file.
    3.15 + *
    3.16 + * This software is provided "AS IS" with no warranty of any kind,
    3.17 + * express or implied, and with no claim as to its suitability for any
    3.18 + * purpose.
    3.19 + *
    3.20 + */
    3.21 +
    3.22 +///\ingroup demos
    3.23 +///\file
    3.24 +///\brief Demonstrating the usage of LEMON's General Flow algorithm
    3.25 +///
    3.26 +/// This demo program reads a general network circulation problem from the
    3.27 +/// file 'circulation-input.lgf', runs the preflow based algorithm for
    3.28 +/// finding a feasible solution and writes the output
    3.29 +/// to 'circulation-input.lgf'
    3.30 +///
    3.31 +/// \include circulation_demo.cc
    3.32 +
    3.33 +#include <iostream>
    3.34 +
    3.35 +#include <lemon/list_graph.h>
    3.36 +#include <lemon/circulation.h>
    3.37 +#include <lemon/lgf_reader.h>
    3.38 +#include <lemon/lgf_writer.h>
    3.39 +
    3.40 +using namespace lemon;
    3.41 +
    3.42 +
    3.43 +int main (int, char*[])
    3.44 +{
    3.45 +
    3.46 +    typedef ListDigraph Digraph;
    3.47 +    typedef Digraph::Node Node;
    3.48 +    typedef Digraph::NodeIt NodeIt;
    3.49 +    typedef Digraph::Arc Arc;
    3.50 +    typedef Digraph::ArcIt ArcIt;
    3.51 +    typedef Digraph::ArcMap<int> ArcMap;
    3.52 +    typedef Digraph::NodeMap<int> NodeMap;
    3.53 +    typedef Digraph::NodeMap<double> DNodeMap;
    3.54 +
    3.55 +    Digraph g;
    3.56 +    ArcMap lo(g);
    3.57 +    ArcMap up(g);
    3.58 +    NodeMap delta(g);
    3.59 +    NodeMap nid(g);
    3.60 +    ArcMap eid(g);
    3.61 +    DNodeMap cx(g);
    3.62 +    DNodeMap cy(g);
    3.63 +
    3.64 +    DigraphReader<Digraph>(g,"circulation-input.lgf").
    3.65 +      arcMap("lo_cap", lo).
    3.66 +      arcMap("up_cap", up).
    3.67 +      nodeMap("delta", delta).
    3.68 +      arcMap("label", eid).
    3.69 +      nodeMap("label", nid).
    3.70 +      nodeMap("coordinates_x", cx).
    3.71 +      nodeMap("coordinates_y", cy).
    3.72 +      run();
    3.73 +
    3.74 +    Circulation<Digraph> gen(g,lo,up,delta);
    3.75 +    bool ret=gen.run();
    3.76 +    if(ret)
    3.77 +      {
    3.78 +        std::cout << "\nA feasible flow has been found.\n";
    3.79 +        if(!gen.checkFlow()) std::cerr << "Oops!!!\n";
    3.80 +        DigraphWriter<Digraph>(g, "circulation-output.lgf").
    3.81 +          arcMap("lo_cap", lo).
    3.82 +          arcMap("up_cap", up).
    3.83 +          arcMap("flow", gen.flowMap()).
    3.84 +          nodeMap("delta", delta).
    3.85 +          arcMap("label", eid).
    3.86 +          nodeMap("label", nid).
    3.87 +          nodeMap("coordinates_x", cx).
    3.88 +          nodeMap("coordinates_y", cy).
    3.89 +          run();
    3.90 +      }
    3.91 +    else {
    3.92 +      std::cout << "\nThere is no such a flow\n";
    3.93 +      Digraph::NodeMap<int> bar(g);
    3.94 +      gen.barrierMap(bar);
    3.95 +      if(!gen.checkBarrier()) std::cerr << "Dual Oops!!!\n";
    3.96 +
    3.97 +      DigraphWriter<Digraph>(g, "circulation-output.lgf").
    3.98 +        arcMap("lo_cap", lo).
    3.99 +        arcMap("up_cap", up).
   3.100 +        nodeMap("barrier", bar).
   3.101 +        nodeMap("delta", delta).
   3.102 +        arcMap("label", eid).
   3.103 +        nodeMap("label", nid).
   3.104 +        nodeMap("coordinates_x", cx).
   3.105 +        nodeMap("coordinates_y", cy).
   3.106 +        run();
   3.107 +    }
   3.108 +  std::cout << "The output is written to file 'circulation-output.lgf'\n\n";
   3.109 +
   3.110 +}
     4.1 --- a/lemon/Makefile.am	Fri Nov 21 10:49:39 2008 +0000
     4.2 +++ b/lemon/Makefile.am	Fri Nov 21 14:42:47 2008 +0000
     4.3 @@ -20,6 +20,7 @@
     4.4  	lemon/assert.h \
     4.5          lemon/bfs.h \
     4.6          lemon/bin_heap.h \
     4.7 +        lemon/circulation.h \
     4.8          lemon/color.h \
     4.9  	lemon/concept_check.h \
    4.10          lemon/counter.h \
     5.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     5.2 +++ b/lemon/circulation.h	Fri Nov 21 14:42:47 2008 +0000
     5.3 @@ -0,0 +1,656 @@
     5.4 +/* -*- mode: C++; indent-tabs-mode: nil; -*-
     5.5 + *
     5.6 + * This file is a part of LEMON, a generic C++ optimization library.
     5.7 + *
     5.8 + * Copyright (C) 2003-2008
     5.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    5.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    5.11 + *
    5.12 + * Permission to use, modify and distribute this software is granted
    5.13 + * provided that this copyright notice appears in all copies. For
    5.14 + * precise terms see the accompanying LICENSE file.
    5.15 + *
    5.16 + * This software is provided "AS IS" with no warranty of any kind,
    5.17 + * express or implied, and with no claim as to its suitability for any
    5.18 + * purpose.
    5.19 + *
    5.20 + */
    5.21 +
    5.22 +#ifndef LEMON_CIRCULATION_H
    5.23 +#define LEMON_CIRCULATION_H
    5.24 +
    5.25 +#include <iostream>
    5.26 +#include <queue>
    5.27 +#include <lemon/tolerance.h>
    5.28 +#include <lemon/elevator.h>
    5.29 +
    5.30 +///\ingroup max_flow
    5.31 +///\file
    5.32 +///\brief Push-prelabel algorithm for finding a feasible circulation.
    5.33 +///
    5.34 +namespace lemon {
    5.35 +
    5.36 +  /// \brief Default traits class of Circulation class.
    5.37 +  ///
    5.38 +  /// Default traits class of Circulation class.
    5.39 +  /// \param _Graph Digraph type.
    5.40 +  /// \param _CapacityMap Type of capacity map.
    5.41 +  template <typename _Graph, typename _LCapMap,
    5.42 +            typename _UCapMap, typename _DeltaMap>
    5.43 +  struct CirculationDefaultTraits {
    5.44 +
    5.45 +    /// \brief The digraph type the algorithm runs on.
    5.46 +    typedef _Graph Digraph;
    5.47 +
    5.48 +    /// \brief The type of the map that stores the circulation lower
    5.49 +    /// bound.
    5.50 +    ///
    5.51 +    /// The type of the map that stores the circulation lower bound.
    5.52 +    /// It must meet the \ref concepts::ReadMap "ReadMap" concept.
    5.53 +    typedef _LCapMap LCapMap;
    5.54 +
    5.55 +    /// \brief The type of the map that stores the circulation upper
    5.56 +    /// bound.
    5.57 +    ///
    5.58 +    /// The type of the map that stores the circulation upper bound.
    5.59 +    /// It must meet the \ref concepts::ReadMap "ReadMap" concept.
    5.60 +    typedef _UCapMap UCapMap;
    5.61 +
    5.62 +    /// \brief The type of the map that stores the upper bound of
    5.63 +    /// node excess.
    5.64 +    ///
    5.65 +    /// The type of the map that stores the lower bound of node
    5.66 +    /// excess. It must meet the \ref concepts::ReadMap "ReadMap"
    5.67 +    /// concept.
    5.68 +    typedef _DeltaMap DeltaMap;
    5.69 +
    5.70 +    /// \brief The type of the length of the arcs.
    5.71 +    typedef typename DeltaMap::Value Value;
    5.72 +
    5.73 +    /// \brief The map type that stores the flow values.
    5.74 +    ///
    5.75 +    /// The map type that stores the flow values.
    5.76 +    /// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept.
    5.77 +    typedef typename Digraph::template ArcMap<Value> FlowMap;
    5.78 +
    5.79 +    /// \brief Instantiates a FlowMap.
    5.80 +    ///
    5.81 +    /// This function instantiates a \ref FlowMap.
    5.82 +    /// \param digraph The digraph, to which we would like to define
    5.83 +    /// the flow map.
    5.84 +    static FlowMap* createFlowMap(const Digraph& digraph) {
    5.85 +      return new FlowMap(digraph);
    5.86 +    }
    5.87 +
    5.88 +    /// \brief The eleavator type used by Circulation algorithm.
    5.89 +    ///
    5.90 +    /// The elevator type used by Circulation algorithm.
    5.91 +    ///
    5.92 +    /// \sa Elevator
    5.93 +    /// \sa LinkedElevator
    5.94 +    typedef lemon::Elevator<Digraph, typename Digraph::Node> Elevator;
    5.95 +
    5.96 +    /// \brief Instantiates an Elevator.
    5.97 +    ///
    5.98 +    /// This function instantiates a \ref Elevator.
    5.99 +    /// \param digraph The digraph, to which we would like to define
   5.100 +    /// the elevator.
   5.101 +    /// \param max_level The maximum level of the elevator.
   5.102 +    static Elevator* createElevator(const Digraph& digraph, int max_level) {
   5.103 +      return new Elevator(digraph, max_level);
   5.104 +    }
   5.105 +
   5.106 +    /// \brief The tolerance used by the algorithm
   5.107 +    ///
   5.108 +    /// The tolerance used by the algorithm to handle inexact computation.
   5.109 +    typedef lemon::Tolerance<Value> Tolerance;
   5.110 +
   5.111 +  };
   5.112 +
   5.113 +  ///Push-relabel algorithm for the Network Circulation Problem.
   5.114 +
   5.115 +  /**
   5.116 +     \ingroup max_flow
   5.117 +     This class implements a push-relabel algorithm
   5.118 +     or the Network Circulation Problem.
   5.119 +     The exact formulation of this problem is the following.
   5.120 +     \f[\sum_{e\in\rho(v)}x(e)-\sum_{e\in\delta(v)}x(e)\leq
   5.121 +     -delta(v)\quad \forall v\in V \f]
   5.122 +     \f[ lo(e)\leq x(e) \leq up(e) \quad \forall e\in E \f]
   5.123 +  */
   5.124 +  template<class _Graph,
   5.125 +           class _LCapMap=typename _Graph::template ArcMap<int>,
   5.126 +           class _UCapMap=_LCapMap,
   5.127 +           class _DeltaMap=typename _Graph::template NodeMap<
   5.128 +             typename _UCapMap::Value>,
   5.129 +           class _Traits=CirculationDefaultTraits<_Graph, _LCapMap,
   5.130 +                                                  _UCapMap, _DeltaMap> >
   5.131 +  class Circulation {
   5.132 +
   5.133 +    typedef _Traits Traits;
   5.134 +    typedef typename Traits::Digraph Digraph;
   5.135 +    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);
   5.136 +
   5.137 +    typedef typename Traits::Value Value;
   5.138 +
   5.139 +    typedef typename Traits::LCapMap LCapMap;
   5.140 +    typedef typename Traits::UCapMap UCapMap;
   5.141 +    typedef typename Traits::DeltaMap DeltaMap;
   5.142 +    typedef typename Traits::FlowMap FlowMap;
   5.143 +    typedef typename Traits::Elevator Elevator;
   5.144 +    typedef typename Traits::Tolerance Tolerance;
   5.145 +
   5.146 +    typedef typename Digraph::template NodeMap<Value> ExcessMap;
   5.147 +
   5.148 +    const Digraph &_g;
   5.149 +    int _node_num;
   5.150 +
   5.151 +    const LCapMap *_lo;
   5.152 +    const UCapMap *_up;
   5.153 +    const DeltaMap *_delta;
   5.154 +
   5.155 +    FlowMap *_flow;
   5.156 +    bool _local_flow;
   5.157 +
   5.158 +    Elevator* _level;
   5.159 +    bool _local_level;
   5.160 +
   5.161 +    ExcessMap* _excess;
   5.162 +
   5.163 +    Tolerance _tol;
   5.164 +    int _el;
   5.165 +
   5.166 +  public:
   5.167 +
   5.168 +    typedef Circulation Create;
   5.169 +
   5.170 +    ///\name Named template parameters
   5.171 +
   5.172 +    ///@{
   5.173 +
   5.174 +    template <typename _FlowMap>
   5.175 +    struct DefFlowMapTraits : public Traits {
   5.176 +      typedef _FlowMap FlowMap;
   5.177 +      static FlowMap *createFlowMap(const Digraph&) {
   5.178 +        LEMON_ASSERT(false, "FlowMap is not initialized");
   5.179 +        return 0; // ignore warnings
   5.180 +      }
   5.181 +    };
   5.182 +
   5.183 +    /// \brief \ref named-templ-param "Named parameter" for setting
   5.184 +    /// FlowMap type
   5.185 +    ///
   5.186 +    /// \ref named-templ-param "Named parameter" for setting FlowMap
   5.187 +    /// type
   5.188 +    template <typename _FlowMap>
   5.189 +    struct DefFlowMap
   5.190 +      : public Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
   5.191 +                           DefFlowMapTraits<_FlowMap> > {
   5.192 +      typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
   5.193 +                          DefFlowMapTraits<_FlowMap> > Create;
   5.194 +    };
   5.195 +
   5.196 +    template <typename _Elevator>
   5.197 +    struct DefElevatorTraits : public Traits {
   5.198 +      typedef _Elevator Elevator;
   5.199 +      static Elevator *createElevator(const Digraph&, int) {
   5.200 +        LEMON_ASSERT(false, "Elevator is not initialized");
   5.201 +        return 0; // ignore warnings
   5.202 +      }
   5.203 +    };
   5.204 +
   5.205 +    /// \brief \ref named-templ-param "Named parameter" for setting
   5.206 +    /// Elevator type
   5.207 +    ///
   5.208 +    /// \ref named-templ-param "Named parameter" for setting Elevator
   5.209 +    /// type
   5.210 +    template <typename _Elevator>
   5.211 +    struct DefElevator
   5.212 +      : public Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
   5.213 +                           DefElevatorTraits<_Elevator> > {
   5.214 +      typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
   5.215 +                          DefElevatorTraits<_Elevator> > Create;
   5.216 +    };
   5.217 +
   5.218 +    template <typename _Elevator>
   5.219 +    struct DefStandardElevatorTraits : public Traits {
   5.220 +      typedef _Elevator Elevator;
   5.221 +      static Elevator *createElevator(const Digraph& digraph, int max_level) {
   5.222 +        return new Elevator(digraph, max_level);
   5.223 +      }
   5.224 +    };
   5.225 +
   5.226 +    /// \brief \ref named-templ-param "Named parameter" for setting
   5.227 +    /// Elevator type
   5.228 +    ///
   5.229 +    /// \ref named-templ-param "Named parameter" for setting Elevator
   5.230 +    /// type. The Elevator should be standard constructor interface, ie.
   5.231 +    /// the digraph and the maximum level should be passed to it.
   5.232 +    template <typename _Elevator>
   5.233 +    struct DefStandardElevator
   5.234 +      : public Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
   5.235 +                       DefStandardElevatorTraits<_Elevator> > {
   5.236 +      typedef Circulation<Digraph, LCapMap, UCapMap, DeltaMap,
   5.237 +                      DefStandardElevatorTraits<_Elevator> > Create;
   5.238 +    };
   5.239 +
   5.240 +    /// @}
   5.241 +
   5.242 +  protected:
   5.243 +
   5.244 +    Circulation() {}
   5.245 +
   5.246 +  public:
   5.247 +
   5.248 +    /// The constructor of the class.
   5.249 +
   5.250 +    /// The constructor of the class.
   5.251 +    /// \param g The digraph the algorithm runs on.
   5.252 +    /// \param lo The lower bound capacity of the arcs.
   5.253 +    /// \param up The upper bound capacity of the arcs.
   5.254 +    /// \param delta The lower bound on node excess.
   5.255 +    Circulation(const Digraph &g,const LCapMap &lo,
   5.256 +                const UCapMap &up,const DeltaMap &delta)
   5.257 +      : _g(g), _node_num(),
   5.258 +        _lo(&lo),_up(&up),_delta(&delta),_flow(0),_local_flow(false),
   5.259 +        _level(0), _local_level(false), _excess(0), _el() {}
   5.260 +
   5.261 +    /// Destrcutor.
   5.262 +    ~Circulation() {
   5.263 +      destroyStructures();
   5.264 +    }
   5.265 +
   5.266 +  private:
   5.267 +
   5.268 +    void createStructures() {
   5.269 +      _node_num = _el = countNodes(_g);
   5.270 +
   5.271 +      if (!_flow) {
   5.272 +        _flow = Traits::createFlowMap(_g);
   5.273 +        _local_flow = true;
   5.274 +      }
   5.275 +      if (!_level) {
   5.276 +        _level = Traits::createElevator(_g, _node_num);
   5.277 +        _local_level = true;
   5.278 +      }
   5.279 +      if (!_excess) {
   5.280 +        _excess = new ExcessMap(_g);
   5.281 +      }
   5.282 +    }
   5.283 +
   5.284 +    void destroyStructures() {
   5.285 +      if (_local_flow) {
   5.286 +        delete _flow;
   5.287 +      }
   5.288 +      if (_local_level) {
   5.289 +        delete _level;
   5.290 +      }
   5.291 +      if (_excess) {
   5.292 +        delete _excess;
   5.293 +      }
   5.294 +    }
   5.295 +
   5.296 +  public:
   5.297 +
   5.298 +    /// Sets the lower bound capacity map.
   5.299 +
   5.300 +    /// Sets the lower bound capacity map.
   5.301 +    /// \return \c (*this)
   5.302 +    Circulation& lowerCapMap(const LCapMap& map) {
   5.303 +      _lo = &map;
   5.304 +      return *this;
   5.305 +    }
   5.306 +
   5.307 +    /// Sets the upper bound capacity map.
   5.308 +
   5.309 +    /// Sets the upper bound capacity map.
   5.310 +    /// \return \c (*this)
   5.311 +    Circulation& upperCapMap(const LCapMap& map) {
   5.312 +      _up = &map;
   5.313 +      return *this;
   5.314 +    }
   5.315 +
   5.316 +    /// Sets the lower bound map on excess.
   5.317 +
   5.318 +    /// Sets the lower bound map on excess.
   5.319 +    /// \return \c (*this)
   5.320 +    Circulation& deltaMap(const DeltaMap& map) {
   5.321 +      _delta = &map;
   5.322 +      return *this;
   5.323 +    }
   5.324 +
   5.325 +    /// Sets the flow map.
   5.326 +
   5.327 +    /// Sets the flow map.
   5.328 +    /// \return \c (*this)
   5.329 +    Circulation& flowMap(FlowMap& map) {
   5.330 +      if (_local_flow) {
   5.331 +        delete _flow;
   5.332 +        _local_flow = false;
   5.333 +      }
   5.334 +      _flow = &map;
   5.335 +      return *this;
   5.336 +    }
   5.337 +
   5.338 +    /// Returns the flow map.
   5.339 +
   5.340 +    /// \return The flow map.
   5.341 +    ///
   5.342 +    const FlowMap& flowMap() {
   5.343 +      return *_flow;
   5.344 +    }
   5.345 +
   5.346 +    /// Sets the elevator.
   5.347 +
   5.348 +    /// Sets the elevator.
   5.349 +    /// \return \c (*this)
   5.350 +    Circulation& elevator(Elevator& elevator) {
   5.351 +      if (_local_level) {
   5.352 +        delete _level;
   5.353 +        _local_level = false;
   5.354 +      }
   5.355 +      _level = &elevator;
   5.356 +      return *this;
   5.357 +    }
   5.358 +
   5.359 +    /// Returns the elevator.
   5.360 +
   5.361 +    /// \return The elevator.
   5.362 +    ///
   5.363 +    const Elevator& elevator() {
   5.364 +      return *_level;
   5.365 +    }
   5.366 +
   5.367 +    /// Sets the tolerance used by algorithm.
   5.368 +
   5.369 +    /// Sets the tolerance used by algorithm.
   5.370 +    ///
   5.371 +    Circulation& tolerance(const Tolerance& tolerance) const {
   5.372 +      _tol = tolerance;
   5.373 +      return *this;
   5.374 +    }
   5.375 +
   5.376 +    /// Returns the tolerance used by algorithm.
   5.377 +
   5.378 +    /// Returns the tolerance used by algorithm.
   5.379 +    ///
   5.380 +    const Tolerance& tolerance() const {
   5.381 +      return tolerance;
   5.382 +    }
   5.383 +
   5.384 +    /// \name Execution control
   5.385 +    /// The simplest way to execute the algorithm is to use one of the
   5.386 +    /// member functions called \c run().
   5.387 +    /// \n
   5.388 +    /// If you need more control on initial solution or execution then
   5.389 +    /// you have to call one \ref init() function and then the start()
   5.390 +    /// function.
   5.391 +
   5.392 +    ///@{
   5.393 +
   5.394 +    /// Initializes the internal data structures.
   5.395 +
   5.396 +    /// Initializes the internal data structures. This function sets
   5.397 +    /// all flow values to the lower bound.
   5.398 +    /// \return This function returns false if the initialization
   5.399 +    /// process found a barrier.
   5.400 +    void init()
   5.401 +    {
   5.402 +      createStructures();
   5.403 +
   5.404 +      for(NodeIt n(_g);n!=INVALID;++n) {
   5.405 +        _excess->set(n, (*_delta)[n]);
   5.406 +      }
   5.407 +
   5.408 +      for (ArcIt e(_g);e!=INVALID;++e) {
   5.409 +        _flow->set(e, (*_lo)[e]);
   5.410 +        _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_flow)[e]);
   5.411 +        _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_flow)[e]);
   5.412 +      }
   5.413 +
   5.414 +      // global relabeling tested, but in general case it provides
   5.415 +      // worse performance for random digraphs
   5.416 +      _level->initStart();
   5.417 +      for(NodeIt n(_g);n!=INVALID;++n)
   5.418 +        _level->initAddItem(n);
   5.419 +      _level->initFinish();
   5.420 +      for(NodeIt n(_g);n!=INVALID;++n)
   5.421 +        if(_tol.positive((*_excess)[n]))
   5.422 +          _level->activate(n);
   5.423 +    }
   5.424 +
   5.425 +    /// Initializes the internal data structures.
   5.426 +
   5.427 +    /// Initializes the internal data structures. This functions uses
   5.428 +    /// greedy approach to construct the initial solution.
   5.429 +    void greedyInit()
   5.430 +    {
   5.431 +      createStructures();
   5.432 +
   5.433 +      for(NodeIt n(_g);n!=INVALID;++n) {
   5.434 +        _excess->set(n, (*_delta)[n]);
   5.435 +      }
   5.436 +
   5.437 +      for (ArcIt e(_g);e!=INVALID;++e) {
   5.438 +        if (!_tol.positive((*_excess)[_g.target(e)] + (*_up)[e])) {
   5.439 +          _flow->set(e, (*_up)[e]);
   5.440 +          _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_up)[e]);
   5.441 +          _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_up)[e]);
   5.442 +        } else if (_tol.positive((*_excess)[_g.target(e)] + (*_lo)[e])) {
   5.443 +          _flow->set(e, (*_lo)[e]);
   5.444 +          _excess->set(_g.target(e), (*_excess)[_g.target(e)] + (*_lo)[e]);
   5.445 +          _excess->set(_g.source(e), (*_excess)[_g.source(e)] - (*_lo)[e]);
   5.446 +        } else {
   5.447 +          Value fc = -(*_excess)[_g.target(e)];
   5.448 +          _flow->set(e, fc);
   5.449 +          _excess->set(_g.target(e), 0);
   5.450 +          _excess->set(_g.source(e), (*_excess)[_g.source(e)] - fc);
   5.451 +        }
   5.452 +      }
   5.453 +
   5.454 +      _level->initStart();
   5.455 +      for(NodeIt n(_g);n!=INVALID;++n)
   5.456 +        _level->initAddItem(n);
   5.457 +      _level->initFinish();
   5.458 +      for(NodeIt n(_g);n!=INVALID;++n)
   5.459 +        if(_tol.positive((*_excess)[n]))
   5.460 +          _level->activate(n);
   5.461 +    }
   5.462 +
   5.463 +    ///Starts the algorithm
   5.464 +
   5.465 +    ///This function starts the algorithm.
   5.466 +    ///\return This function returns true if it found a feasible circulation.
   5.467 +    ///
   5.468 +    ///\sa barrier()
   5.469 +    bool start()
   5.470 +    {
   5.471 +
   5.472 +      Node act;
   5.473 +      Node bact=INVALID;
   5.474 +      Node last_activated=INVALID;
   5.475 +      while((act=_level->highestActive())!=INVALID) {
   5.476 +        int actlevel=(*_level)[act];
   5.477 +        int mlevel=_node_num;
   5.478 +        Value exc=(*_excess)[act];
   5.479 +
   5.480 +        for(OutArcIt e(_g,act);e!=INVALID; ++e) {
   5.481 +          Node v = _g.target(e);
   5.482 +          Value fc=(*_up)[e]-(*_flow)[e];
   5.483 +          if(!_tol.positive(fc)) continue;
   5.484 +          if((*_level)[v]<actlevel) {
   5.485 +            if(!_tol.less(fc, exc)) {
   5.486 +              _flow->set(e, (*_flow)[e] + exc);
   5.487 +              _excess->set(v, (*_excess)[v] + exc);
   5.488 +              if(!_level->active(v) && _tol.positive((*_excess)[v]))
   5.489 +                _level->activate(v);
   5.490 +              _excess->set(act,0);
   5.491 +              _level->deactivate(act);
   5.492 +              goto next_l;
   5.493 +            }
   5.494 +            else {
   5.495 +              _flow->set(e, (*_up)[e]);
   5.496 +              _excess->set(v, (*_excess)[v] + fc);
   5.497 +              if(!_level->active(v) && _tol.positive((*_excess)[v]))
   5.498 +                _level->activate(v);
   5.499 +              exc-=fc;
   5.500 +            }
   5.501 +          }
   5.502 +          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
   5.503 +        }
   5.504 +        for(InArcIt e(_g,act);e!=INVALID; ++e) {
   5.505 +          Node v = _g.source(e);
   5.506 +          Value fc=(*_flow)[e]-(*_lo)[e];
   5.507 +          if(!_tol.positive(fc)) continue;
   5.508 +          if((*_level)[v]<actlevel) {
   5.509 +            if(!_tol.less(fc, exc)) {
   5.510 +              _flow->set(e, (*_flow)[e] - exc);
   5.511 +              _excess->set(v, (*_excess)[v] + exc);
   5.512 +              if(!_level->active(v) && _tol.positive((*_excess)[v]))
   5.513 +                _level->activate(v);
   5.514 +              _excess->set(act,0);
   5.515 +              _level->deactivate(act);
   5.516 +              goto next_l;
   5.517 +            }
   5.518 +            else {
   5.519 +              _flow->set(e, (*_lo)[e]);
   5.520 +              _excess->set(v, (*_excess)[v] + fc);
   5.521 +              if(!_level->active(v) && _tol.positive((*_excess)[v]))
   5.522 +                _level->activate(v);
   5.523 +              exc-=fc;
   5.524 +            }
   5.525 +          }
   5.526 +          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
   5.527 +        }
   5.528 +
   5.529 +        _excess->set(act, exc);
   5.530 +        if(!_tol.positive(exc)) _level->deactivate(act);
   5.531 +        else if(mlevel==_node_num) {
   5.532 +          _level->liftHighestActiveToTop();
   5.533 +          _el = _node_num;
   5.534 +          return false;
   5.535 +        }
   5.536 +        else {
   5.537 +          _level->liftHighestActive(mlevel+1);
   5.538 +          if(_level->onLevel(actlevel)==0) {
   5.539 +            _el = actlevel;
   5.540 +            return false;
   5.541 +          }
   5.542 +        }
   5.543 +      next_l:
   5.544 +        ;
   5.545 +      }
   5.546 +      return true;
   5.547 +    }
   5.548 +
   5.549 +    /// Runs the circulation algorithm.
   5.550 +
   5.551 +    /// Runs the circulation algorithm.
   5.552 +    /// \note fc.run() is just a shortcut of the following code.
   5.553 +    /// \code
   5.554 +    ///   fc.greedyInit();
   5.555 +    ///   return fc.start();
   5.556 +    /// \endcode
   5.557 +    bool run() {
   5.558 +      greedyInit();
   5.559 +      return start();
   5.560 +    }
   5.561 +
   5.562 +    /// @}
   5.563 +
   5.564 +    /// \name Query Functions
   5.565 +    /// The result of the %Circulation algorithm can be obtained using
   5.566 +    /// these functions.
   5.567 +    /// \n
   5.568 +    /// Before the use of these functions,
   5.569 +    /// either run() or start() must be called.
   5.570 +
   5.571 +    ///@{
   5.572 +
   5.573 +    /**
   5.574 +       \brief Returns a barrier
   5.575 +       
   5.576 +       Barrier is a set \e B of nodes for which
   5.577 +       \f[ \sum_{v\in B}-delta(v)<
   5.578 +       \sum_{e\in\rho(B)}lo(e)-\sum_{e\in\delta(B)}up(e) \f]
   5.579 +       holds. The existence of a set with this property prooves that a feasible
   5.580 +       flow cannot exists.
   5.581 +       \sa checkBarrier()
   5.582 +       \sa run()
   5.583 +    */
   5.584 +    template<class GT>
   5.585 +    void barrierMap(GT &bar)
   5.586 +    {
   5.587 +      for(NodeIt n(_g);n!=INVALID;++n)
   5.588 +        bar.set(n, (*_level)[n] >= _el);
   5.589 +    }
   5.590 +
   5.591 +    ///Returns true if the node is in the barrier
   5.592 +
   5.593 +    ///Returns true if the node is in the barrier
   5.594 +    ///\sa barrierMap()
   5.595 +    bool barrier(const Node& node)
   5.596 +    {
   5.597 +      return (*_level)[node] >= _el;
   5.598 +    }
   5.599 +
   5.600 +    /// \brief Returns the flow on the arc.
   5.601 +    ///
   5.602 +    /// Sets the \c flowMap to the flow on the arcs. This method can
   5.603 +    /// be called after the second phase of algorithm.
   5.604 +    Value flow(const Arc& arc) const {
   5.605 +      return (*_flow)[arc];
   5.606 +    }
   5.607 +
   5.608 +    /// @}
   5.609 +
   5.610 +    /// \name Checker Functions
   5.611 +    /// The feasibility  of the results can be checked using
   5.612 +    /// these functions.
   5.613 +    /// \n
   5.614 +    /// Before the use of these functions,
   5.615 +    /// either run() or start() must be called.
   5.616 +
   5.617 +    ///@{
   5.618 +
   5.619 +    ///Check if the  \c flow is a feasible circulation
   5.620 +    bool checkFlow() {
   5.621 +      for(ArcIt e(_g);e!=INVALID;++e)
   5.622 +        if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false;
   5.623 +      for(NodeIt n(_g);n!=INVALID;++n)
   5.624 +        {
   5.625 +          Value dif=-(*_delta)[n];
   5.626 +          for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e];
   5.627 +          for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e];
   5.628 +          if(_tol.negative(dif)) return false;
   5.629 +        }
   5.630 +      return true;
   5.631 +    }
   5.632 +
   5.633 +    ///Check whether or not the last execution provides a barrier
   5.634 +
   5.635 +    ///Check whether or not the last execution provides a barrier
   5.636 +    ///\sa barrier()
   5.637 +    bool checkBarrier()
   5.638 +    {
   5.639 +      Value delta=0;
   5.640 +      for(NodeIt n(_g);n!=INVALID;++n)
   5.641 +        if(barrier(n))
   5.642 +          delta-=(*_delta)[n];
   5.643 +      for(ArcIt e(_g);e!=INVALID;++e)
   5.644 +        {
   5.645 +          Node s=_g.source(e);
   5.646 +          Node t=_g.target(e);
   5.647 +          if(barrier(s)&&!barrier(t)) delta+=(*_up)[e];
   5.648 +          else if(barrier(t)&&!barrier(s)) delta-=(*_lo)[e];
   5.649 +        }
   5.650 +      return _tol.negative(delta);
   5.651 +    }
   5.652 +
   5.653 +    /// @}
   5.654 +
   5.655 +  };
   5.656 +
   5.657 +}
   5.658 +
   5.659 +#endif