COIN-OR::LEMON - Graph Library

source: lemon-0.x/src/work/jacint/preflow_push_hl.h @ 101:d2ac583ed195

Last change on this file since 101:d2ac583ed195 was 101:d2ac583ed195, checked in by jacint, 17 years ago

another heuristic

File size: 8.1 KB
Line 
1// -*- C++ -*-
2
3//kerdesek: nem tudom lehet-e a
4//kieleket csak a legf n szintu pontokra nezni.
5
6/*
7preflow_push_hl.h
8by jacint.
9Runs the highest label variant of the preflow push algorithm with
10running time O(n^2\sqrt(m)), and with the 'empty level' heuristic.
11
12'A' is a parameter for the empty_level heuristic
13
14Member functions:
15
16void run() : runs the algorithm
17
18 The following functions should be used after run() was already run.
19
20T maxflow() : returns the value of a maximum flow
21
22T flowonedge(EdgeIt e) : for a fixed maximum flow x it returns x(e)
23
24FlowMap allflow() : returns the fixed maximum flow x
25
26void mincut(CutMap& M) : sets M to the characteristic vector of a
27     minimum cut. M should be a map of bools initialized to false.
28
29void min_mincut(CutMap& M) : sets M to the characteristic vector of the
30     minimum min cut. M should be a map of bools initialized to false.
31
32void max_mincut(CutMap& M) : sets M to the characteristic vector of the
33     maximum min cut. M should be a map of bools initialized to false.
34
35*/
36
37#ifndef PREFLOW_PUSH_HL_H
38#define PREFLOW_PUSH_HL_H
39
40#define A 1
41
42#include <vector>
43#include <stack>
44#include <queue>
45
46namespace marci {
47
48  template <typename Graph, typename T,
49    typename FlowMap=typename Graph::EdgeMap<T>, typename CapMap=typename Graph::EdgeMap<T>,
50    typename IntMap=typename Graph::NodeMap<int>, typename TMap=typename Graph::NodeMap<T> >
51  class preflow_push_hl {
52   
53    typedef typename Graph::NodeIt NodeIt;
54    typedef typename Graph::EdgeIt EdgeIt;
55    typedef typename Graph::EachNodeIt EachNodeIt;
56    typedef typename Graph::OutEdgeIt OutEdgeIt;
57    typedef typename Graph::InEdgeIt InEdgeIt;
58   
59    Graph& G;
60    NodeIt s;
61    NodeIt t;
62    FlowMap flow;
63    CapMap& capacity; 
64    T value;
65   
66  public:
67
68    preflow_push_hl(Graph& _G, NodeIt _s, NodeIt _t, CapMap& _capacity) :
69      G(_G), s(_s), t(_t), flow(_G, 0), capacity(_capacity) { }
70
71
72   
73
74    void run() {
75 
76      int n=G.nodeNum();
77      int b=n-2;
78      /*
79        b is a bound on the highest level of an active node.
80      */
81
82      IntMap level(G,n);     
83      TMap excess(G);
84
85      std::vector<int> numb(n);   
86      /*
87        The number of nodes on level i < n. It is
88        initialized to n+1, because of the reverse_bfs-part.
89      */
90
91      std::vector<std::stack<NodeIt> > stack(2*n-1);   
92      //Stack of the active nodes in level i.
93
94
95      /*Reverse_bfs from t, to find the starting level.*/
96      level.set(t,0);
97      std::queue<NodeIt> bfs_queue;
98      bfs_queue.push(t);
99
100      while (!bfs_queue.empty()) {
101
102        NodeIt v=bfs_queue.front();     
103        bfs_queue.pop();
104        int l=level.get(v)+1;
105
106        for(InEdgeIt e=G.template first<InEdgeIt>(v); e.valid(); ++e) {
107          NodeIt w=G.tail(e);
108          if ( level.get(w) == n ) {
109            bfs_queue.push(w);
110            ++numb[l];
111            level.set(w, l);
112          }
113        }
114      }
115       
116      level.set(s,n);
117
118
119
120      /* Starting flow. It is everywhere 0 at the moment. */     
121      for(OutEdgeIt e=G.template first<OutEdgeIt>(s); e.valid(); ++e)
122        {
123          T c=capacity.get(e);
124          if ( c == 0 ) continue;
125          NodeIt w=G.head(e);
126          if ( w!=s ) {   
127            if ( excess.get(w) == 0 && w!=t ) stack[level.get(w)].push(w);
128            flow.set(e, c);
129            excess.set(w, excess.get(w)+c);
130          }
131        }
132     
133      /*
134         End of preprocessing
135      */
136
137
138      /*
139        Push/relabel on the highest level active nodes.
140      */
141      /*While there exists an active node.*/
142      while (b) {
143        if ( stack[b].empty() ) {
144          --b;
145          continue;
146        }
147       
148        NodeIt w=stack[b].top();        //w is a highest label active node.
149        stack[b].pop();           
150        int lev=level.get(w);
151        int exc=excess.get(w);
152        int newlevel=2*n;      //In newlevel we bound the next level of w.
153        //vagy MAXINT   
154
155        //  if ( level.get(w) < n ) { //Nem tudom ez mukodik-e
156          for(OutEdgeIt e=G.template first<OutEdgeIt>(w); e.valid(); ++e) {
157           
158            if ( flow.get(e) == capacity.get(e) ) continue;
159            NodeIt v=G.head(e);           
160            //e=wv         
161           
162            if( lev > level.get(v) ) {     
163              /*Push is allowed now*/
164             
165              if ( excess.get(v)==0 && v != s && v !=t )
166                stack[level.get(v)].push(v);
167              /*v becomes active.*/
168             
169              int cap=capacity.get(e);
170              int flo=flow.get(e);
171              int remcap=cap-flo;
172             
173              if ( remcap >= exc ) {       
174                /*A nonsaturating push.*/
175               
176                flow.set(e, flo+exc);
177                excess.set(v, excess.get(v)+exc);
178                exc=0;
179                break;
180               
181              } else {
182                /*A saturating push.*/
183               
184                flow.set(e, cap );
185                excess.set(v, excess.get(v)+remcap);
186                exc-=remcap;
187              }
188            } else if ( newlevel > level.get(v) ) newlevel = level.get(v);
189           
190          } //for out edges wv
191       
192       
193        if ( exc > 0 ) {       
194          for( InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) {
195           
196            if( flow.get(e) == 0 ) continue;
197            NodeIt v=G.tail(e); 
198            //e=vw
199           
200            if( lev > level.get(v) ) { 
201              /*Push is allowed now*/
202             
203              if ( excess.get(v)==0 && v != s && v !=t)
204                stack[level.get(v)].push(v);
205              /*v becomes active.*/
206             
207              int flo=flow.get(e);
208             
209              if ( flo >= exc ) {
210                /*A nonsaturating push.*/
211               
212                flow.set(e, flo-exc);
213                excess.set(v, excess.get(v)+exc);
214                exc=0;
215                break;
216              } else {                                               
217                /*A saturating push.*/
218               
219                excess.set(v, excess.get(v)+flo);
220                exc-=flo;
221                flow.set(e,0);
222              } 
223            } else if ( newlevel > level.get(v) ) newlevel = level.get(v);
224           
225          } //for in edges vw
226         
227        } // if w still has excess after the out edge for cycle
228       
229        excess.set(w, exc);
230       
231
232       
233
234        /*
235          Relabel
236        */
237       
238        if ( exc > 0 ) {
239          //now 'lev' is the old level of w
240          level.set(w,++newlevel);
241         
242          if ( lev < n ) {
243            --numb[lev];
244           
245            if ( !numb[lev] && lev < A*n ) {  //If the level of w gets empty.
246             
247              for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid() ; ++v) {
248                if (level.get(v) > lev && level.get(v) < n ) level.set(v,n); 
249              }
250              for (int i=lev+1 ; i!=n ; ++i) numb[i]=0;
251              if ( newlevel < n ) newlevel=n;
252            } else {
253              if ( newlevel < n ) ++numb[newlevel];
254            }
255          }
256         
257          stack[newlevel].push(w);
258          b=newlevel;
259         
260        }
261       
262      } // while(b)
263     
264     
265      value = excess.get(t);
266      /*Max flow value.*/
267
268
269    } //void run()
270
271
272
273
274
275    /*
276      Returns the maximum value of a flow.
277     */
278
279    T maxflow() {
280      return value;
281    }
282
283
284
285    /*
286      For the maximum flow x found by the algorithm, it returns the flow value on Edge e, i.e. x(e).
287    */
288
289    T flowonedge(const EdgeIt e) {
290      return flow.get(e);
291    }
292
293
294
295    /*
296      Returns the maximum flow x found by the algorithm.
297    */
298
299    FlowMap allflow() {
300      return flow;
301    }
302
303
304
305
306    /*
307      Returns the minimum min cut, by a bfs from s in the residual graph.
308    */
309   
310    template<typename CutMap>
311    void mincut(CutMap& M) {
312   
313      std::queue<NodeIt> queue;
314     
315      M.set(s,true);     
316      queue.push(s);
317
318      while (!queue.empty()) {
319        NodeIt w=queue.front();
320        queue.pop();
321       
322        for(OutEdgeIt e=G.template first<OutEdgeIt>(w) ; e.valid(); ++e) {
323          NodeIt v=G.head(e);
324          if (!M.get(v) && flow.get(e) < capacity.get(e) ) {
325            queue.push(v);
326            M.set(v, true);
327          }
328        }
329
330        for(InEdgeIt e=G.template first<InEdgeIt>(w) ; e.valid(); ++e) {
331          NodeIt v=G.tail(e);
332          if (!M.get(v) && flow.get(e) > 0 ) {
333            queue.push(v);
334            M.set(v, true);
335          }
336        }
337
338      }
339    }
340
341
342
343    /*
344      Returns the maximum min cut, by a reverse bfs
345      from t in the residual graph.
346    */
347   
348    template<typename CutMap>
349    void max_mincut(CutMap& M) {
350   
351      std::queue<NodeIt> queue;
352     
353      M.set(t,true);       
354      queue.push(t);
355
356      while (!queue.empty()) {
357        NodeIt w=queue.front();
358        queue.pop();
359
360        for(InEdgeIt e=G.template first<InEdgeIt>(w) ; e.valid(); ++e) {
361          NodeIt v=G.tail(e);
362          if (!M.get(v) && flow.get(e) < capacity.get(e) ) {
363            queue.push(v);
364            M.set(v, true);
365          }
366        }
367
368        for(OutEdgeIt e=G.template first<OutEdgeIt>(w) ; e.valid(); ++e) {
369          NodeIt v=G.head(e);
370          if (!M.get(v) && flow.get(e) > 0 ) {
371            queue.push(v);
372            M.set(v, true);
373          }
374        }
375      }
376
377      for(EachNodeIt v=G.template first<EachNodeIt>() ; v.valid(); ++v) {
378        M.set(v, !M.get(v));
379      }
380
381    }
382
383
384
385    template<typename CutMap>
386    void min_mincut(CutMap& M) {
387      mincut(M);
388    }
389
390
391
392  };
393}//namespace marci
394#endif
395
396
397
398
Note: See TracBrowser for help on using the repository browser.