source: lemon-0.x/src/work/athos/preflow_push.hh @ 631:26819ef1611f

#ifndef HUGO_PREFLOW_PUSH_HH
#define HUGO_PREFLOW_PUSH_HH

//#include <algorithm>
#include <list>
#include <vector>
#include <queue>
//#include "pf_hiba.hh"
//#include <marci_list_graph.hh>
//#include <marci_graph_traits.hh>
#include <invalid.h>
#include <graph_wrapper.h>
//#include <reverse_bfs.hh>

using namespace std;

namespace hugo {

  template <typename Graph, typename T>
  class preflow_push {

    //Useful typedefs
    typedef typename Graph::Node Node;
    typedef typename Graph::NodeIt NodeIt;
    typedef typename Graph::Edge Edge;
    typedef typename Graph::OutEdgeIt OutEdgeIt;
    typedef typename Graph::InEdgeIt InEdgeIt;
    typedef typename Graph::EdgeMap<T> CapacityType;

    typedef ResGraphWrapper<const Graph,int,CapacityType,CapacityType> ResGraphType;


    //---------------------------------------------
    //Parameters of the algorithm
    //---------------------------------------------
    //Fully examine an active node until excess becomes 0
    enum node_examination_t {examine_full, examine_to_relabel};
    //No more implemented yet:, examine_only_one_edge};
    node_examination_t node_examination;
    //Which implementation to be used
    enum implementation_t {impl_fifo, impl_highest_label};
    //No more implemented yet:};
    implementation_t implementation;
    //---------------------------------------------
    //Parameters of the algorithm
    //---------------------------------------------
 
  private:
    //input
    Graph& G;
    Node s;
    Node t;
    CapacityType &capacity;

    //output
    CapacityType preflow;
    T maxflow_value;
  
    //auxiliary variables for computation
    //The number of the nodes
    int number_of_nodes;
    //A nodemap for the level
    typename Graph::NodeMap<int> level;
    //A nodemap for the excess
    typename Graph::NodeMap<T> excess;
    
    //Number of nodes on each level
    vector<int> num_of_nodes_on_level;
    
    //For the FIFO implementation
    list<Node> fifo_nodes;
    //For 'highest label' implementation
    int highest_active;
    //int second_highest_active;
    vector< list<Node> > active_nodes;

  public:
  
    //Constructing the object using the graph, source, sink and capacity vector
    preflow_push(
                      Graph& _G, 
                      Node _s, 
                      Node _t, 
                      typename Graph::EdgeMap<T> & _capacity)
      : G(_G), s(_s), t(_t), 
        capacity(_capacity), 
        preflow(_G),
        //Counting the number of nodes
        //number_of_nodes(count(G.first<EachNodeIt>())),
        number_of_nodes(G.nodeNum()),

        level(_G),
        excess(_G)//,
        // Default constructor: active_nodes()
    { 
      //Simplest parameter settings
      node_examination = examine_full;//examine_to_relabel;//
      //Which implementation to be usedexamine_full
      implementation = impl_highest_label;//impl_fifo;
 
      //
      num_of_nodes_on_level.resize(2*number_of_nodes-1);
      num_of_nodes_on_level.clear();

      switch(implementation){
      case impl_highest_label :{
        active_nodes.clear();
        active_nodes.resize(2*number_of_nodes-1);
        
        break;
      }
      default:
        break;
      }

    }

    //Returns the value of a maximal flow 
    T run();
  
    typename Graph::EdgeMap<T>  getmaxflow(){
      return preflow;
    }


