COIN-OR::LEMON - Graph Library

Ignore:
Timestamp:
02/18/04 15:42:38 (16 years ago)
Author:
jacint
Branch:
default
Phase:
public
Convert:
svn:c9d7d8f5-90d6-0310-b91f-818b3a526b0e/lemon/trunk@126
Message:

javitott valtozat

File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/work/jacint/preflow_push_max_flow.h

    r85 r97  
     1// -*- C++ -*-
    12/*
    23preflow_push_max_flow_h
     
    1617T maxflow() : returns the value of a maximum flow
    1718
    18 NodeMap<bool> mincut(): returns a
    19      characteristic vector of a minimum cut.
     19void mincut(CutMap& M) : sets M to the characteristic vector of a
     20     minimum cut. M should be a map of bools initialized to false.
     21
    2022*/
    2123
    2224#ifndef PREFLOW_PUSH_MAX_FLOW_H
    2325#define PREFLOW_PUSH_MAX_FLOW_H
     26
     27#define A 1
    2428
    2529#include <algorithm>
     
    3236namespace marci {
    3337
    34   template <typename Graph, typename T>
     38  template <typename Graph, typename T, 
     39    typename FlowMap=typename Graph::EdgeMap<T>, typename CapMap=typename Graph::EdgeMap<T>,
     40    typename IntMap=typename Graph::NodeMap<int>, typename TMap=typename Graph::NodeMap<T> >
    3541  class preflow_push_max_flow {
    3642   
     
    4349    NodeIt s;
    4450    NodeIt t;
    45     typename Graph::EdgeMap<T>& capacity;
    46     T value;
    47     typename Graph::NodeMap<bool> mincutvector;   
    48 
    49 
    50      
     51    IntMap level;
     52    CapMap& capacity; 
     53    int empty_level;    //an empty level in the end of run()
     54    T value;
     55   
    5156  public:
    52        
    53     preflow_push_max_flow ( Graph& _G, NodeIt _s, NodeIt _t,
    54                             typename Graph::EdgeMap<T>& _capacity) :
    55       G(_G), s(_s), t(_t), capacity(_capacity), mincutvector(_G, false) { }
     57     
     58    preflow_push_max_flow(Graph& _G, NodeIt _s, NodeIt _t, CapMap& _capacity) :
     59      G(_G), s(_s), t(_t), level(_G), capacity(_capacity) { }
    5660
    5761
     
    6367    void run() {
    6468 
    65       typename Graph::EdgeMap<T> flow(G, 0);
    66       typename Graph::NodeMap<int> level(G);   
    67       typename Graph::NodeMap<T> excess(G);   
    68            
    69       int n=G.nodeNum();                       
     69      int n=G.nodeNum();
    7070      int b=n-2;
    7171      /*
    72         b is a bound on the highest level of an active Node.
    73         In the beginning it is at most n-2.
    74       */
    75      
    76       std::vector<int> numb(n);     //The number of Nodes on level i < n.
    77       std::vector<std::stack<NodeIt> > stack(2*n-1);   
    78       //Stack of the active Nodes in level i.
     72        b is a bound on the highest level of an active node.
     73      */
     74
     75      IntMap level(G,n);     
     76      TMap excess(G);
     77      FlowMap flow(G,0);
     78
     79      std::vector<int> numb(n);   
     80      /*
     81        The number of nodes on level i < n. It is
     82        initialized to n+1, because of the reverse_bfs-part.
     83      */
     84
     85      std::vector<std::stack<NodeIt> > stack(n);   
     86      //Stack of the active nodes in level i.
     87
    7988
    8089      /*Reverse_bfs from t, to find the starting level.*/
    81       reverse_bfs<Graph> bfs(G, t);
    82       bfs.run();
    83       for(EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v)
    84         {
    85           int dist=bfs.dist(v);
    86           level.set(v, dist);
    87           ++numb[dist];
     90      level.set(t,0);
     91      std::queue<NodeIt> bfs_queue;
     92      bfs_queue.push(t);
     93
     94      while (!bfs_queue.empty()) {
     95
     96        NodeIt v=bfs_queue.front();     
     97        bfs_queue.pop();
     98        int l=level.get(v)+1;
     99
     100        for(InEdgeIt e=G.template first<InEdgeIt>(v); e.valid(); ++e) {
     101          NodeIt w=G.tail(e);
     102          if ( level.get(w) == n ) {
     103            bfs_queue.push(w);
     104            ++numb[l];
     105            level.set(w, l);
     106          }
    88107        }
    89 
     108      }
     109       
    90110      level.set(s,n);
    91111
    92       /* Starting flow. It is everywhere 0 at the moment. */
     112
     113
     114      /* Starting flow. It is everywhere 0 at the moment. */     
    93115      for(OutEdgeIt e=G.template first<OutEdgeIt>(s); e.valid(); ++e)
    94116        {
    95           if ( capacity.get(e) > 0 ) {
    96             NodeIt w=G.head(e);
     117          if ( capacity.get(e) == 0 ) continue;
     118          NodeIt w=G.head(e);
     119          if ( level.get(w) < n ) {       
     120            if ( excess.get(w) == 0 && w!=t ) stack[level.get(w)].push(w);
    97121            flow.set(e, capacity.get(e));
    98             stack[level.get(w)].push(w);
    99122            excess.set(w, excess.get(w)+capacity.get(e));
    100123          }
    101124        }
    102 
     125     
    103126      /*
    104127         End of preprocessing
     
    106129
    107130
    108 
    109131      /*
    110         Push/relabel on the highest level active Nodes.
    111       */
    112        
    113       /*While there exists an active Node.*/
     132        Push/relabel on the highest level active nodes.
     133      */
     134      /*While there exists an active node.*/
    114135      while (b) {
    115 
    116         /*We decrease the bound if there is no active node of level b.*/
    117         if (stack[b].empty()) {
     136        if ( stack[b].empty() ) {
    118137          --b;
    119         } else {
    120 
    121           NodeIt w=stack[b].top();    //w is the highest label active Node.
    122           stack[b].pop();                    //We delete w from the stack.
    123        
    124           int newlevel=2*n-2;                //In newlevel we maintain the next level of w.
    125        
     138          continue;
     139        }
     140       
     141        NodeIt w=stack[b].top();        //w is a highest label active node.
     142        stack[b].pop();           
     143        int lev=level.get(w);
     144        int exc=excess.get(w);
     145        int newlevel=2*n-2;      //In newlevel we bound the next level of w.
     146       
     147        //  if ( level.get(w) < n ) { //Nem tudom ez mukodik-e
    126148          for(OutEdgeIt e=G.template first<OutEdgeIt>(w); e.valid(); ++e) {
    127             NodeIt v=G.head(e);
    128             /*e is the Edge wv.*/
    129 
    130             if (flow.get(e)<capacity.get(e)) {             
    131               /*e is an Edge of the residual graph */
    132 
    133               if(level.get(w)==level.get(v)+1) {     
    134                 /*Push is allowed now*/
    135 
    136                 if (capacity.get(e)-flow.get(e) > excess.get(w)) {       
    137                   /*A nonsaturating push.*/
    138                  
    139                   if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
    140                   /*v becomes active.*/
    141                  
    142                   flow.set(e, flow.get(e)+excess.get(w));
    143                   excess.set(v, excess.get(v)+excess.get(w));
    144                   excess.set(w,0);
    145                   //std::cout << w << " " << v <<" elore elen nonsat pump "  << std::endl;
    146                   break;
    147                 } else {
    148                   /*A saturating push.*/
    149 
    150                   if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
    151                   /*v becomes active.*/
    152 
    153                   excess.set(v, excess.get(v)+capacity.get(e)-flow.get(e));
    154                   excess.set(w, excess.get(w)-capacity.get(e)+flow.get(e));
    155                   flow.set(e, capacity.get(e));
    156                   //std::cout << w <<" " << v <<" elore elen sat pump "   << std::endl;
    157                   if (excess.get(w)==0) break;
    158                   /*If w is not active any more, then we go on to the next Node.*/
    159                  
    160                 } // if (capacity.get(e)-flow.get(e) > excess.get(w))
    161               } // if (level.get(w)==level.get(v)+1)
    162            
    163               else {newlevel = newlevel < level.get(v) ? newlevel : level.get(v);}
    164            
    165             } //if (flow.get(e)<capacity.get(e))
    166          
    167           } //for(OutEdgeIt e=G.first_OutEdge(w); e.valid(); ++e)
     149           
     150            if ( flow.get(e) == capacity.get(e) ) continue;
     151            NodeIt v=G.head(e);           
     152            //e=wv         
     153           
     154            if( lev > level.get(v) ) {     
     155              /*Push is allowed now*/
     156             
     157              if ( excess.get(v)==0 && v != s && v !=t )
     158                stack[level.get(v)].push(v);
     159              /*v becomes active.*/
     160             
     161              int cap=capacity.get(e);
     162              int flo=flow.get(e);
     163              int remcap=cap-flo;
     164             
     165              if ( remcap >= exc ) {       
     166                /*A nonsaturating push.*/
     167               
     168                flow.set(e, flo+exc);
     169                excess.set(v, excess.get(v)+exc);
     170                exc=0;
     171                break;
     172               
     173              } else {
     174                /*A saturating push.*/
     175               
     176                flow.set(e, cap );
     177                excess.set(v, excess.get(v)+remcap);
     178                exc-=remcap;
     179              }
     180            } else if ( newlevel > level.get(v) ) newlevel = level.get(v);
     181           
     182          } //for out edges wv
     183       
     184       
     185        if ( exc > 0 ) {       
     186          for( InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) {
     187           
     188            if( flow.get(e) == 0 ) continue;
     189            NodeIt v=G.tail(e); 
     190            //e=vw
     191           
     192            if( lev > level.get(v) ) { 
     193              /*Push is allowed now*/
     194             
     195              if ( excess.get(v)==0 && v != s && v !=t)
     196                stack[level.get(v)].push(v);
     197              /*v becomes active.*/
     198             
     199              int flo=flow.get(e);
     200             
     201              if ( flo >= exc ) {
     202                /*A nonsaturating push.*/
     203               
     204                flow.set(e, flo-exc);
     205                excess.set(v, excess.get(v)+exc);
     206                exc=0;
     207                break;
     208              } else {                                               
     209                /*A saturating push.*/
     210               
     211                excess.set(v, excess.get(v)+flo);
     212                exc-=flo;
     213                flow.set(e,0);
     214              } 
     215            } else if ( newlevel > level.get(v) ) newlevel = level.get(v);
     216           
     217          } //for in edges vw
    168218         
    169 
    170 
    171           for(InEdgeIt e=G.template first<InEdgeIt>(w); e.valid(); ++e) {
    172             NodeIt v=G.tail(e);
    173             /*e is the Edge vw.*/
    174 
    175             if (excess.get(w)==0) break;
    176             /*It may happen, that w became inactive in the first 'for' cycle.*/         
    177  
    178             if(flow.get(e)>0) {             
    179               /*e is an Edge of the residual graph */
    180 
    181               if(level.get(w)==level.get(v)+1) { 
    182                 /*Push is allowed now*/
    183                
    184                 if (flow.get(e) > excess.get(w)) {
    185                   /*A nonsaturating push.*/
    186                  
    187                   if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
    188                   /*v becomes active.*/
    189 
    190                   flow.set(e, flow.get(e)-excess.get(w));
    191                   excess.set(v, excess.get(v)+excess.get(w));
    192                   excess.set(w,0);
    193                   //std::cout << v << " " << w << " vissza elen nonsat pump "     << std::endl;
    194                   break;
    195                 } else {                                               
    196                   /*A saturating push.*/
    197                  
    198                   if (excess.get(v)==0 && v != s) stack[level.get(v)].push(v);
    199                   /*v becomes active.*/
    200                  
    201                   flow.set(e,0);
    202                   excess.set(v, excess.get(v)+flow.get(e));
    203                   excess.set(w, excess.get(w)-flow.get(e));
    204                   //std::cout << v <<" " << w << " vissza elen sat pump "     << std::endl;
    205                   if (excess.get(w)==0) { break;}
    206                 } //if (flow.get(e) > excess.get(v))
    207               } //if(level.get(w)==level.get(v)+1)
    208              
    209               else {newlevel = newlevel < level.get(v) ? newlevel : level.get(v);}
    210               //std::cout << "Leveldecrease of Node " << w << " to " << newlevel << std::endl;
    211 
    212             } //if (flow.get(e)>0)
    213 
    214           } //for in-Edge
    215 
    216 
    217 
    218 
    219           /*
    220             Relabel
    221           */
    222           if (excess.get(w)>0) {
    223             /*Now newlevel <= n*/
    224 
    225             int l=level.get(w);         //l is the old level of w.
    226             --numb[l];
    227            
    228             if (newlevel == n) {
    229               level.set(w,n);
    230              
    231             } else {
    232              
    233               if (numb[l]) {
    234                 /*If the level of w remains nonempty.*/
    235                
    236                 level.set(w,++newlevel);
    237                 ++numb[newlevel];
    238                 stack[newlevel].push(w);
    239                 b=newlevel;
    240               } else {
    241                 /*If the level of w gets empty.*/
    242              
    243                 for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid() ; ++v) {
    244                   if (level.get(v) >= l ) {
    245                     level.set(v,n); 
    246                   }
    247                 }
    248                
    249                 for (int i=l+1 ; i!=n ; ++i) numb[i]=0;
    250               } //if (numb[l])
    251        
    252             } // if (newlevel = n)
    253          
    254           } // if (excess.get(w)>0)
    255 
    256 
    257         } //else
    258        
     219        } // if w still has excess after the out edge for cycle
     220       
     221        excess.set(w, exc);
     222       
     223       
     224        /*
     225          Relabel
     226        */
     227         
     228        if ( exc > 0 ) {
     229          //now 'lev' is the old level of w
     230          level.set(w,++newlevel);
     231          --numb[lev];
     232           
     233          if ( !numb[lev] && lev < A*n ) {  //If the level of w gets empty.
     234             
     235            for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid() ; ++v) {
     236              if (level.get(v) > lev ) level.set(v,n); 
     237            }
     238            for (int i=lev+1 ; i!=n ; ++i) numb[i]=0;
     239            if ( newlevel < n ) newlevel=n;
     240          } else if ( newlevel < n ) {
     241            ++numb[newlevel];
     242            stack[newlevel].push(w);
     243            b=newlevel;
     244          }
     245        }
     246
     247
     248
    259249      } //while(b)
    260250
     
    263253     
    264254
    265 
    266       /*
    267         We find an empty level, e. The Nodes above this level give
    268         a minimum cut.
    269       */
    270      
    271       int e=1;
    272      
    273       while(e) {
    274         if(numb[e]) ++e;
     255      /*
     256         We count empty_level. The nodes above this level is a mincut.
     257      */
     258      while(true) {
     259        if(numb[empty_level]) ++empty_level;
    275260        else break;
    276261      }
    277       for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v) {
    278         if (level.get(v) > e) mincutvector.set(v, true);
    279       }
    280      
    281 
     262     
    282263    } // void run()
    283264
     
    296277    /*
    297278      Returns a minimum cut.
    298     */
    299    
    300     typename Graph::NodeMap<bool> mincut() {
    301       return mincutvector;
     279    */   
     280    template<typename CutMap>
     281    void mincut(CutMap& M) {
     282
     283      for (EachNodeIt v=G.template first<EachNodeIt>(); v.valid(); ++v) {
     284        if ( level.get(v) > empty_level ) M.set(v, true);
     285      }
    302286    }
    303    
     287
    304288
    305289  };
Note: See TracChangeset for help on using the changeset viewer.