COIN-OR::LEMON - Graph Library

source: lemon-0.x/src/work/athos/preflow_push.hh @ 547:50184b822370

Last change on this file since 547:50184b822370 was 505:8589c0658839, checked in by athos, 20 years ago

I changed it to correspond changing requirements

File size: 11.4 KB
RevLine 
[331]1#ifndef HUGO_PREFLOW_PUSH_HH
2#define HUGO_PREFLOW_PUSH_HH
[36]3
[331]4//#include <algorithm>
[36]5#include <list>
6#include <vector>
[331]7#include <queue>
[36]8//#include "pf_hiba.hh"
9//#include <marci_list_graph.hh>
[77]10//#include <marci_graph_traits.hh>
[331]11#include <invalid.h>
12#include <graph_wrapper.h>
13//#include <reverse_bfs.hh>
[36]14
15using namespace std;
16
[105]17namespace hugo {
[36]18
[331]19  template <typename Graph, typename T>
[36]20  class preflow_push {
21
[331]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;
[505]28    typedef typename Graph::EdgeMap<T> CapacityType;
29
30    typedef ResGraphWrapper<const Graph,int,CapacityType,CapacityType> ResGraphType;
[77]31
32
[36]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
[331]50    Graph& G;
51    Node s;
52    Node t;
[505]53    CapacityType &capacity;
[331]54
[36]55    //output
[505]56    CapacityType preflow;
[36]57    T maxflow_value;
58 
59    //auxiliary variables for computation
[331]60    //The number of the nodes
[36]61    int number_of_nodes;
[331]62    //A nodemap for the level
63    typename Graph::NodeMap<int> level;
64    //A nodemap for the excess
65    typename Graph::NodeMap<T> excess;
[36]66   
67    //Number of nodes on each level
68    vector<int> num_of_nodes_on_level;
69   
70    //For the FIFO implementation
[331]71    list<Node> fifo_nodes;
[36]72    //For 'highest label' implementation
73    int highest_active;
74    //int second_highest_active;
[331]75    vector< list<Node> > active_nodes;
[36]76
77  public:
78 
79    //Constructing the object using the graph, source, sink and capacity vector
80    preflow_push(
[331]81                      Graph& _G,
82                      Node _s,
83                      Node _t,
84                      typename Graph::EdgeMap<T> & _capacity)
[36]85      : G(_G), s(_s), t(_t),
86        capacity(_capacity),
87        preflow(_G),
88        //Counting the number of nodes
[77]89        //number_of_nodes(count(G.first<EachNodeIt>())),
90        number_of_nodes(G.nodeNum()),
91
[36]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 :{
[331]107        active_nodes.clear();
[36]108        active_nodes.resize(2*number_of_nodes-1);
[331]109       
[36]110        break;
111      }
112      default:
113        break;
114      }
115
116    }
117
118    //Returns the value of a maximal flow
119    T run();
120 
[331]121    typename Graph::EdgeMap<T>  getmaxflow(){
[36]122      return preflow;
123    }
124
125
126  private:
127    //For testing purposes only
128    //Lists the node_properties
[331]129    void write_property_vector(typename Graph::NodeMap<T> a,
130                               //node_property_vector<Graph, T> a,
[36]131                               char* prop_name="property"){
[331]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;
[36]134      }
135      cout<<endl;
136    }
[505]137    /*
[36]138    //Modifies the excess of the node and makes sufficient changes
[331]139    void modify_excess(const Node& a ,T v){
140      //T old_value=excess[a];
141      excess[a] += v;
[36]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
[331]148    void modify_preflow(Edge j, const T& v){
[36]149
150      //Modifiyng the edge
[331]151      preflow[j] += v;
[36]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    }
[505]161    */
[36]162    //Gives the active node to work with
163    //(depending on the implementation to be used)
[331]164    Node get_active_node(){
[119]165     
[36]166
167      switch(implementation) {
168      case impl_highest_label : {
169
[331]170        //First need to find the highest label for which there's an active node
[36]171        while( highest_active>=0 && active_nodes[highest_active].empty() ){
172          --highest_active;
173        }
174
175        if( highest_active>=0) {
[119]176         
177
[331]178          Node a=active_nodes[highest_active].front();
[36]179          active_nodes[highest_active].pop_front();
[119]180         
[36]181          return a;
182        }
183        else {
[331]184          return INVALID;
[36]185        }
186       
187        break;
188       
189      }
190      case impl_fifo : {
191
192        if( ! fifo_nodes.empty() ) {
[331]193          Node a=fifo_nodes.front();
[36]194          fifo_nodes.pop_front();
195          return a;
196        }
197        else {
[331]198          return INVALID;
[36]199        }
200        break;
201      }
202      }
203      //
[331]204      return INVALID;
[36]205    }
206
207    //Puts node 'a' among the active nodes
[331]208    void make_active(const Node& a){
[36]209      //s and t never become active
210      if (a!=s && a!= t){
211        switch(implementation){
212        case impl_highest_label :
[331]213          active_nodes[level[a]].push_back(a);
[36]214          break;
215        case impl_fifo :
216          fifo_nodes.push_back(a);
217          break;
218        }
219
220      }
221
222      //Update highest_active label
[331]223      if (highest_active<level[a]){
224        highest_active=level[a];
[36]225      }
226
227    }
228
229    //Changes the level of node a and make sufficent changes
[331]230    void change_level_to(Node a, int new_value){
231      int seged = level[a];
[77]232      level.set(a,new_value);
[36]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
[77]238
[36]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
[331]247      for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
[77]248        level.set(i,0);
249        excess.set(i,0);
[331]250        for(OutEdgeIt j=G.template first<OutEdgeIt>(i); G.valid(j); G.next(j))
[77]251          preflow.set(j, 0);
[36]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      //------------------------------------
[331]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      /*
[36]295      //Setting starting level values using reverse bfs
[331]296      reverse_bfs<Graph> rev_bfs(G,t);
[36]297      rev_bfs.run();
298      //write_property_vector(rev_bfs.dist,"rev_bfs");
[331]299      for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
[36]300        change_level_to(i,rev_bfs.dist(i));
301        //level.put(i,rev_bfs.dist.get(i));
302      }
[331]303      */
[36]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
[331]315      for(OutEdgeIt j=G.template first<OutEdgeIt>(s); G.valid(j); G.next(j)){
316        modify_preflow(j,capacity[j] );
[36]317        make_active(G.head(j));
[331]318        int lev=level[G.head(j)];
[36]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
[331]329    bool is_admissible_forward_edge(Edge j, int& new_level){
[119]330
[331]331      if (capacity[j]>preflow[j]){
332        if(level[G.tail(j)]==level[G.head(j)]+1){
[36]333          return true;
334        }
335        else{
[331]336          if (level[G.head(j)] < new_level)
337            new_level=level[G.head(j)];
[36]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
[331]345    bool is_admissible_backward_edge(Edge j, int& new_level){
[119]346     
[331]347      if (0<preflow[j]){
348        if(level[G.tail(j)]==level[G.head(j)]-1){
[119]349         
[36]350          return true;
351        }
352        else{
[331]353          if (level[G.tail(j)] < new_level)
354            new_level=level[G.tail(j)];
[36]355        }
356       
357      }
358      return false;
359    }
360
361 
362  };  //class preflow_push 
363
[331]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);
[36]369   
370    preprocess();
[119]371    //write_property_vector(level,"level");
[36]372    T e,v;
[331]373    Node a;
374    while (a=get_active_node(), G.valid(a)){
[119]375     
[331]376      bool go_to_next_node=false;
377      e = excess[a];
378      while (!go_to_next_node){
[36]379
[77]380        //Initial value for the new level for the active node we are dealing with
381        int new_level=2*number_of_nodes;
[331]382
383
[36]384        //Out edges from node a
385        {
[505]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        {
[77]406          OutEdgeIt j=G.template first<OutEdgeIt>(a);
[331]407          while (G.valid(j) && e){
[36]408
409            if (is_admissible_forward_edge(j,new_level)){
[331]410              v=min(e,capacity[j] - preflow[j]);
[36]411              e -= v;
412              //New node might become active
[331]413              if (excess[G.head(j)]==0){
[36]414                make_active(G.head(j));
415              }
416              modify_preflow(j,v);
417            }
[331]418            G.next(j);
[36]419          }
420        }
421        //In edges to node a
422        {
[77]423          InEdgeIt j=G.template first<InEdgeIt>(a);
[331]424          while (G.valid(j) && e){
[36]425            if (is_admissible_backward_edge(j,new_level)){
[331]426              v=min(e,preflow[j]);
[36]427              e -= v;
428              //New node might become active
[331]429              if (excess[G.tail(j)]==0){
[36]430                make_active(G.tail(j));
431              }
432              modify_preflow(j,-v);
433            }
[331]434            G.next(j);
[36]435          }
436        }
[505]437        */
[36]438
[119]439        //if (G.id(a)==999)
440        //cout<<new_level<<" e: "<<e<<endl;
[36]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
[77]448         
449          //change_level_to(a,new_level+1);
450         
[36]451          //Level remains empty
[331]452          if (num_of_nodes_on_level[level[a]]==1){
[36]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          }
[77]460         
[36]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    }
[331]479    maxflow_value = excess[t];
[36]480    return maxflow_value;
481  }//run
482
483
[105]484}//namespace hugo
[36]485
486#endif //PREFLOW_PUSH_HH
Note: See TracBrowser for help on using the repository browser.