athos@331: #ifndef HUGO_PREFLOW_PUSH_HH
athos@331: #define HUGO_PREFLOW_PUSH_HH
athos@331: 
athos@331: //#include <algorithm>
athos@331: #include <list>
athos@331: #include <vector>
athos@331: #include <queue>
athos@331: //#include "pf_hiba.hh"
athos@331: //#include <marci_list_graph.hh>
athos@331: //#include <marci_graph_traits.hh>
athos@331: #include <invalid.h>
athos@331: //#include <reverse_bfs.hh>
athos@331: 
athos@331: using namespace std;
athos@331: 
athos@331: namespace hugo {
athos@331: 
athos@331:   template <typename Graph, typename T>
athos@331:   class preflow_push {
athos@331: 
athos@331:     //Useful typedefs
athos@331:     typedef typename Graph::Node Node;
athos@331:     typedef typename Graph::NodeIt NodeIt;
athos@331:     typedef typename Graph::Edge Edge;
athos@331:     typedef typename Graph::OutEdgeIt OutEdgeIt;
athos@331:     typedef typename Graph::InEdgeIt InEdgeIt;
athos@331: 
athos@331: 
athos@331:     //---------------------------------------------
athos@331:     //Parameters of the algorithm
athos@331:     //---------------------------------------------
athos@331:     //Fully examine an active node until excess becomes 0
athos@331:     enum node_examination_t {examine_full, examine_to_relabel};
athos@331:     //No more implemented yet:, examine_only_one_edge};
athos@331:     node_examination_t node_examination;
athos@331:     //Which implementation to be used
athos@331:     enum implementation_t {impl_fifo, impl_highest_label};
athos@331:     //No more implemented yet:};
athos@331:     implementation_t implementation;
athos@331:     //---------------------------------------------
athos@331:     //Parameters of the algorithm
athos@331:     //---------------------------------------------
athos@331:  
athos@331:   private:
athos@331:     //input
athos@331:     Graph& G;
athos@331:     Node s;
athos@331:     Node t;
athos@331:     typename Graph::EdgeMap<T> &capacity;
athos@331: 
athos@331:     //output
athos@331:     typename Graph::EdgeMap<T> preflow;
athos@331:     T maxflow_value;
athos@331:   
athos@331:     //auxiliary variables for computation
athos@331:     //The number of the nodes
athos@331:     int number_of_nodes;
athos@331:     //A nodemap for the level
athos@331:     typename Graph::NodeMap<int> level;
athos@331:     //A nodemap for the excess
athos@331:     typename Graph::NodeMap<T> excess;
athos@331:     
athos@331:     //Number of nodes on each level
athos@331:     vector<int> num_of_nodes_on_level;
athos@331:     
athos@331:     //For the FIFO implementation
athos@331:     list<Node> fifo_nodes;
athos@331:     //For 'highest label' implementation
athos@331:     int highest_active;
athos@331:     //int second_highest_active;
athos@331:     vector< list<Node> > active_nodes;
athos@331: 
athos@331:   public:
athos@331:   
athos@331:     //Constructing the object using the graph, source, sink and capacity vector
athos@331:     preflow_push(
athos@331: 		      Graph& _G, 
athos@331: 		      Node _s, 
athos@331: 		      Node _t, 
athos@331: 		      typename Graph::EdgeMap<T> & _capacity)
athos@331:       : G(_G), s(_s), t(_t), 
athos@331: 	capacity(_capacity), 
athos@331: 	preflow(_G),
athos@331: 	//Counting the number of nodes
athos@331: 	//number_of_nodes(count(G.first<EachNodeIt>())),
athos@331: 	number_of_nodes(G.nodeNum()),
athos@331: 
athos@331: 	level(_G),
athos@331: 	excess(_G)//,
athos@331:         // Default constructor: active_nodes()
athos@331:     { 
athos@331:       //Simplest parameter settings
athos@331:       node_examination = examine_full;//examine_to_relabel;//
athos@331:       //Which implementation to be usedexamine_full
athos@331:       implementation = impl_highest_label;//impl_fifo;
athos@331:  
athos@331:       //
athos@331:       num_of_nodes_on_level.resize(2*number_of_nodes-1);
athos@331:       num_of_nodes_on_level.clear();
athos@331: 
athos@331:       switch(implementation){
athos@331:       case impl_highest_label :{
athos@331: 	active_nodes.clear();
athos@331: 	active_nodes.resize(2*number_of_nodes-1);
athos@331: 	
athos@331: 	break;
athos@331:       }
athos@331:       default:
athos@331: 	break;
athos@331:       }
athos@331: 
athos@331:     }
athos@331: 
athos@331:     //Returns the value of a maximal flow 
athos@331:     T run();
athos@331:   
athos@331:     typename Graph::EdgeMap<T>  getmaxflow(){
athos@331:       return preflow;
athos@331:     }
athos@331: 
athos@331: 
athos@331:   private:
athos@331:     //For testing purposes only
athos@331:     //Lists the node_properties
athos@331:     void write_property_vector(typename Graph::NodeMap<T> a,
athos@331: 			       //node_property_vector<Graph, T> a, 
athos@331: 			       char* prop_name="property"){
athos@331:       for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
athos@331: 	cout<<"Node id.: "<<G.id(i)<<", "<<prop_name<<" value: "<<a[i]<<endl;
athos@331:       }
athos@331:       cout<<endl;
athos@331:     }
athos@331: 
athos@331:     //Modifies the excess of the node and makes sufficient changes
athos@331:     void modify_excess(const Node& a ,T v){
athos@331:       //T old_value=excess[a];
athos@331:       excess[a] += v;
athos@331:     }
athos@331:   
athos@331:     //This private procedure is supposed to modify the preflow on edge j
athos@331:     //by value v (which can be positive or negative as well) 
athos@331:     //and maintain the excess on the head and tail
athos@331:     //Here we do not check whether this is possible or not
athos@331:     void modify_preflow(Edge j, const T& v){
athos@331: 
athos@331:       //Modifiyng the edge
athos@331:       preflow[j] += v;
athos@331: 
athos@331: 
athos@331:       //Modifiyng the head
athos@331:       modify_excess(G.head(j),v);
athos@331: 	
athos@331:       //Modifiyng the tail
athos@331:       modify_excess(G.tail(j),-v);
athos@331: 
athos@331:     }
athos@331: 
athos@331:     //Gives the active node to work with 
athos@331:     //(depending on the implementation to be used)
athos@331:     Node get_active_node(){
athos@331:       
athos@331: 
athos@331:       switch(implementation) {
athos@331:       case impl_highest_label : {
athos@331: 
athos@331: 	//First need to find the highest label for which there's an active node
athos@331: 	while( highest_active>=0 && active_nodes[highest_active].empty() ){ 
athos@331: 	  --highest_active;
athos@331: 	}
athos@331: 
athos@331: 	if( highest_active>=0) {
athos@331: 	  
athos@331: 
athos@331: 	  Node a=active_nodes[highest_active].front();
athos@331: 	  active_nodes[highest_active].pop_front();
athos@331: 	  
athos@331: 	  return a;
athos@331: 	}
athos@331: 	else {
athos@331: 	  return INVALID;
athos@331: 	}
athos@331: 	
athos@331: 	break;
athos@331: 	
athos@331:       }
athos@331:       case impl_fifo : {
athos@331: 
athos@331: 	if( ! fifo_nodes.empty() ) {
athos@331: 	  Node a=fifo_nodes.front();
athos@331: 	  fifo_nodes.pop_front();
athos@331: 	  return a;
athos@331: 	}
athos@331: 	else {
athos@331: 	  return INVALID;
athos@331: 	}
athos@331: 	break;
athos@331:       }
athos@331:       }
athos@331:       //
athos@331:       return INVALID;
athos@331:     }
athos@331: 
athos@331:     //Puts node 'a' among the active nodes
athos@331:     void make_active(const Node& a){
athos@331:       //s and t never become active
athos@331:       if (a!=s && a!= t){
athos@331: 	switch(implementation){
athos@331: 	case impl_highest_label :
athos@331: 	  active_nodes[level[a]].push_back(a);
athos@331: 	  break;
athos@331: 	case impl_fifo :
athos@331: 	  fifo_nodes.push_back(a);
athos@331: 	  break;
athos@331: 	}
athos@331: 
athos@331:       }
athos@331: 
athos@331:       //Update highest_active label
athos@331:       if (highest_active<level[a]){
athos@331: 	highest_active=level[a];
athos@331:       }
athos@331: 
athos@331:     }
athos@331: 
athos@331:     //Changes the level of node a and make sufficent changes
athos@331:     void change_level_to(Node a, int new_value){
athos@331:       int seged = level[a];
athos@331:       level.set(a,new_value);
athos@331:       --num_of_nodes_on_level[seged];
athos@331:       ++num_of_nodes_on_level[new_value];
athos@331:     }
athos@331: 
athos@331:     //Collection of things useful (or necessary) to do before running
athos@331: 
athos@331:     void preprocess(){
athos@331: 
athos@331:       //---------------------------------------
athos@331:       //Initialize parameters
athos@331:       //---------------------------------------
athos@331: 
athos@331:       //Setting starting preflow, level and excess values to zero
athos@331:       //This can be important, if the algorithm is run more then once
athos@331:       for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
athos@331:         level.set(i,0);
athos@331:         excess.set(i,0);
athos@331: 	for(OutEdgeIt j=G.template first<OutEdgeIt>(i); G.valid(j); G.next(j)) 
athos@331: 	  preflow.set(j, 0);
athos@331:       }
athos@331:       num_of_nodes_on_level[0]=number_of_nodes;
athos@331:       highest_active=0;
athos@331:       //---------------------------------------
athos@331:       //Initialize parameters
athos@331:       //---------------------------------------
athos@331: 
athos@331:       
athos@331:       //------------------------------------
athos@331:       //This is the only part that uses BFS
athos@331:       //------------------------------------
athos@331: 
athos@331:       /*Reverse_bfs from t, to find the starting level.*/
athos@331:       //Copyright: Jacint
athos@331:       change_level_to(t,0);
athos@331: 
athos@331:       std::queue<Node> bfs_queue;
athos@331:       bfs_queue.push(t);
athos@331: 
athos@331:       while (!bfs_queue.empty()) {
athos@331: 
athos@331: 	Node v=bfs_queue.front();	
athos@331: 	bfs_queue.pop();
athos@331: 	int l=level[v]+1;
athos@331: 
athos@331: 	InEdgeIt e;
athos@331: 	for(G.first(e,v); G.valid(e); G.next(e)) {
athos@331: 	  Node w=G.tail(e);
athos@331: 	  if ( level[w] == number_of_nodes && w != s ) {
athos@331: 	    bfs_queue.push(w);
athos@331: 	    //Node first=level_list[l];
athos@331: 	    //if ( G.valid(first) ) left.set(first,w);
athos@331: 	    //right.set(w,first);
athos@331: 	    //level_list[l]=w;
athos@331: 	    change_level_to(w, l);
athos@331: 	    //level.set(w, l);
athos@331: 	  }
athos@331: 	}
athos@331:       }
athos@331:       change_level_to(s,number_of_nodes);
athos@331:       //level.set(s,number_of_nodes);
athos@331: 
athos@331:       /*
athos@331:       //Setting starting level values using reverse bfs
athos@331:       reverse_bfs<Graph> rev_bfs(G,t);
athos@331:       rev_bfs.run();
athos@331:       //write_property_vector(rev_bfs.dist,"rev_bfs");
athos@331:       for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
athos@331:         change_level_to(i,rev_bfs.dist(i));
athos@331: 	//level.put(i,rev_bfs.dist.get(i));
athos@331:       }
athos@331:       */
athos@331:       //------------------------------------
athos@331:       //This is the only part that uses BFS
athos@331:       //------------------------------------
athos@331:       
athos@331:       
athos@331:       //Starting level of s
athos@331:       change_level_to(s,number_of_nodes);
athos@331:       //level.put(s,number_of_nodes);
athos@331:       
athos@331:       
athos@331:       //we push as much preflow from s as possible to start with
athos@331:       for(OutEdgeIt j=G.template first<OutEdgeIt>(s); G.valid(j); G.next(j)){ 
athos@331: 	modify_preflow(j,capacity[j] );
athos@331: 	make_active(G.head(j));
athos@331: 	int lev=level[G.head(j)];
athos@331: 	if(highest_active<lev){
athos@331: 	  highest_active=lev;
athos@331: 	}
athos@331:       }
athos@331:       //cout<<highest_active<<endl;
athos@331:     } 
athos@331: 
athos@331:     
athos@331:     //If the preflow is less than the capacity on the given edge
athos@331:     //then it is an edge in the residual graph
athos@331:     bool is_admissible_forward_edge(Edge j, int& new_level){
athos@331: 
athos@331:       if (capacity[j]>preflow[j]){
athos@331: 	if(level[G.tail(j)]==level[G.head(j)]+1){
athos@331: 	  return true;
athos@331: 	}
athos@331: 	else{
athos@331: 	  if (level[G.head(j)] < new_level)
athos@331: 	    new_level=level[G.head(j)];
athos@331: 	}
athos@331:       }
athos@331:       return false;
athos@331:     }
athos@331: 
athos@331:     //If the preflow is greater than 0 on the given edge
athos@331:     //then the edge reversd is an edge in the residual graph
athos@331:     bool is_admissible_backward_edge(Edge j, int& new_level){
athos@331:       
athos@331:       if (0<preflow[j]){
athos@331: 	if(level[G.tail(j)]==level[G.head(j)]-1){
athos@331: 	 
athos@331: 	  return true;
athos@331: 	}
athos@331: 	else{
athos@331: 	  if (level[G.tail(j)] < new_level)
athos@331: 	    new_level=level[G.tail(j)];
athos@331: 	}
athos@331: 	
athos@331:       }
athos@331:       return false;
athos@331:     }
athos@331: 
athos@331:  
athos@331:   };  //class preflow_push  
athos@331: 
athos@331:   template<typename Graph, typename T>
athos@331:     T preflow_push<Graph, T>::run() {
athos@331:     
athos@331:     preprocess();
athos@331:     //write_property_vector(level,"level");
athos@331:     T e,v;
athos@331:     Node a;
athos@331:     while (a=get_active_node(), G.valid(a)){
athos@331:       
athos@331:       //cout<<G.id(a)<<endl;
athos@331:       //write_property_vector(excess,"excess");
athos@331:       //write_property_vector(level,"level");
athos@331: 
athos@331: 
athos@331:       bool go_to_next_node=false;
athos@331:       e = excess[a];
athos@331:       while (!go_to_next_node){
athos@331: 	//Initial value for the new level for the active node we are dealing with
athos@331: 	int new_level=2*number_of_nodes;
athos@331: 	//write_property_vector(excess,"excess");
athos@331: 	//write_property_vector(level,"level");
athos@331: 	//cout<<G.id(a)<<endl;
athos@331: 	//Out edges from node a
athos@331: 	{
athos@331: 	  OutEdgeIt j=G.template first<OutEdgeIt>(a);
athos@331: 	  while (G.valid(j) && e){
athos@331: 
athos@331: 	    if (is_admissible_forward_edge(j,new_level)){
athos@331: 	      v=min(e,capacity[j] - preflow[j]);
athos@331: 	      e -= v;
athos@331: 	      //New node might become active
athos@331: 	      if (excess[G.head(j)]==0){
athos@331: 		make_active(G.head(j));
athos@331: 	      }
athos@331: 	      modify_preflow(j,v);
athos@331: 	    }
athos@331: 	    G.next(j);
athos@331: 	  }
athos@331: 	}
athos@331: 	//In edges to node a
athos@331: 	{
athos@331: 	  InEdgeIt j=G.template first<InEdgeIt>(a);
athos@331: 	  while (G.valid(j) && e){
athos@331: 	    if (is_admissible_backward_edge(j,new_level)){
athos@331: 	      v=min(e,preflow[j]);
athos@331: 	      e -= v;
athos@331: 	      //New node might become active
athos@331: 	      if (excess[G.tail(j)]==0){
athos@331: 		make_active(G.tail(j));
athos@331: 	      }
athos@331: 	      modify_preflow(j,-v);
athos@331: 	    }
athos@331: 	    G.next(j);
athos@331: 	  }
athos@331: 	}
athos@331: 
athos@331: 	//if (G.id(a)==999)
athos@331: 	//cout<<new_level<<" e: "<<e<<endl;
athos@331: 	//cout<<G.id(a)<<" "<<new_level<<endl;
athos@331: 
athos@331: 	if (0==e){
athos@331: 	  //Saturating push
athos@331: 	  go_to_next_node=true;
athos@331: 	}
athos@331: 	else{//If there is still excess in node a
athos@331: 	  
athos@331: 	  //change_level_to(a,new_level+1);
athos@331: 	  
athos@331: 	  //Level remains empty
athos@331: 	  if (num_of_nodes_on_level[level[a]]==1){
athos@331: 	    change_level_to(a,number_of_nodes);
athos@331: 	    //go_to_next_node=True;
athos@331: 	  }
athos@331: 	  else{
athos@331: 	    change_level_to(a,new_level+1);
athos@331: 	    //increase_level(a);
athos@331: 	  }
athos@331: 	  
athos@331:     
athos@331: 	  
athos@331: 
athos@331: 	  switch(node_examination){
athos@331: 	  case examine_to_relabel:
athos@331: 	    make_active(a);
athos@331: 
athos@331: 	    go_to_next_node = true;
athos@331: 	    break;
athos@331: 	  default:
athos@331: 	    break;
athos@331: 	  }
athos@331: 	  
athos@331:     
athos@331: 	
athos@331: 	}//if (0==e)
athos@331:       }
athos@331:     }
athos@331:     maxflow_value = excess[t];
athos@331:     return maxflow_value;
athos@331:   }//run
athos@331: 
athos@331: 
athos@331: }//namespace hugo
athos@331: 
athos@331: #endif //PREFLOW_PUSH_HH