| 1 | #ifndef HUGO_PREFLOW_PUSH_HH | 
|---|
| 2 | #define HUGO_PREFLOW_PUSH_HH | 
|---|
| 3 |  | 
|---|
| 4 | //#include <algorithm> | 
|---|
| 5 | #include <list> | 
|---|
| 6 | #include <vector> | 
|---|
| 7 | #include <queue> | 
|---|
| 8 | //#include "pf_hiba.hh" | 
|---|
| 9 | //#include <marci_list_graph.hh> | 
|---|
| 10 | //#include <marci_graph_traits.hh> | 
|---|
| 11 | #include <invalid.h> | 
|---|
| 12 | //#include <reverse_bfs.hh> | 
|---|
| 13 |  | 
|---|
| 14 | using namespace std; | 
|---|
| 15 |  | 
|---|
| 16 | namespace hugo { | 
|---|
| 17 |  | 
|---|
| 18 | template <typename Graph, typename T> | 
|---|
| 19 | class preflow_push { | 
|---|
| 20 |  | 
|---|
| 21 | //Useful typedefs | 
|---|
| 22 | typedef typename Graph::Node Node; | 
|---|
| 23 | typedef typename Graph::NodeIt NodeIt; | 
|---|
| 24 | typedef typename Graph::Edge Edge; | 
|---|
| 25 | typedef typename Graph::OutEdgeIt OutEdgeIt; | 
|---|
| 26 | typedef typename Graph::InEdgeIt InEdgeIt; | 
|---|
| 27 |  | 
|---|
| 28 |  | 
|---|
| 29 | //--------------------------------------------- | 
|---|
| 30 | //Parameters of the algorithm | 
|---|
| 31 | //--------------------------------------------- | 
|---|
| 32 | //Fully examine an active node until excess becomes 0 | 
|---|
| 33 | enum node_examination_t {examine_full, examine_to_relabel}; | 
|---|
| 34 | //No more implemented yet:, examine_only_one_edge}; | 
|---|
| 35 | node_examination_t node_examination; | 
|---|
| 36 | //Which implementation to be used | 
|---|
| 37 | enum implementation_t {impl_fifo, impl_highest_label}; | 
|---|
| 38 | //No more implemented yet:}; | 
|---|
| 39 | implementation_t implementation; | 
|---|
| 40 | //--------------------------------------------- | 
|---|
| 41 | //Parameters of the algorithm | 
|---|
| 42 | //--------------------------------------------- | 
|---|
| 43 |  | 
|---|
| 44 | private: | 
|---|
| 45 | //input | 
|---|
| 46 | Graph& G; | 
|---|
| 47 | Node s; | 
|---|
| 48 | Node t; | 
|---|
| 49 | typename Graph::EdgeMap<T> &capacity; | 
|---|
| 50 |  | 
|---|
| 51 | //output | 
|---|
| 52 | typename Graph::EdgeMap<T> preflow; | 
|---|
| 53 | T maxflow_value; | 
|---|
| 54 |  | 
|---|
| 55 | //auxiliary variables for computation | 
|---|
| 56 | //The number of the nodes | 
|---|
| 57 | int number_of_nodes; | 
|---|
| 58 | //A nodemap for the level | 
|---|
| 59 | typename Graph::NodeMap<int> level; | 
|---|
| 60 | //A nodemap for the excess | 
|---|
| 61 | typename Graph::NodeMap<T> excess; | 
|---|
| 62 |  | 
|---|
| 63 | //Number of nodes on each level | 
|---|
| 64 | vector<int> num_of_nodes_on_level; | 
|---|
| 65 |  | 
|---|
| 66 | //For the FIFO implementation | 
|---|
| 67 | list<Node> fifo_nodes; | 
|---|
| 68 | //For 'highest label' implementation | 
|---|
| 69 | int highest_active; | 
|---|
| 70 | //int second_highest_active; | 
|---|
| 71 | vector< list<Node> > active_nodes; | 
|---|
| 72 |  | 
|---|
| 73 | public: | 
|---|
| 74 |  | 
|---|
| 75 | //Constructing the object using the graph, source, sink and capacity vector | 
|---|
| 76 | preflow_push( | 
|---|
| 77 | Graph& _G, | 
|---|
| 78 | Node _s, | 
|---|
| 79 | Node _t, | 
|---|
| 80 | typename Graph::EdgeMap<T> & _capacity) | 
|---|
| 81 | : G(_G), s(_s), t(_t), | 
|---|
| 82 | capacity(_capacity), | 
|---|
| 83 | preflow(_G), | 
|---|
| 84 | //Counting the number of nodes | 
|---|
| 85 | //number_of_nodes(count(G.first<EachNodeIt>())), | 
|---|
| 86 | number_of_nodes(G.nodeNum()), | 
|---|
| 87 |  | 
|---|
| 88 | level(_G), | 
|---|
| 89 | excess(_G)//, | 
|---|
| 90 | // Default constructor: active_nodes() | 
|---|
| 91 | { | 
|---|
| 92 | //Simplest parameter settings | 
|---|
| 93 | node_examination = examine_full;//examine_to_relabel;// | 
|---|
| 94 | //Which implementation to be usedexamine_full | 
|---|
| 95 | implementation = impl_highest_label;//impl_fifo; | 
|---|
| 96 |  | 
|---|
| 97 | // | 
|---|
| 98 | num_of_nodes_on_level.resize(2*number_of_nodes-1); | 
|---|
| 99 | num_of_nodes_on_level.clear(); | 
|---|
| 100 |  | 
|---|
| 101 | switch(implementation){ | 
|---|
| 102 | case impl_highest_label :{ | 
|---|
| 103 | active_nodes.clear(); | 
|---|
| 104 | active_nodes.resize(2*number_of_nodes-1); | 
|---|
| 105 |  | 
|---|
| 106 | break; | 
|---|
| 107 | } | 
|---|
| 108 | default: | 
|---|
| 109 | break; | 
|---|
| 110 | } | 
|---|
| 111 |  | 
|---|
| 112 | } | 
|---|
| 113 |  | 
|---|
| 114 | //Returns the value of a maximal flow | 
|---|
| 115 | T run(); | 
|---|
| 116 |  | 
|---|
| 117 | typename Graph::EdgeMap<T>  getmaxflow(){ | 
|---|
| 118 | return preflow; | 
|---|
| 119 | } | 
|---|
| 120 |  | 
|---|
| 121 |  | 
|---|
| 122 | private: | 
|---|
| 123 | //For testing purposes only | 
|---|
| 124 | //Lists the node_properties | 
|---|
| 125 | void write_property_vector(typename Graph::NodeMap<T> a, | 
|---|
| 126 | //node_property_vector<Graph, T> a, | 
|---|
| 127 | char* prop_name="property"){ | 
|---|
| 128 | for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) { | 
|---|
| 129 | cout<<"Node id.: "<<G.id(i)<<", "<<prop_name<<" value: "<<a[i]<<endl; | 
|---|
| 130 | } | 
|---|
| 131 | cout<<endl; | 
|---|
| 132 | } | 
|---|
| 133 |  | 
|---|
| 134 | //Modifies the excess of the node and makes sufficient changes | 
|---|
| 135 | void modify_excess(const Node& a ,T v){ | 
|---|
| 136 | //T old_value=excess[a]; | 
|---|
| 137 | excess[a] += v; | 
|---|
| 138 | } | 
|---|
| 139 |  | 
|---|
| 140 | //This private procedure is supposed to modify the preflow on edge j | 
|---|
| 141 | //by value v (which can be positive or negative as well) | 
|---|
| 142 | //and maintain the excess on the head and tail | 
|---|
| 143 | //Here we do not check whether this is possible or not | 
|---|
| 144 | void modify_preflow(Edge j, const T& v){ | 
|---|
| 145 |  | 
|---|
| 146 | //Modifiyng the edge | 
|---|
| 147 | preflow[j] += v; | 
|---|
| 148 |  | 
|---|
| 149 |  | 
|---|
| 150 | //Modifiyng the head | 
|---|
| 151 | modify_excess(G.head(j),v); | 
|---|
| 152 |  | 
|---|
| 153 | //Modifiyng the tail | 
|---|
| 154 | modify_excess(G.tail(j),-v); | 
|---|
| 155 |  | 
|---|
| 156 | } | 
|---|
| 157 |  | 
|---|
| 158 | //Gives the active node to work with | 
|---|
| 159 | //(depending on the implementation to be used) | 
|---|
| 160 | Node get_active_node(){ | 
|---|
| 161 |  | 
|---|
| 162 |  | 
|---|
| 163 | switch(implementation) { | 
|---|
| 164 | case impl_highest_label : { | 
|---|
| 165 |  | 
|---|
| 166 | //First need to find the highest label for which there's an active node | 
|---|
| 167 | while( highest_active>=0 && active_nodes[highest_active].empty() ){ | 
|---|
| 168 | --highest_active; | 
|---|
| 169 | } | 
|---|
| 170 |  | 
|---|
| 171 | if( highest_active>=0) { | 
|---|
| 172 |  | 
|---|
| 173 |  | 
|---|
| 174 | Node a=active_nodes[highest_active].front(); | 
|---|
| 175 | active_nodes[highest_active].pop_front(); | 
|---|
| 176 |  | 
|---|
| 177 | return a; | 
|---|
| 178 | } | 
|---|
| 179 | else { | 
|---|
| 180 | return INVALID; | 
|---|
| 181 | } | 
|---|
| 182 |  | 
|---|
| 183 | break; | 
|---|
| 184 |  | 
|---|
| 185 | } | 
|---|
| 186 | case impl_fifo : { | 
|---|
| 187 |  | 
|---|
| 188 | if( ! fifo_nodes.empty() ) { | 
|---|
| 189 | Node a=fifo_nodes.front(); | 
|---|
| 190 | fifo_nodes.pop_front(); | 
|---|
| 191 | return a; | 
|---|
| 192 | } | 
|---|
| 193 | else { | 
|---|
| 194 | return INVALID; | 
|---|
| 195 | } | 
|---|
| 196 | break; | 
|---|
| 197 | } | 
|---|
| 198 | } | 
|---|
| 199 | // | 
|---|
| 200 | return INVALID; | 
|---|
| 201 | } | 
|---|
| 202 |  | 
|---|
| 203 | //Puts node 'a' among the active nodes | 
|---|
| 204 | void make_active(const Node& a){ | 
|---|
| 205 | //s and t never become active | 
|---|
| 206 | if (a!=s && a!= t){ | 
|---|
| 207 | switch(implementation){ | 
|---|
| 208 | case impl_highest_label : | 
|---|
| 209 | active_nodes[level[a]].push_back(a); | 
|---|
| 210 | break; | 
|---|
| 211 | case impl_fifo : | 
|---|
| 212 | fifo_nodes.push_back(a); | 
|---|
| 213 | break; | 
|---|
| 214 | } | 
|---|
| 215 |  | 
|---|
| 216 | } | 
|---|
| 217 |  | 
|---|
| 218 | //Update highest_active label | 
|---|
| 219 | if (highest_active<level[a]){ | 
|---|
| 220 | highest_active=level[a]; | 
|---|
| 221 | } | 
|---|
| 222 |  | 
|---|
| 223 | } | 
|---|
| 224 |  | 
|---|
| 225 | //Changes the level of node a and make sufficent changes | 
|---|
| 226 | void change_level_to(Node a, int new_value){ | 
|---|
| 227 | int seged = level[a]; | 
|---|
| 228 | level.set(a,new_value); | 
|---|
| 229 | --num_of_nodes_on_level[seged]; | 
|---|
| 230 | ++num_of_nodes_on_level[new_value]; | 
|---|
| 231 | } | 
|---|
| 232 |  | 
|---|
| 233 | //Collection of things useful (or necessary) to do before running | 
|---|
| 234 |  | 
|---|
| 235 | void preprocess(){ | 
|---|
| 236 |  | 
|---|
| 237 | //--------------------------------------- | 
|---|
| 238 | //Initialize parameters | 
|---|
| 239 | //--------------------------------------- | 
|---|
| 240 |  | 
|---|
| 241 | //Setting starting preflow, level and excess values to zero | 
|---|
| 242 | //This can be important, if the algorithm is run more then once | 
|---|
| 243 | for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) { | 
|---|
| 244 | level.set(i,0); | 
|---|
| 245 | excess.set(i,0); | 
|---|
| 246 | for(OutEdgeIt j=G.template first<OutEdgeIt>(i); G.valid(j); G.next(j)) | 
|---|
| 247 | preflow.set(j, 0); | 
|---|
| 248 | } | 
|---|
| 249 | num_of_nodes_on_level[0]=number_of_nodes; | 
|---|
| 250 | highest_active=0; | 
|---|
| 251 | //--------------------------------------- | 
|---|
| 252 | //Initialize parameters | 
|---|
| 253 | //--------------------------------------- | 
|---|
| 254 |  | 
|---|
| 255 |  | 
|---|
| 256 | //------------------------------------ | 
|---|
| 257 | //This is the only part that uses BFS | 
|---|
| 258 | //------------------------------------ | 
|---|
| 259 |  | 
|---|
| 260 | /*Reverse_bfs from t, to find the starting level.*/ | 
|---|
| 261 | //Copyright: Jacint | 
|---|
| 262 | change_level_to(t,0); | 
|---|
| 263 |  | 
|---|
| 264 | std::queue<Node> bfs_queue; | 
|---|
| 265 | bfs_queue.push(t); | 
|---|
| 266 |  | 
|---|
| 267 | while (!bfs_queue.empty()) { | 
|---|
| 268 |  | 
|---|
| 269 | Node v=bfs_queue.front(); | 
|---|
| 270 | bfs_queue.pop(); | 
|---|
| 271 | int l=level[v]+1; | 
|---|
| 272 |  | 
|---|
| 273 | InEdgeIt e; | 
|---|
| 274 | for(G.first(e,v); G.valid(e); G.next(e)) { | 
|---|
| 275 | Node w=G.tail(e); | 
|---|
| 276 | if ( level[w] == number_of_nodes && w != s ) { | 
|---|
| 277 | bfs_queue.push(w); | 
|---|
| 278 | //Node first=level_list[l]; | 
|---|
| 279 | //if ( G.valid(first) ) left.set(first,w); | 
|---|
| 280 | //right.set(w,first); | 
|---|
| 281 | //level_list[l]=w; | 
|---|
| 282 | change_level_to(w, l); | 
|---|
| 283 | //level.set(w, l); | 
|---|
| 284 | } | 
|---|
| 285 | } | 
|---|
| 286 | } | 
|---|
| 287 | change_level_to(s,number_of_nodes); | 
|---|
| 288 | //level.set(s,number_of_nodes); | 
|---|
| 289 |  | 
|---|
| 290 | /* | 
|---|
| 291 | //Setting starting level values using reverse bfs | 
|---|
| 292 | reverse_bfs<Graph> rev_bfs(G,t); | 
|---|
| 293 | rev_bfs.run(); | 
|---|
| 294 | //write_property_vector(rev_bfs.dist,"rev_bfs"); | 
|---|
| 295 | for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) { | 
|---|
| 296 | change_level_to(i,rev_bfs.dist(i)); | 
|---|
| 297 | //level.put(i,rev_bfs.dist.get(i)); | 
|---|
| 298 | } | 
|---|
| 299 | */ | 
|---|
| 300 | //------------------------------------ | 
|---|
| 301 | //This is the only part that uses BFS | 
|---|
| 302 | //------------------------------------ | 
|---|
| 303 |  | 
|---|
| 304 |  | 
|---|
| 305 | //Starting level of s | 
|---|
| 306 | change_level_to(s,number_of_nodes); | 
|---|
| 307 | //level.put(s,number_of_nodes); | 
|---|
| 308 |  | 
|---|
| 309 |  | 
|---|
| 310 | //we push as much preflow from s as possible to start with | 
|---|
| 311 | for(OutEdgeIt j=G.template first<OutEdgeIt>(s); G.valid(j); G.next(j)){ | 
|---|
| 312 | modify_preflow(j,capacity[j] ); | 
|---|
| 313 | make_active(G.head(j)); | 
|---|
| 314 | int lev=level[G.head(j)]; | 
|---|
| 315 | if(highest_active<lev){ | 
|---|
| 316 | highest_active=lev; | 
|---|
| 317 | } | 
|---|
| 318 | } | 
|---|
| 319 | //cout<<highest_active<<endl; | 
|---|
| 320 | } | 
|---|
| 321 |  | 
|---|
| 322 |  | 
|---|
| 323 | //If the preflow is less than the capacity on the given edge | 
|---|
| 324 | //then it is an edge in the residual graph | 
|---|
| 325 | bool is_admissible_forward_edge(Edge j, int& new_level){ | 
|---|
| 326 |  | 
|---|
| 327 | if (capacity[j]>preflow[j]){ | 
|---|
| 328 | if(level[G.tail(j)]==level[G.head(j)]+1){ | 
|---|
| 329 | return true; | 
|---|
| 330 | } | 
|---|
| 331 | else{ | 
|---|
| 332 | if (level[G.head(j)] < new_level) | 
|---|
| 333 | new_level=level[G.head(j)]; | 
|---|
| 334 | } | 
|---|
| 335 | } | 
|---|
| 336 | return false; | 
|---|
| 337 | } | 
|---|
| 338 |  | 
|---|
| 339 | //If the preflow is greater than 0 on the given edge | 
|---|
| 340 | //then the edge reversd is an edge in the residual graph | 
|---|
| 341 | bool is_admissible_backward_edge(Edge j, int& new_level){ | 
|---|
| 342 |  | 
|---|
| 343 | if (0<preflow[j]){ | 
|---|
| 344 | if(level[G.tail(j)]==level[G.head(j)]-1){ | 
|---|
| 345 |  | 
|---|
| 346 | return true; | 
|---|
| 347 | } | 
|---|
| 348 | else{ | 
|---|
| 349 | if (level[G.tail(j)] < new_level) | 
|---|
| 350 | new_level=level[G.tail(j)]; | 
|---|
| 351 | } | 
|---|
| 352 |  | 
|---|
| 353 | } | 
|---|
| 354 | return false; | 
|---|
| 355 | } | 
|---|
| 356 |  | 
|---|
| 357 |  | 
|---|
| 358 | };  //class preflow_push | 
|---|
| 359 |  | 
|---|
| 360 | template<typename Graph, typename T> | 
|---|
| 361 | T preflow_push<Graph, T>::run() { | 
|---|
| 362 |  | 
|---|
| 363 | preprocess(); | 
|---|
| 364 | //write_property_vector(level,"level"); | 
|---|
| 365 | T e,v; | 
|---|
| 366 | Node a; | 
|---|
| 367 | while (a=get_active_node(), G.valid(a)){ | 
|---|
| 368 |  | 
|---|
| 369 | //cout<<G.id(a)<<endl; | 
|---|
| 370 | //write_property_vector(excess,"excess"); | 
|---|
| 371 | //write_property_vector(level,"level"); | 
|---|
| 372 |  | 
|---|
| 373 |  | 
|---|
| 374 | bool go_to_next_node=false; | 
|---|
| 375 | e = excess[a]; | 
|---|
| 376 | while (!go_to_next_node){ | 
|---|
| 377 | //Initial value for the new level for the active node we are dealing with | 
|---|
| 378 | int new_level=2*number_of_nodes; | 
|---|
| 379 | //write_property_vector(excess,"excess"); | 
|---|
| 380 | //write_property_vector(level,"level"); | 
|---|
| 381 | //cout<<G.id(a)<<endl; | 
|---|
| 382 | //Out edges from node a | 
|---|
| 383 | { | 
|---|
| 384 | OutEdgeIt j=G.template first<OutEdgeIt>(a); | 
|---|
| 385 | while (G.valid(j) && e){ | 
|---|
| 386 |  | 
|---|
| 387 | if (is_admissible_forward_edge(j,new_level)){ | 
|---|
| 388 | v=min(e,capacity[j] - preflow[j]); | 
|---|
| 389 | e -= v; | 
|---|
| 390 | //New node might become active | 
|---|
| 391 | if (excess[G.head(j)]==0){ | 
|---|
| 392 | make_active(G.head(j)); | 
|---|
| 393 | } | 
|---|
| 394 | modify_preflow(j,v); | 
|---|
| 395 | } | 
|---|
| 396 | G.next(j); | 
|---|
| 397 | } | 
|---|
| 398 | } | 
|---|
| 399 | //In edges to node a | 
|---|
| 400 | { | 
|---|
| 401 | InEdgeIt j=G.template first<InEdgeIt>(a); | 
|---|
| 402 | while (G.valid(j) && e){ | 
|---|
| 403 | if (is_admissible_backward_edge(j,new_level)){ | 
|---|
| 404 | v=min(e,preflow[j]); | 
|---|
| 405 | e -= v; | 
|---|
| 406 | //New node might become active | 
|---|
| 407 | if (excess[G.tail(j)]==0){ | 
|---|
| 408 | make_active(G.tail(j)); | 
|---|
| 409 | } | 
|---|
| 410 | modify_preflow(j,-v); | 
|---|
| 411 | } | 
|---|
| 412 | G.next(j); | 
|---|
| 413 | } | 
|---|
| 414 | } | 
|---|
| 415 |  | 
|---|
| 416 | //if (G.id(a)==999) | 
|---|
| 417 | //cout<<new_level<<" e: "<<e<<endl; | 
|---|
| 418 | //cout<<G.id(a)<<" "<<new_level<<endl; | 
|---|
| 419 |  | 
|---|
| 420 | if (0==e){ | 
|---|
| 421 | //Saturating push | 
|---|
| 422 | go_to_next_node=true; | 
|---|
| 423 | } | 
|---|
| 424 | else{//If there is still excess in node a | 
|---|
| 425 |  | 
|---|
| 426 | //change_level_to(a,new_level+1); | 
|---|
| 427 |  | 
|---|
| 428 | //Level remains empty | 
|---|
| 429 | if (num_of_nodes_on_level[level[a]]==1){ | 
|---|
| 430 | change_level_to(a,number_of_nodes); | 
|---|
| 431 | //go_to_next_node=True; | 
|---|
| 432 | } | 
|---|
| 433 | else{ | 
|---|
| 434 | change_level_to(a,new_level+1); | 
|---|
| 435 | //increase_level(a); | 
|---|
| 436 | } | 
|---|
| 437 |  | 
|---|
| 438 |  | 
|---|
| 439 |  | 
|---|
| 440 |  | 
|---|
| 441 | switch(node_examination){ | 
|---|
| 442 | case examine_to_relabel: | 
|---|
| 443 | make_active(a); | 
|---|
| 444 |  | 
|---|
| 445 | go_to_next_node = true; | 
|---|
| 446 | break; | 
|---|
| 447 | default: | 
|---|
| 448 | break; | 
|---|
| 449 | } | 
|---|
| 450 |  | 
|---|
| 451 |  | 
|---|
| 452 |  | 
|---|
| 453 | }//if (0==e) | 
|---|
| 454 | } | 
|---|
| 455 | } | 
|---|
| 456 | maxflow_value = excess[t]; | 
|---|
| 457 | return maxflow_value; | 
|---|
| 458 | }//run | 
|---|
| 459 |  | 
|---|
| 460 |  | 
|---|
| 461 | }//namespace hugo | 
|---|
| 462 |  | 
|---|
| 463 | #endif //PREFLOW_PUSH_HH | 
|---|