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