  private:
    //For testing purposes only
    //Lists the node_properties
    void write_property_vector(typename Graph::NodeMap<T> a,
                               //node_property_vector<Graph, T> a, 
                               char* prop_name="property"){
      for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
        cout<<"Node id.: "<<G.id(i)<<", "<<prop_name<<" value: "<<a[i]<<endl;
      }
      cout<<endl;
    }
    /*
    //Modifies the excess of the node and makes sufficient changes
    void modify_excess(const Node& a ,T v){
      //T old_value=excess[a];
      excess[a] += v;
    }
  
    //This private procedure is supposed to modify the preflow on edge j
    //by value v (which can be positive or negative as well) 
    //and maintain the excess on the head and tail
    //Here we do not check whether this is possible or not
    void modify_preflow(Edge j, const T& v){

      //Modifiyng the edge
      preflow[j] += v;


      //Modifiyng the head
      modify_excess(G.head(j),v);
        
      //Modifiyng the tail
      modify_excess(G.tail(j),-v);

    }
    */
    //Gives the active node to work with 
    //(depending on the implementation to be used)
    Node get_active_node(){
      

      switch(implementation) {
      case impl_highest_label : {

        //First need to find the highest label for which there's an active node
        while( highest_active>=0 && active_nodes[highest_active].empty() ){ 
          --highest_active;
        }

        if( highest_active>=0) {
          

          Node a=active_nodes[highest_active].front();
          active_nodes[highest_active].pop_front();
          
          return a;
        }
        else {
          return INVALID;
        }
        
        break;
        
      }
      case impl_fifo : {

        if( ! fifo_nodes.empty() ) {
          Node a=fifo_nodes.front();
          fifo_nodes.pop_front();
          return a;
        }
        else {
          return INVALID;
        }
        break;
      }
      }
      //
      return INVALID;
    }

    //Puts node 'a' among the active nodes
    void make_active(const Node& a){
      //s and t never become active
      if (a!=s && a!= t){
        switch(implementation){
        case impl_highest_label :
          active_nodes[level[a]].push_back(a);
          break;
        case impl_fifo :
          fifo_nodes.push_back(a);
          break;
        }

      }

      //Update highest_active label
      if (highest_active<level[a]){
        highest_active=level[a];
      }

    }

    //Changes the level of node a and make sufficent changes
    void change_level_to(Node a, int new_value){
      int seged = level[a];
      level.set(a,new_value);
      --num_of_nodes_on_level[seged];
      ++num_of_nodes_on_level[new_value];
    }

    //Collection of things useful (or necessary) to do before running

    void preprocess(){

      //---------------------------------------
      //Initialize parameters
      //---------------------------------------

      //Setting starting preflow, level and excess values to zero
      //This can be important, if the algorithm is run more then once
      for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
        level.set(i,0);
        excess.set(i,0);
        for(OutEdgeIt j=G.template first<OutEdgeIt>(i); G.valid(j); G.next(j)) 
          preflow.set(j, 0);
      }
      num_of_nodes_on_level[0]=number_of_nodes;
      highest_active=0;
      //---------------------------------------
      //Initialize parameters
      //---------------------------------------

      
      //------------------------------------
      //This is the only part that uses BFS
      //------------------------------------

      /*Reverse_bfs from t, to find the starting level.*/
      //Copyright: Jacint
      change_level_to(t,0);

      std::queue<Node> bfs_queue;
      bfs_queue.push(t);

      while (!bfs_queue.empty()) {

        Node v=bfs_queue.front();       
        bfs_queue.pop();
        int l=level[v]+1;

        InEdgeIt e;
        for(G.first(e,v); G.valid(e); G.next(e)) {
          Node w=G.tail(e);
          if ( level[w] == number_of_nodes && w != s ) {
            bfs_queue.push(w);
            //Node first=level_list[l];
            //if ( G.valid(first) ) left.set(first,w);
            //right.set(w,first);
            //level_list[l]=w;
            change_level_to(w, l);
            //level.set(w, l);
          }
        }
      }
      change_level_to(s,number_of_nodes);
      //level.set(s,number_of_nodes);

