# source:lemon-0.x/lemon/circulation.h@2391:14a343be7a5a

Last change on this file since 2391:14a343be7a5a was 2391:14a343be7a5a, checked in by Alpar Juttner, 12 years ago

Happy New Year to all source files!

File size: 7.8 KB
Line
1/* -*- C++ -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library
4 *
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/graph_utils.h>
23#include <iostream>
24#include <queue>
25#include <lemon/tolerance.h>
26#include <lemon/elevator.h>
27
28///\ingroup max_flow
29///\file
30///\brief Push-prelabel algorithm for finding a feasible circulation.
31///
32namespace lemon {
33
34  ///Preflow algorithm for the Network Circulation Problem.
35
36  ///\ingroup max_flow
37  ///This class implements a preflow algorithm
38  ///for the Network Circulation Problem.
39  ///The exact formulation of this problem is the following.
40  /// \f[\sum_{e\in\rho(v)}x(e)-\sum_{e\in\delta(v)}x(e)\leq delta(v)\quad \forall v\in V \f]
41  /// \f[ lo(e)\leq x(e) \leq up(e) \quad \forall e\in E \f]
42  ///
43  template<class Graph,
44           class Value,
45           class FlowMap=typename Graph::template EdgeMap<Value>,
46           class LCapMap=typename Graph::template EdgeMap<Value>,
47           class UCapMap=LCapMap,
48           class DeltaMap=typename Graph::template NodeMap<Value>
49            >
50  class Circulation {
51    typedef typename Graph::Node Node;
52    typedef typename Graph::NodeIt NodeIt;
53    typedef typename Graph::Edge Edge;
54    typedef typename Graph::EdgeIt EdgeIt;
55    typedef typename Graph::InEdgeIt InEdgeIt;
56    typedef typename Graph::OutEdgeIt OutEdgeIt;
57
58
59    const Graph &_g;
60    int _node_num;
61
62    const LCapMap &_lo;
63    const UCapMap &_up;
64    const DeltaMap &_delta;
65    FlowMap &_x;
66    Tolerance<Value> _tol;
67    Elevator<Graph,typename Graph::Node> _levels;
68    typename Graph::template NodeMap<Value> _excess;
69
70  public:
71    ///\e
72    Circulation(const Graph &g,const LCapMap &lo,const UCapMap &up,
73                const DeltaMap &delta,FlowMap &x) :
74      _g(g),
75      _node_num(countNodes(g)),
76      _lo(lo),_up(up),_delta(delta),_x(x),
77      _levels(g,_node_num),
78      _excess(g)
79    {
80    }
81
82  private:
83
85    {
86      if(_tol.positive(_excess[n]+=v))
87        {
88          if(!_levels.active(n)) _levels.activate(n);
89        }
90      else if(_levels.active(n)) _levels.deactivate(n);
91    }
92
93    void init()
94    {
95
96      _x=_lo;
97
98      for(NodeIt n(_g);n!=INVALID;++n) _excess[n]=-_delta[n];
99
100      for(EdgeIt e(_g);e!=INVALID;++e)
101        {
102          _excess[_g.target(e)]+=_x[e];
103          _excess[_g.source(e)]-=_x[e];
104        }
105
106      _levels.initStart();
107      for(NodeIt n(_g);n!=INVALID;++n)
109      _levels.initFinish();
110      for(NodeIt n(_g);n!=INVALID;++n)
111        if(_tol.positive(_excess[n]))
112          _levels.activate(n);
113    }
114
115  public:
116    ///Check if \c x is a feasible circulation
117    template<class FT>
118    bool checkX(FT &x) {
119      for(EdgeIt e(_g);e!=INVALID;++e)
120        if(x[e]<_lo[e]||x[e]>_up[e]) return false;
121      for(NodeIt n(_g);n!=INVALID;++n)
122        {
123          Value dif=_delta[n];
124          for(InEdgeIt e(_g,n);e!=INVALID;++e) dif-=x[e];
125          for(OutEdgeIt e(_g,n);e!=INVALID;++e) dif+=x[e];
126          if(_tol.negative(dif)) return false;
127        }
128      return true;
129    };
130    ///Check if the default \c x is a feasible circulation
131    bool checkX() { return checkX(_x); }
132
133    ///Chech if \c bar is a real barrier
134
135    ///Chech if \c bar is a real barrier
136    ///\sa barrier()
137    template<class GT>
138    bool checkBarrier(GT &bar)
139    {
140      Value delta=0;
141      for(NodeIt n(_g);n!=INVALID;++n)
142        if(bar[n])
143          delta+=_delta[n];
144      for(EdgeIt e(_g);e!=INVALID;++e)
145        {
146          Node s=_g.source(e);
147          Node t=_g.target(e);
148          if(bar[s]&&!bar[t]) delta+=_up[e];
149          else if(bar[t]&&!bar[s]) delta-=_lo[e];
150        }
151      return _tol.negative(delta);
152    }
153    ///Chech whether or not the last execution provides a barrier
154
155    ///Chech whether or not the last execution provides a barrier
156    ///\sa barrier()
157    bool checkBarrier()
158    {
159      typename Graph:: template NodeMap<bool> bar(_g);
160      barrier(bar);
161      return checkBarrier(bar);
162    }
163    ///Run the algorithm
164
165    ///This function runs the algorithm.
166    ///\return This function returns -1 if it found a feasible circulation.
167    /// nonnegative values (including 0) mean that no feasible solution is
168    /// found. In this case the return value means an "empty level".
169    ///
170    ///\sa barrier()
171    int run()
172    {
173      init();
174
175#ifdef LEMON_CIRCULATION_DEBUG
176      for(NodeIt n(_g);n!=INVALID;++n)
177        std::cerr<< _levels[n] << ' ';
178      std::cerr << std::endl;
179#endif
180      Node act;
181      Node bact=INVALID;
182      Node last_activated=INVALID;
183      while((act=_levels.highestActive())!=INVALID) {
184        int actlevel=_levels[act];
185        int tlevel;
186        int mlevel=_node_num;
187        Value exc=_excess[act];
188        Value fc;
189
190#ifdef LEMON_CIRCULATION_DEBUG
191        for(NodeIt n(_g);n!=INVALID;++n)
192          std::cerr<< _levels[n] << ' ';
193        std::cerr << std::endl;
194        std::cerr << "Process node " << _g.id(act)
195                  << " on level " << actlevel
196                  << " with excess " << exc
197                  << std::endl;
198#endif
199        for(OutEdgeIt e(_g,act);e!=INVALID; ++e)
200          if((fc=_up[e]-_x[e])>0)
201            if((tlevel=_levels[_g.target(e)])<actlevel)
202              if(fc<=exc) {
203                _x[e]=_up[e];
205                exc-=fc;
206#ifdef LEMON_CIRCULATION_DEBUG
207                std::cerr << "  Push " << fc
208                          << " toward " << _g.id(_g.target(e)) << std::endl;
209#endif
210              }
211              else {
212                _x[e]+=exc;
214                //exc=0;
215                _excess[act]=0;
216                _levels.deactivate(act);
217#ifdef LEMON_CIRCULATION_DEBUG
218                std::cerr << "  Push " << exc
219                          << " toward " << _g.id(_g.target(e)) << std::endl;
220                std::cerr << "  Deactivate\n";
221#endif
222                goto next_l;
223              }
224            else if(tlevel<mlevel) mlevel=tlevel;
225
226        for(InEdgeIt e(_g,act);e!=INVALID; ++e)
227          if((fc=_x[e]-_lo[e])>0)
228            if((tlevel=_levels[_g.source(e)])<actlevel)
229              if(fc<=exc) {
230                _x[e]=_lo[e];
232                exc-=fc;
233#ifdef LEMON_CIRCULATION_DEBUG
234                std::cerr << "  Push " << fc
235                          << " toward " << _g.id(_g.source(e)) << std::endl;
236#endif
237              }
238              else {
239                _x[e]-=exc;
241                //exc=0;
242                _excess[act]=0;
243                _levels.deactivate(act);
244#ifdef LEMON_CIRCULATION_DEBUG
245                std::cerr << "  Push " << exc
246                          << " toward " << _g.id(_g.source(e)) << std::endl;
247                std::cerr << "  Deactivate\n";
248#endif
249                goto next_l;
250              }
251            else if(tlevel<mlevel) mlevel=tlevel;
252
253        _excess[act]=exc;
254        if(!_tol.positive(exc)) _levels.deactivate(act);
255        else if(mlevel==_node_num) {
256          _levels.liftHighestActiveTo(_node_num);
257#ifdef LEMON_CIRCULATION_DEBUG
258          std::cerr << "  Lift to level " << _node_num << std::endl;
259#endif
260          return _levels.onLevel(_node_num-1)==0?_node_num-1:actlevel;
261        }
262        else {
263          _levels.liftHighestActiveTo(mlevel+1);
264#ifdef LEMON_CIRCULATION_DEBUG
265          std::cerr << "  Lift to level " << mlevel+1 << std::endl;
266#endif
267          if(_levels.onLevel(actlevel)==0)
268            return actlevel;
269        }
270      next_l:
271        ;
272      }
273#ifdef LEMON_CIRCULATION_DEBUG
274      std::cerr << "Feasible flow found.\n";
275#endif
276      return -1;
277    }
278
279    ///Return a barrier
280
281    ///Barrier is a set \e B of nodes for which
282    /// \f[ \sum_{v\in B}delta(v)<\sum_{e\in\rho(B)}lo(e)-\sum_{e\in\delta(B)}up(e) \f]
283    ///holds. The existence of a set with this property prooves that a feasible
284    ///flow cannot exists.
285    ///\pre The run() must have been executed, and its return value was -1.
286    ///\sa checkBarrier()
287    ///\sa run()
288    template<class GT>
289    void barrier(GT &bar,int empty_level=-1)
290    {
291      if(empty_level==-1)
292        for(empty_level=0;_levels.onLevel(empty_level);empty_level++) ;
293      for(NodeIt n(_g);n!=INVALID;++n)
294        bar[n] = _levels[n]>empty_level;
295    }
296
297  };
298
299}
300
301#endif
Note: See TracBrowser for help on using the repository browser.