Changeset 97:a5127ecb2914 in lemon0.x for src/work/jacint/preflow_push_max_flow.h
 Timestamp:
 02/18/04 15:42:38 (16 years ago)
 Branch:
 default
 Phase:
 public
 Convert:
 svn:c9d7d8f590d60310b91f818b3a526b0e/lemon/trunk@126
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

src/work/jacint/preflow_push_max_flow.h
r85 r97 1 // * C++ * 1 2 /* 2 3 preflow_push_max_flow_h … … 16 17 T maxflow() : returns the value of a maximum flow 17 18 18 NodeMap<bool> mincut(): returns a 19 characteristic vector of a minimum cut. 19 void mincut(CutMap& M) : sets M to the characteristic vector of a 20 minimum cut. M should be a map of bools initialized to false. 21 20 22 */ 21 23 22 24 #ifndef PREFLOW_PUSH_MAX_FLOW_H 23 25 #define PREFLOW_PUSH_MAX_FLOW_H 26 27 #define A 1 24 28 25 29 #include <algorithm> … … 32 36 namespace marci { 33 37 34 template <typename Graph, typename T> 38 template <typename Graph, typename T, 39 typename FlowMap=typename Graph::EdgeMap<T>, typename CapMap=typename Graph::EdgeMap<T>, 40 typename IntMap=typename Graph::NodeMap<int>, typename TMap=typename Graph::NodeMap<T> > 35 41 class preflow_push_max_flow { 36 42 … … 43 49 NodeIt s; 44 50 NodeIt t; 45 typename Graph::EdgeMap<T>& capacity; 46 T value; 47 typename Graph::NodeMap<bool> mincutvector; 48 49 50 51 IntMap level; 52 CapMap& capacity; 53 int empty_level; //an empty level in the end of run() 54 T value; 55 51 56 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) { } 57 58 preflow_push_max_flow(Graph& _G, NodeIt _s, NodeIt _t, CapMap& _capacity) : 59 G(_G), s(_s), t(_t), level(_G), capacity(_capacity) { } 56 60 57 61 … … 63 67 void run() { 64 68 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(); 69 int n=G.nodeNum(); 70 70 int b=n2; 71 71 /* 72 b is a bound on the highest level of an active Node. 73 In the beginning it is at most n2. 74 */ 75 76 std::vector<int> numb(n); //The number of Nodes on level i < n. 77 std::vector<std::stack<NodeIt> > stack(2*n1); 78 //Stack of the active Nodes in level i. 72 b is a bound on the highest level of an active node. 73 */ 74 75 IntMap level(G,n); 76 TMap excess(G); 77 FlowMap flow(G,0); 78 79 std::vector<int> numb(n); 80 /* 81 The number of nodes on level i < n. It is 82 initialized to n+1, because of the reverse_bfspart. 83 */ 84 85 std::vector<std::stack<NodeIt> > stack(n); 86 //Stack of the active nodes in level i. 87 79 88 80 89 /*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]; 90 level.set(t,0); 91 std::queue<NodeIt> bfs_queue; 92 bfs_queue.push(t); 93 94 while (!bfs_queue.empty()) { 95 96 NodeIt v=bfs_queue.front(); 97 bfs_queue.pop(); 98 int l=level.get(v)+1; 99 100 for(InEdgeIt e=G.template first<InEdgeIt>(v); e.valid(); ++e) { 101 NodeIt w=G.tail(e); 102 if ( level.get(w) == n ) { 103 bfs_queue.push(w); 104 ++numb[l]; 105 level.set(w, l); 106 } 88 107 } 89 108 } 109 90 110 level.set(s,n); 91 111 92 /* Starting flow. It is everywhere 0 at the moment. */ 112 113 114 /* Starting flow. It is everywhere 0 at the moment. */ 93 115 for(OutEdgeIt e=G.template first<OutEdgeIt>(s); e.valid(); ++e) 94 116 { 95 if ( capacity.get(e) > 0 ) { 96 NodeIt w=G.head(e); 117 if ( capacity.get(e) == 0 ) continue; 118 NodeIt w=G.head(e); 119 if ( level.get(w) < n ) { 120 if ( excess.get(w) == 0 && w!=t ) stack[level.get(w)].push(w); 97 121 flow.set(e, capacity.get(e)); 98 stack[level.get(w)].push(w);99 122 excess.set(w, excess.get(w)+capacity.get(e)); 100 123 } 101 124 } 102 125 103 126 /* 104 127 End of preprocessing … … 106 129 107 130 108 109 131 /* 110 Push/relabel on the highest level active Nodes. 111 */ 112 113 /*While there exists an active Node.*/ 132 Push/relabel on the highest level active nodes. 133 */ 134 /*While there exists an active node.*/ 114 135 while (b) { 115 116 /*We decrease the bound if there is no active node of level b.*/ 117 if (stack[b].empty()) { 136 if ( stack[b].empty() ) { 118 137 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*n2; //In newlevel we maintain the next level of w. 125 138 continue; 139 } 140 141 NodeIt w=stack[b].top(); //w is a highest label active node. 142 stack[b].pop(); 143 int lev=level.get(w); 144 int exc=excess.get(w); 145 int newlevel=2*n2; //In newlevel we bound the next level of w. 146 147 // if ( level.get(w) < n ) { //Nem tudom ez mukodike 126 148 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) 149 150 if ( flow.get(e) == capacity.get(e) ) continue; 151 NodeIt v=G.head(e); 152 //e=wv 153 154 if( lev > level.get(v) ) { 155 /*Push is allowed now*/ 156 157 if ( excess.get(v)==0 && v != s && v !=t ) 158 stack[level.get(v)].push(v); 159 /*v becomes active.*/ 160 161 int cap=capacity.get(e); 162 int flo=flow.get(e); 163 int remcap=capflo; 164 165 if ( remcap >= exc ) { 166 /*A nonsaturating push.*/ 167 168 flow.set(e, flo+exc); 169 excess.set(v, excess.get(v)+exc); 170 exc=0; 171 break; 172 173 } else { 174 /*A saturating push.*/ 175 176 flow.set(e, cap ); 177 excess.set(v, excess.get(v)+remcap); 178 exc=remcap; 179 } 180 } else if ( newlevel > level.get(v) ) newlevel = level.get(v); 181 182 } //for out edges wv 183 184 185 if ( exc > 0 ) { 186 for( InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) { 187 188 if( flow.get(e) == 0 ) continue; 189 NodeIt v=G.tail(e); 190 //e=vw 191 192 if( lev > level.get(v) ) { 193 /*Push is allowed now*/ 194 195 if ( excess.get(v)==0 && v != s && v !=t) 196 stack[level.get(v)].push(v); 197 /*v becomes active.*/ 198 199 int flo=flow.get(e); 200 201 if ( flo >= exc ) { 202 /*A nonsaturating push.*/ 203 204 flow.set(e, floexc); 205 excess.set(v, excess.get(v)+exc); 206 exc=0; 207 break; 208 } else { 209 /*A saturating push.*/ 210 211 excess.set(v, excess.get(v)+flo); 212 exc=flo; 213 flow.set(e,0); 214 } 215 } else if ( newlevel > level.get(v) ) newlevel = level.get(v); 216 217 } //for in edges vw 168 218 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 inEdge 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 219 } // if w still has excess after the out edge for cycle 220 221 excess.set(w, exc); 222 223 224 /* 225 Relabel 226 */ 227 228 if ( exc > 0 ) { 229 //now 'lev' is the old level of w 230 level.set(w,++newlevel); 231 numb[lev]; 232 233 if ( !numb[lev] && lev < A*n ) { //If the level of w gets empty. 234 235 for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid() ; ++v) { 236 if (level.get(v) > lev ) level.set(v,n); 237 } 238 for (int i=lev+1 ; i!=n ; ++i) numb[i]=0; 239 if ( newlevel < n ) newlevel=n; 240 } else if ( newlevel < n ) { 241 ++numb[newlevel]; 242 stack[newlevel].push(w); 243 b=newlevel; 244 } 245 } 246 247 248 259 249 } //while(b) 260 250 … … 263 253 264 254 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; 255 /* 256 We count empty_level. The nodes above this level is a mincut. 257 */ 258 while(true) { 259 if(numb[empty_level]) ++empty_level; 275 260 else break; 276 261 } 277 for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v) { 278 if (level.get(v) > e) mincutvector.set(v, true); 279 } 280 281 262 282 263 } // void run() 283 264 … … 296 277 /* 297 278 Returns a minimum cut. 298 */ 299 300 typename Graph::NodeMap<bool> mincut() { 301 return mincutvector; 279 */ 280 template<typename CutMap> 281 void mincut(CutMap& M) { 282 283 for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v) { 284 if ( level.get(v) > empty_level ) M.set(v, true); 285 } 302 286 } 303 287 304 288 305 289 };
Note: See TracChangeset
for help on using the changeset viewer.