COIN-OR::LEMON - Graph Library

source: lemon-0.x/src/work/jacint/preflow_push_max_flow.h @ 96:e2e18eb0fd10

Last change on this file since 96:e2e18eb0fd10 was 85:15362fafaf1a, checked in by jacint, 17 years ago

* empty log message *

File size: 7.5 KB
Line 
1/*
2preflow_push_max_flow_h
3by jacint.
4Runs a preflow push algorithm with the modification,
5that we do not push on nodes with level at least n.
6Moreover, if a level gets empty, we set all nodes above that
7level to level n. Hence, in the end, we arrive at a maximum preflow
8with value of a max flow value. An empty level gives a minimum cut.
9
10Member functions:
11
12void run() : runs the algorithm
13
14  The following functions should be used after run() was already run.
15
16T maxflow() : returns the value of a maximum flow
17
18NodeMap<bool> mincut(): returns a
19     characteristic vector of a minimum cut.
20*/
21
22#ifndef PREFLOW_PUSH_MAX_FLOW_H
23#define PREFLOW_PUSH_MAX_FLOW_H
24
25#include <algorithm>
26#include <vector>
27#include <stack>
28
29#include <reverse_bfs.h>
30
31
32namespace marci {
33
34  template <typename Graph, typename T>
35  class preflow_push_max_flow {
36   
37    typedef typename Graph::NodeIt NodeIt;
38    typedef typename Graph::EachNodeIt EachNodeIt;
39    typedef typename Graph::OutEdgeIt OutEdgeIt;
40    typedef typename Graph::InEdgeIt InEdgeIt;
41   
42    Graph& G;
43    NodeIt s;
44    NodeIt t;
45    typename Graph::EdgeMap<T>& capacity;
46    T value;
47    typename Graph::NodeMap<bool> mincutvector;   
48
49
50     
51  public:
52       
53    preflow_push_max_flow ( Graph& _G, NodeIt _s, NodeIt _t,
54                            typename Graph::EdgeMap<T>& _capacity) :
55      G(_G), s(_s), t(_t), capacity(_capacity), mincutvector(_G, false) { }
56
57
58    /*
59      The run() function runs a modified version of the
60      highest label preflow-push, which only
61      finds a maximum preflow, hence giving the value of a maximum flow.
62    */
63    void run() {
64 
65      typename Graph::EdgeMap<T> flow(G, 0);
66      typename Graph::NodeMap<int> level(G);   
67      typename Graph::NodeMap<T> excess(G);   
68           
69      int n=G.nodeNum();                       
70      int b=n-2;
71      /*
72        b is a bound on the highest level of an active Node.
73        In the beginning it is at most n-2.
74      */
75     
76      std::vector<int> numb(n);     //The number of Nodes on level i < n.
77      std::vector<std::stack<NodeIt> > stack(2*n-1);   
78      //Stack of the active Nodes in level i.
79
80      /*Reverse_bfs from t, to find the starting level.*/
81      reverse_bfs<Graph> bfs(G, t);
82      bfs.run();
83      for(EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v)
84        {
85          int dist=bfs.dist(v);
86          level.set(v, dist);
87          ++numb[dist];
88        }
89
90      level.set(s,n);
91
92      /* Starting flow. It is everywhere 0 at the moment. */
93      for(OutEdgeIt e=G.template first<OutEdgeIt>(s); e.valid(); ++e)
94        {
95          if ( capacity.get(e) > 0 ) {
96            NodeIt w=G.head(e);
97            flow.set(e, capacity.get(e));
98            stack[level.get(w)].push(w);
99            excess.set(w, excess.get(w)+capacity.get(e));
100          }
101        }
102
103      /*
104         End of preprocessing
105      */
106
107
108
109      /*
110        Push/relabel on the highest level active Nodes.
111      */
112       
113      /*While there exists an active Node.*/
114      while (b) {
115
116        /*We decrease the bound if there is no active node of level b.*/
117        if (stack[b].empty()) {
118          --b;
119        } else {
120
121          NodeIt w=stack[b].top();    //w is the highest label active Node.
122          stack[b].pop();                    //We delete w from the stack.
123       
124          int newlevel=2*n-2;                //In newlevel we maintain the next level of w.
125       
126          for(OutEdgeIt e=G.template first<OutEdgeIt>(w); e.valid(); ++e) {
127            NodeIt v=G.head(e);
128            /*e is the Edge wv.*/
129
130            if (flow.get(e)<capacity.get(e)) {             
131              /*e is an Edge of the residual graph */
132
133              if(level.get(w)==level.get(v)+1) {     
134                /*Push is allowed now*/
135
136                if (capacity.get(e)-flow.get(e) > excess.get(w)) {       
137                  /*A nonsaturating push.*/
138                 
139                  if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
140                  /*v becomes active.*/
141                 
142                  flow.set(e, flow.get(e)+excess.get(w));
143                  excess.set(v, excess.get(v)+excess.get(w));
144                  excess.set(w,0);
145                  //std::cout << w << " " << v <<" elore elen nonsat pump "  << std::endl;
146                  break;
147                } else {
148                  /*A saturating push.*/
149
150                  if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
151                  /*v becomes active.*/
152
153                  excess.set(v, excess.get(v)+capacity.get(e)-flow.get(e));
154                  excess.set(w, excess.get(w)-capacity.get(e)+flow.get(e));
155                  flow.set(e, capacity.get(e));
156                  //std::cout << w <<" " << v <<" elore elen sat pump "   << std::endl;
157                  if (excess.get(w)==0) break;
158                  /*If w is not active any more, then we go on to the next Node.*/
159                 
160                } // if (capacity.get(e)-flow.get(e) > excess.get(w))
161              } // if (level.get(w)==level.get(v)+1)
162           
163              else {newlevel = newlevel < level.get(v) ? newlevel : level.get(v);}
164           
165            } //if (flow.get(e)<capacity.get(e))
166         
167          } //for(OutEdgeIt e=G.first_OutEdge(w); e.valid(); ++e)
168         
169
170
171          for(InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) {
172            NodeIt v=G.tail(e);
173            /*e is the Edge vw.*/
174
175            if (excess.get(w)==0) break;
176            /*It may happen, that w became inactive in the first 'for' cycle.*/         
177 
178            if(flow.get(e)>0) {             
179              /*e is an Edge of the residual graph */
180
181              if(level.get(w)==level.get(v)+1) { 
182                /*Push is allowed now*/
183               
184                if (flow.get(e) > excess.get(w)) {
185                  /*A nonsaturating push.*/
186                 
187                  if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
188                  /*v becomes active.*/
189
190                  flow.set(e, flow.get(e)-excess.get(w));
191                  excess.set(v, excess.get(v)+excess.get(w));
192                  excess.set(w,0);
193                  //std::cout << v << " " << w << " vissza elen nonsat pump "     << std::endl;
194                  break;
195                } else {                                               
196                  /*A saturating push.*/
197                 
198                  if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
199                  /*v becomes active.*/
200                 
201                  flow.set(e,0);
202                  excess.set(v, excess.get(v)+flow.get(e));
203                  excess.set(w, excess.get(w)-flow.get(e));
204                  //std::cout << v <<" " << w << " vissza elen sat pump "     << std::endl;
205                  if (excess.get(w)==0) { break;}
206                } //if (flow.get(e) > excess.get(v))
207              } //if(level.get(w)==level.get(v)+1)
208             
209              else {newlevel = newlevel < level.get(v) ? newlevel : level.get(v);}
210              //std::cout << "Leveldecrease of Node " << w << " to " << newlevel << std::endl;
211
212            } //if (flow.get(e)>0)
213
214          } //for in-Edge
215
216
217
218
219          /*
220            Relabel
221          */
222          if (excess.get(w)>0) {
223            /*Now newlevel <= n*/
224
225            int l=level.get(w);         //l is the old level of w.
226            --numb[l];
227           
228            if (newlevel == n) {
229              level.set(w,n);
230             
231            } else {
232             
233              if (numb[l]) {
234                /*If the level of w remains nonempty.*/
235               
236                level.set(w,++newlevel);
237                ++numb[newlevel];
238                stack[newlevel].push(w);
239                b=newlevel;
240              } else {
241                /*If the level of w gets empty.*/
242             
243                for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid() ; ++v) {
244                  if (level.get(v) >= l ) {
245                    level.set(v,n); 
246                  }
247                }
248               
249                for (int i=l+1 ; i!=n ; ++i) numb[i]=0;
250              } //if (numb[l])
251       
252            } // if (newlevel = n)
253         
254          } // if (excess.get(w)>0)
255
256
257        } //else
258       
259      } //while(b)
260
261      value=excess.get(t);
262      /*Max flow value.*/
263     
264
265
266      /*
267        We find an empty level, e. The Nodes above this level give
268        a minimum cut.
269      */
270     
271      int e=1;
272     
273      while(e) {
274        if(numb[e]) ++e;
275        else break;
276      }
277      for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v) {
278        if (level.get(v) > e) mincutvector.set(v, true);
279      }
280     
281
282    } // void run()
283
284
285
286    /*
287      Returns the maximum value of a flow.
288     */
289
290    T maxflow() {
291      return value;
292    }
293
294
295
296    /*
297      Returns a minimum cut.
298    */
299   
300    typename Graph::NodeMap<bool> mincut() {
301      return mincutvector;
302    }
303   
304
305  };
306}//namespace marci
307#endif
308
309
310
311
312
Note: See TracBrowser for help on using the repository browser.