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