      /*
      //Setting starting level values using reverse bfs
      reverse_bfs<Graph> rev_bfs(G,t);
      rev_bfs.run();
      //write_property_vector(rev_bfs.dist,"rev_bfs");
      for(NodeIt i=G.template first<NodeIt>(); G.valid(i); G.next(i)) {
        change_level_to(i,rev_bfs.dist(i));
        //level.put(i,rev_bfs.dist.get(i));
      }
      */
      //------------------------------------
      //This is the only part that uses BFS
      //------------------------------------
      
      
      //Starting level of s
      change_level_to(s,number_of_nodes);
      //level.put(s,number_of_nodes);
      
      
      //we push as much preflow from s as possible to start with
      for(OutEdgeIt j=G.template first<OutEdgeIt>(s); G.valid(j); G.next(j)){ 
        modify_preflow(j,capacity[j] );
        make_active(G.head(j));
        int lev=level[G.head(j)];
        if(highest_active<lev){
          highest_active=lev;
        }
      }
      //cout<<highest_active<<endl;
    } 

    
    //If the preflow is less than the capacity on the given edge
    //then it is an edge in the residual graph
    bool is_admissible_forward_edge(Edge j, int& new_level){

      if (capacity[j]>preflow[j]){
        if(level[G.tail(j)]==level[G.head(j)]+1){
          return true;
        }
        else{
          if (level[G.head(j)] < new_level)
            new_level=level[G.head(j)];
        }
      }
      return false;
    }

    //If the preflow is greater than 0 on the given edge
    //then the edge reversd is an edge in the residual graph
    bool is_admissible_backward_edge(Edge j, int& new_level){
      
      if (0<preflow[j]){
        if(level[G.tail(j)]==level[G.head(j)]-1){
         
          return true;
        }
        else{
          if (level[G.tail(j)] < new_level)
            new_level=level[G.tail(j)];
        }
        
      }
      return false;
    }

 
  };  //class preflow_push  

  template<typename Graph, typename T>
    T preflow_push<Graph, T>::run() {
    
    //We need a residual graph
    ResGraphType res_graph(G, preflow, capacity);
    
    preprocess();
    //write_property_vector(level,"level");
    T e,v;
    Node a;
    while (a=get_active_node(), G.valid(a)){
      
      bool go_to_next_node=false;
      e = excess[a];
      while (!go_to_next_node){

        //Initial value for the new level for the active node we are dealing with
        int new_level=2*number_of_nodes;


        //Out edges from node a
        {
          ResGraphType::OutEdgeIt j=res_graph.first(j,a);
          while (res_graph.valid(j) && e){
            if (is_admissible_forward_edge(j,new_level)){
              v=min(e,res_graph.resCap(j));
              e -= v;
              //New node might become active
              if (excess[res_graph.head(j)]==0){
                make_active(res_graph.head(j));
              }
              res_graph.augment(j,v);
              excess[res_graph.tail(j)] -= v;
              excess[res_graph.head(j)] += v;
            }
            res_graph.next(j);
          }
        }

        /*
        //Out edges from node a
        {
          OutEdgeIt j=G.template first<OutEdgeIt>(a);
          while (G.valid(j) && e){

            if (is_admissible_forward_edge(j,new_level)){
              v=min(e,capacity[j] - preflow[j]);
              e -= v;
              //New node might become active
              if (excess[G.head(j)]==0){
                make_active(G.head(j));
              }
              modify_preflow(j,v);
            }
            G.next(j);
          }
        }
        //In edges to node a
        {
          InEdgeIt j=G.template first<InEdgeIt>(a);
          while (G.valid(j) && e){
            if (is_admissible_backward_edge(j,new_level)){
              v=min(e,preflow[j]);
              e -= v;
              //New node might become active
              if (excess[G.tail(j)]==0){
                make_active(G.tail(j));
              }
              modify_preflow(j,-v);
            }
            G.next(j);
          }
        }
        */

        //if (G.id(a)==999)
        //cout<<new_level<<" e: "<<e<<endl;
        //cout<<G.id(a)<<" "<<new_level<<endl;

        if (0==e){
          //Saturating push
          go_to_next_node=true;
        }
        else{//If there is still excess in node a
          
          //change_level_to(a,new_level+1);
          
          //Level remains empty
          if (num_of_nodes_on_level[level[a]]==1){
            change_level_to(a,number_of_nodes);
            //go_to_next_node=True;
          }
          else{
            change_level_to(a,new_level+1);
            //increase_level(a);
          }
          
    
          

          switch(node_examination){
          case examine_to_relabel:
            make_active(a);

            go_to_next_node = true;
            break;
          default:
            break;
          }
          
    
        
        }//if (0==e)
      }
    }
    maxflow_value = excess[t];
    return maxflow_value;
  }//run


}//namespace hugo

#endif //PREFLOW_PUSH_HH