65 */ |
64 */ |
66 void run() { |
65 void run() { |
67 |
66 |
68 typename Graph::NodeMap<int> level(G); |
67 typename Graph::NodeMap<int> level(G); |
69 typename Graph::NodeMap<T> excess(G); |
68 typename Graph::NodeMap<T> excess(G); |
70 std::cout <<"excess s:"<<excess.get(s); //delme |
69 |
71 |
|
72 int n=G.nodeNum(); |
70 int n=G.nodeNum(); |
73 int b=n-2; |
71 int b=n-2; |
74 /* |
72 /* |
75 b is a bound on the highest level of an active node. |
73 b is a bound on the highest level of an active node. |
76 In the beginning it is at most n-2. |
74 In the beginning it is at most n-2. |
77 */ |
75 */ |
78 |
76 |
|
77 std::vector<int> numb(n); //The number of nodes on level i < n. |
79 std::vector<std::stack<NodeIt> > stack(2*n-1); |
78 std::vector<std::stack<NodeIt> > stack(2*n-1); |
80 //Stack of the active nodes in level i. |
79 //Stack of the active nodes in level i. |
81 |
80 |
82 |
81 |
83 /*Reverse_bfs from t, to find the starting level.*/ |
82 /*Reverse_bfs from t, to find the starting level.*/ |
84 reverse_bfs<Graph> bfs(G, t); |
83 reverse_bfs<Graph> bfs(G, t); |
85 bfs.run(); |
84 bfs.run(); |
86 for(EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v) |
85 for(EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v) |
87 { |
86 { |
88 level.set(v, bfs.dist(v)); |
87 int dist=bfs.dist(v); |
|
88 level.set(v, dist); |
|
89 ++numb[dist]; |
89 } |
90 } |
90 |
91 |
91 level.set(s,n); |
92 level.set(s,n); |
92 |
93 |
93 |
94 |
94 /* Starting flow. It is everywhere 0 at the moment. */ |
95 /* Starting flow. It is everywhere 0 at the moment. */ |
95 for(OutEdgeIt e=G.template first<OutEdgeIt>(s); e.valid(); ++e) |
96 for(OutEdgeIt e=G.template first<OutEdgeIt>(s); e.valid(); ++e) |
96 { |
97 { |
97 if ( capacity.get(e) > 0 ) { |
98 if ( capacity.get(e) > 0 ) { |
98 NodeIt w=G.head(e); |
99 NodeIt w=G.head(e); |
99 if ( excess.get(w) == 0 && w!=t && w!=s ) stack[level.get(w)].push(w); |
100 if ( w!=s ) { |
100 flow.set(e, capacity.get(e)); |
101 if ( excess.get(w) == 0 && w!=t ) stack[level.get(w)].push(w); |
101 excess.set(w, excess.get(w)+capacity.get(e)); |
102 flow.set(e, capacity.get(e)); |
|
103 excess.set(w, excess.get(w)+capacity.get(e)); |
|
104 } |
102 } |
105 } |
103 } |
106 } |
104 |
107 |
105 /* |
108 /* |
106 End of preprocessing |
109 End of preprocessing |
130 if ( flow.get(e) < capacity.get(e) ) { |
133 if ( flow.get(e) < capacity.get(e) ) { |
131 /*e is an edge of the residual graph */ |
134 /*e is an edge of the residual graph */ |
132 |
135 |
133 NodeIt v=G.head(e); /*e is the edge wv.*/ |
136 NodeIt v=G.head(e); /*e is the edge wv.*/ |
134 |
137 |
135 if( level.get(w) > level.get(v)+1 ) { std::cout<<"HIBA1!";} //delme |
|
136 |
|
137 if( level.get(w) == level.get(v)+1 ) { |
138 if( level.get(w) == level.get(v)+1 ) { |
138 /*Push is allowed now*/ |
139 /*Push is allowed now*/ |
139 |
140 |
140 if (capacity.get(e)-flow.get(e) > excess.get(w)) { |
141 if ( excess.get(v)==0 && v != s && v !=t ) stack[level.get(v)].push(v); |
|
142 /*v becomes active.*/ |
|
143 |
|
144 if ( capacity.get(e)-flow.get(e) > excess.get(w) ) { |
141 /*A nonsaturating push.*/ |
145 /*A nonsaturating push.*/ |
142 |
146 |
143 if (excess.get(v)==0 && v != s && v !=t) stack[level.get(v)].push(v); |
|
144 /*v becomes active.*/ |
|
145 |
|
146 flow.set(e, flow.get(e)+excess.get(w)); |
147 flow.set(e, flow.get(e)+excess.get(w)); |
147 excess.set(v, excess.get(v)+excess.get(w)); |
148 excess.set(v, excess.get(v)+excess.get(w)); |
148 excess.set(w,0); |
149 excess.set(w,0); |
149 //std::cout << w << " " << v <<" elore elen nonsat pump " << std::endl; |
|
150 break; |
150 break; |
|
151 |
151 } else { |
152 } else { |
152 /*A saturating push.*/ |
153 /*A saturating push.*/ |
153 |
|
154 if (excess.get(v)==0 && v != s&& v !=t) stack[level.get(v)].push(v); |
|
155 /*v becomes active.*/ |
|
156 |
154 |
157 excess.set(v, excess.get(v)+capacity.get(e)-flow.get(e)); |
155 excess.set(v, excess.get(v)+capacity.get(e)-flow.get(e)); |
158 excess.set(w, excess.get(w)-capacity.get(e)+flow.get(e)); |
156 excess.set(w, excess.get(w)-capacity.get(e)+flow.get(e)); |
159 flow.set(e, capacity.get(e)); |
157 flow.set(e, capacity.get(e)); |
160 //std::cout << w<<" " <<v<<" elore elen sat pump " << std::endl; |
158 if ( excess.get(w)==0 ) break; |
161 if (excess.get(w)==0) break; |
159 /*If w is not active any more, then we go on to the next node.*/ |
162 /*If w is not active any more, then we go on to the next Node.*/ |
|
163 |
160 |
164 //std::cout<<++i; |
161 } |
|
162 } else { |
|
163 newlevel = newlevel < level.get(v) ? newlevel : level.get(v); |
|
164 } |
|
165 |
|
166 } //if the out edge wv is in the res graph |
|
167 |
|
168 } //for out edges wv |
|
169 |
|
170 |
|
171 if ( excess.get(w) > 0 ) { |
|
172 |
|
173 for( InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) { |
|
174 NodeIt v=G.tail(e); /*e is the edge vw.*/ |
|
175 |
|
176 if( flow.get(e) > 0 ) { |
|
177 /*e is an edge of the residual graph */ |
|
178 |
|
179 if( level.get(w)==level.get(v)+1 ) { |
|
180 /*Push is allowed now*/ |
|
181 |
|
182 if ( excess.get(v)==0 && v != s && v !=t) stack[level.get(v)].push(v); |
|
183 /*v becomes active.*/ |
|
184 |
|
185 if ( flow.get(e) > excess.get(w) ) { |
|
186 /*A nonsaturating push.*/ |
165 |
187 |
166 } // if (capacity.get(e)-flow.get(e) > excess.get(w)) |
188 flow.set(e, flow.get(e)-excess.get(w)); |
167 } // if(level.get(w)==level.get(v)+1) |
189 excess.set(v, excess.get(v)+excess.get(w)); |
168 |
190 excess.set(w,0); |
169 else {newlevel = newlevel < level.get(v) ? newlevel : level.get(v);} |
191 break; |
170 |
192 } else { |
171 } //if (flow.get(e)<capacity.get(e)) |
193 /*A saturating push.*/ |
172 |
194 |
173 } //for(OutEdgeIt e=G.first_OutEdge(w); e.valid(); ++e) |
195 excess.set(v, excess.get(v)+flow.get(e)); |
|
196 excess.set(w, excess.get(w)-flow.get(e)); |
|
197 flow.set(e,0); |
|
198 if ( excess.get(w)==0 ) break; |
|
199 } |
|
200 } else { |
|
201 newlevel = newlevel < level.get(v) ? newlevel : level.get(v); |
|
202 } |
|
203 |
|
204 } //if in edge vw is in the res graph |
|
205 |
|
206 } //for in edges vw |
|
207 |
|
208 } // if w still has excess after the out edge for cycle |
|
209 |
|
210 |
|
211 /* |
|
212 Relabel |
|
213 */ |
174 |
214 |
175 |
215 if ( excess.get(w) > 0 ) { |
176 |
216 |
177 for(InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) { |
217 int oldlevel=level.get(w); |
178 NodeIt v=G.tail(e); |
218 level.set(w,++newlevel); |
179 /*e is the Edge vw.*/ |
219 |
180 |
220 if ( oldlevel < n ) { |
181 if (excess.get(w)==0) break; |
221 --numb[oldlevel]; |
182 /*It may happen, that w became inactive in the first for cycle.*/ |
222 |
183 if(flow.get(e)>0) { |
223 if ( !numb[oldlevel] ) { //If the level of w gets empty. |
184 /*e is an Edge of the residual graph */ |
|
185 |
|
186 if(level.get(w)==level.get(v)+1) { |
|
187 /*Push is allowed now*/ |
|
188 |
224 |
189 if (flow.get(e) > excess.get(w)) { |
225 for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid() ; ++v) { |
190 /*A nonsaturating push.*/ |
226 if (level.get(v) > oldlevel && level.get(v) < n ) level.set(v,n); |
191 |
227 } |
192 if (excess.get(v)==0 && v != s&& v !=t) stack[level.get(v)].push(v); |
228 for (int i=oldlevel+1 ; i!=n ; ++i) numb[i]=0; |
193 /*v becomes active.*/ |
229 if ( newlevel < n ) newlevel=n; |
194 |
230 } else { |
195 flow.set(e, flow.get(e)-excess.get(w)); |
231 if ( newlevel < n ) ++numb[newlevel]; |
196 excess.set(v, excess.get(v)+excess.get(w)); |
232 } |
197 excess.set(w,0); |
233 } else { |
198 //std::cout << v << " " << w << " vissza elen nonsat pump " << std::endl; |
234 if ( newlevel < n ) ++numb[newlevel]; |
199 break; |
235 } |
200 } else { |
236 |
201 /*A saturating push.*/ |
|
202 |
|
203 if (excess.get(v)==0 && v != s&& v !=t) stack[level.get(v)].push(v); |
|
204 /*v becomes active.*/ |
|
205 |
|
206 excess.set(v, excess.get(v)+flow.get(e)); |
|
207 excess.set(w, excess.get(w)-flow.get(e)); |
|
208 flow.set(e,0); |
|
209 //std::cout << v <<" " << w << " vissza elen sat pump " << std::endl; |
|
210 if (excess.get(w)==0) { break;} |
|
211 } //if (flow.get(e) > excess.get(v)) |
|
212 } //if(level.get(w)==level.get(v)+1) |
|
213 |
|
214 else {newlevel = newlevel < level.get(v) ? newlevel : level.get(v);} |
|
215 |
|
216 |
|
217 } //if (flow.get(e)>0) |
|
218 |
|
219 } //for |
|
220 |
|
221 |
|
222 if (excess.get(w)>0) { |
|
223 level.set(w,++newlevel); |
|
224 stack[newlevel].push(w); |
237 stack[newlevel].push(w); |
225 b=newlevel; |
238 b=newlevel; |
226 //std::cout << "The new level of " << w << " is "<< newlevel <<std::endl; |
239 |
227 } |
240 } |
228 |
241 |
229 |
242 } // if stack[b] is nonempty |
230 } //else |
243 |
231 |
244 } // while(b) |
232 } //while(b) |
245 |
233 |
246 |
234 value = excess.get(t); |
247 value = excess.get(t); |
235 /*Max flow value.*/ |
248 /*Max flow value.*/ |
236 |
|
237 |
|
238 |
249 |
239 |
250 |
240 } //void run() |
251 } //void run() |
241 |
252 |
242 |
253 |