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