src/work/edmonds_karp.hh
changeset 133 0631992fe7a1
parent 127 dcace15b1874
child 135 1e5060d1fa1d
     1.1 --- a/src/work/edmonds_karp.hh	Thu Feb 26 16:07:40 2004 +0000
     1.2 +++ b/src/work/edmonds_karp.hh	Fri Feb 27 12:39:15 2004 +0000
     1.3 @@ -6,7 +6,7 @@
     1.4  #include <iterator>
     1.5  
     1.6  #include <bfs_iterator.hh>
     1.7 -#include <time_measure.h>
     1.8 +//#include <time_measure.h>
     1.9  
    1.10  namespace hugo {
    1.11  
    1.12 @@ -281,6 +281,7 @@
    1.13        bool out_or_in; //1, iff out
    1.14      public:
    1.15        EdgeIt() : out_or_in(1) { } 
    1.16 +      //EdgeIt(const EdgeIt& e) : G(e.G), flow(e.flow), capacity(e.capacity), out(e.out), in(e.in), out_or_in(e.out_or_in) { }
    1.17        Number free() const { 
    1.18  	if (out_or_in) { 
    1.19  	  return (/*resG->*/capacity->get(out)-/*resG->*/flow->get(out)); 
    1.20 @@ -297,6 +298,23 @@
    1.21  	  /*resG->*/flow->set(in, /*resG->*/flow->get(in)-a);
    1.22  	}
    1.23        }
    1.24 +      void print() { 
    1.25 +	if (out_or_in) {
    1.26 +	  std::cout << "out "; 
    1.27 +	  if (out.valid()) 
    1.28 +	    std::cout << G->id(G->tail(out)) << "--"<< G->id(out) <<"->"<< G->id(G->head(out)); 
    1.29 +	  else 
    1.30 +	    std::cout << "invalid"; 
    1.31 +	}
    1.32 +	else {
    1.33 +	  std::cout << "in "; 
    1.34 +	  if (in.valid()) 
    1.35 +	    std::cout << G->id(G->head(in)) << "<-"<< G->id(in) <<"--"<< G->id(G->tail(in)); 
    1.36 +	  else 
    1.37 +	    std::cout << "invalid"; 
    1.38 +	}
    1.39 +	std::cout << std::endl;
    1.40 +      }
    1.41      };
    1.42  
    1.43      class OutEdgeIt : public EdgeIt {
    1.44 @@ -308,11 +326,13 @@
    1.45        	G=&_G;
    1.46  	flow=&_flow;
    1.47  	capacity=&_capacity;
    1.48 -	out=/*resG->*/G->template first<OldOutEdgeIt>(v);
    1.49 +	//out=/*resG->*/G->template first<OldOutEdgeIt>(v);
    1.50 +	G->getFirst(out, v);
    1.51  	while( out.valid() && !(free()>0) ) { ++out; }
    1.52  	if (!out.valid()) {
    1.53  	  out_or_in=0;
    1.54 -	  in=/*resG->*/G->template first<OldInEdgeIt>(v);
    1.55 +	  //in=/*resG->*/G->template first<OldInEdgeIt>(v);
    1.56 +	  G->getFirst(in, v);
    1.57  	  while( in.valid() && !(free()>0) ) { ++in; }
    1.58  	}
    1.59        }
    1.60 @@ -324,7 +344,7 @@
    1.61  	  while( out.valid() && !(free()>0) ) { ++out; }
    1.62  	  if (!out.valid()) {
    1.63  	    out_or_in=0;
    1.64 -	    in=/*resG->*/G->template first<OldInEdgeIt>(v);
    1.65 +	    G->getFirst(in, v); //=/*resG->*/G->template first<OldInEdgeIt>(v);
    1.66  	    while( in.valid() && !(free()>0) ) { ++in; }
    1.67  	  }
    1.68  	} else {
    1.69 @@ -335,9 +355,75 @@
    1.70        }
    1.71      };
    1.72  
    1.73 +    class EachEdgeIt : public EdgeIt {
    1.74 +      typename Graph::EachNodeIt v;
    1.75 +    public:
    1.76 +      EachEdgeIt() { }
    1.77 +      //EachEdgeIt(const EachEdgeIt& e) : EdgeIt(e), v(e.v) { }
    1.78 +      EachEdgeIt(const Graph& _G, FlowMap& _flow, const CapacityMap& _capacity) { 
    1.79 +	G=&_G;
    1.80 +	flow=&_flow;
    1.81 +	capacity=&_capacity;
    1.82 +	out_or_in=1;
    1.83 +	G->getFirst(v);
    1.84 +	if (v.valid()) G->getFirst(out, v); else out=OldOutEdgeIt();
    1.85 +	while (out.valid() && !(free()>0) ) { ++out; }
    1.86 +	while (v.valid() && !out.valid()) { 
    1.87 +	  ++v; 
    1.88 +	  if (v.valid()) G->getFirst(out, v); 
    1.89 +	  while (out.valid() && !(free()>0) ) { ++out; }
    1.90 +	}
    1.91 +	if (!out.valid()) {
    1.92 +	  out_or_in=0;
    1.93 +	  G->getFirst(v);
    1.94 +	  if (v.valid()) G->getFirst(in, v); else in=OldInEdgeIt();
    1.95 +	  while (in.valid() && !(free()>0) ) { ++in; }
    1.96 +	  while (v.valid() && !in.valid()) { 
    1.97 +	    ++v; 
    1.98 +	    if (v.valid()) G->getFirst(in, v); 
    1.99 +	    while (in.valid() && !(free()>0) ) { ++in; }
   1.100 +	  }
   1.101 +	}
   1.102 +      }
   1.103 +      EachEdgeIt& operator++() { 
   1.104 +	if (out_or_in) {
   1.105 +	  ++out;
   1.106 +	  while (out.valid() && !(free()>0) ) { ++out; }
   1.107 +	  while (v.valid() && !out.valid()) { 
   1.108 +	    ++v; 
   1.109 +	    if (v.valid()) G->getFirst(out, v); 
   1.110 +	    while (out.valid() && !(free()>0) ) { ++out; }
   1.111 +	  }
   1.112 +	  if (!out.valid()) {
   1.113 +	    out_or_in=0;
   1.114 +	    G->getFirst(v);
   1.115 +	    if (v.valid()) G->getFirst(in, v); else in=OldInEdgeIt();
   1.116 +	    while (in.valid() && !(free()>0) ) { ++in; }
   1.117 +	    while (v.valid() && !in.valid()) { 
   1.118 +	      ++v; 
   1.119 +	      if (v.valid()) G->getFirst(in, v); 
   1.120 +	      while (in.valid() && !(free()>0) ) { ++in; }
   1.121 +	    }  
   1.122 +	  }
   1.123 +	} else {
   1.124 +	  ++in;
   1.125 +	  while (in.valid() && !(free()>0) ) { ++in; }
   1.126 +	  while (v.valid() && !in.valid()) { 
   1.127 +	    ++v; 
   1.128 +	    if (v.valid()) G->getFirst(in, v); 
   1.129 +	    while (in.valid() && !(free()>0) ) { ++in; }
   1.130 +	  }
   1.131 +	}
   1.132 +	return *this;
   1.133 +      }
   1.134 +    };
   1.135 +
   1.136      void getFirst(OutEdgeIt& e, NodeIt v) const { 
   1.137        e=OutEdgeIt(G, v, flow, capacity); 
   1.138      }
   1.139 +    void getFirst(EachEdgeIt& e) const { 
   1.140 +      e=EachEdgeIt(G, flow, capacity); 
   1.141 +    }
   1.142      void getFirst(EachNodeIt& v) const { G.getFirst(v); }
   1.143      
   1.144      template< typename It >
   1.145 @@ -401,13 +487,13 @@
   1.146      //AugGraph res_graph;    
   1.147      //typedef typename AugGraph::NodeMap<bool> ReachedMap;
   1.148      //typename AugGraph::NodeMap<AugEdgeIt> pred; 
   1.149 -    //typename AugGraph::NodeMap<int> free;
   1.150 +    //typename AugGraph::NodeMap<Number> free;
   1.151    public:
   1.152      MaxFlow(const Graph& _G, NodeIt _s, NodeIt _t, FlowMap& _flow, const CapacityMap& _capacity) : 
   1.153        G(_G), s(_s), t(_t), flow(_flow), capacity(_capacity) //,  
   1.154        //res_graph(G, flow, capacity), pred(res_graph), free(res_graph) 
   1.155        { }
   1.156 -    bool augment() {
   1.157 +    bool augmentOnShortestPath() {
   1.158        AugGraph res_graph(G, flow, capacity);
   1.159        bool _augment=false;
   1.160        
   1.161 @@ -419,7 +505,7 @@
   1.162        //filled up with invalid iterators
   1.163        //pred.set(s, AugEdgeIt());
   1.164        
   1.165 -      typename AugGraph::NodeMap<int> free(res_graph);
   1.166 +      typename AugGraph::NodeMap<Number> free(res_graph);
   1.167  	
   1.168        //searching for augmenting path
   1.169        while ( !res_bfs.finished() ) { 
   1.170 @@ -451,48 +537,130 @@
   1.171  
   1.172        return _augment;
   1.173      }
   1.174 -    bool augmentWithBlockingFlow() {
   1.175 -      BfsIterator4< Graph, OutEdgeIt, typename Graph::NodeMap<bool> > bfs(G);
   1.176 +    template<typename MutableGraph> bool augmentOnBlockingFlow() {
   1.177 +      bool _augment=false;
   1.178 +
   1.179 +      AugGraph res_graph(G, flow, capacity);
   1.180 +
   1.181 +      typedef typename AugGraph::NodeMap<bool> ReachedMap;
   1.182 +      BfsIterator4< AugGraph, AugOutEdgeIt, ReachedMap > bfs(res_graph);
   1.183 +
   1.184        bfs.pushAndSetReached(s);
   1.185 -      typename Graph::NodeMap<int> dist(G); //filled up with 0's
   1.186 +      typename AugGraph::NodeMap<int> dist(res_graph); //filled up with 0's
   1.187        while ( !bfs.finished() ) { 
   1.188 -	OutEdgeIt e=OutEdgeIt(bfs);
   1.189 +	AugOutEdgeIt e=AugOutEdgeIt(bfs);
   1.190  	if (e.valid() && bfs.isBNodeNewlyReached()) {
   1.191 -	  dist.set(G.head(e), dist.get(G.tail(e))+1);
   1.192 -	  //NodeIt v=res_graph.tail(e);
   1.193 -	  //NodeIt w=res_graph.head(e);
   1.194 -	  //pred.set(w, e);
   1.195 -	  //if (pred.get(v).valid()) {
   1.196 -	  //  free.set(w, std::min(free.get(v), e.free()));
   1.197 -	  //} else {
   1.198 -	  //  free.set(w, e.free()); 
   1.199 -	  //}
   1.200 -	  //if (res_graph.head(e)==t) { _augment=true; break; }
   1.201 +	  dist.set(res_graph.head(e), dist.get(res_graph.tail(e))+1);
   1.202  	}
   1.203  	
   1.204  	++bfs;
   1.205 -      } //end of searching augmenting path
   1.206 +      } //computing distances from s in the residual graph
   1.207  
   1.208 -      double pre_time_copy=currTime();
   1.209 -      typedef Graph MutableGraph;
   1.210        MutableGraph F;
   1.211 -      typename Graph::NodeMap<NodeIt> G_to_F(G);
   1.212 -      for(typename Graph::EachNodeIt n=G.template first<typename Graph::EachNodeIt>(); n.valid(); ++n) {
   1.213 -	G_to_F.set(n, F.addNode());
   1.214 +      typename AugGraph::NodeMap<NodeIt> res_graph_to_F(res_graph);
   1.215 +      for(typename AugGraph::EachNodeIt n=res_graph.template first<typename AugGraph::EachNodeIt>(); n.valid(); ++n) {
   1.216 +	res_graph_to_F.set(n, F.addNode());
   1.217        }
   1.218 -      for(typename Graph::EachEdgeIt e=G.template first<typename Graph::EachEdgeIt>(); e.valid(); ++e) {
   1.219 -	if (dist.get(G.head(e))==dist.get(G.tail(e))+1) {
   1.220 -	  F.addEdge(G_to_F.get(G.tail(e)), G_to_F.get(G.head(e)));
   1.221 +      
   1.222 +      typename MutableGraph::NodeIt sF=res_graph_to_F.get(s);
   1.223 +      typename MutableGraph::NodeIt tF=res_graph_to_F.get(t);
   1.224 +
   1.225 +      typename MutableGraph::EdgeMap<AugEdgeIt> original_edge(F);
   1.226 +      typename MutableGraph::EdgeMap<Number> free_on_edge(F);
   1.227 +
   1.228 +      //Making F to the graph containing the edges of the residual graph 
   1.229 +      //which are in some shortest paths
   1.230 +      for(typename AugGraph::EachEdgeIt e=res_graph.template first<typename AugGraph::EachEdgeIt>(); e.valid(); ++e) {
   1.231 +	if (dist.get(res_graph.head(e))==dist.get(res_graph.tail(e))+1) {
   1.232 +	  typename MutableGraph::EdgeIt f=F.addEdge(res_graph_to_F.get(res_graph.tail(e)), res_graph_to_F.get(res_graph.head(e)));
   1.233 +	  original_edge.update();
   1.234 +	  original_edge.set(f, e);
   1.235 +	  free_on_edge.update();
   1.236 +	  free_on_edge.set(f, e.free());
   1.237 +	} 
   1.238 +      }
   1.239 +
   1.240 +      bool __augment=true;
   1.241 +
   1.242 +      while (__augment) {
   1.243 +	__augment=false;
   1.244 +	//computing blocking flow with dfs
   1.245 +	typedef typename MutableGraph::NodeMap<bool> BlockingReachedMap;
   1.246 +	DfsIterator4< MutableGraph, typename MutableGraph::OutEdgeIt, BlockingReachedMap > dfs(F);
   1.247 +	typename MutableGraph::NodeMap<EdgeIt> pred(F); //invalid iterators
   1.248 +	typename MutableGraph::NodeMap<Number> free(F);
   1.249 +
   1.250 +	dfs.pushAndSetReached(sF);      
   1.251 +	while (!dfs.finished()) {
   1.252 +	  ++dfs;
   1.253 +	  if (typename MutableGraph::OutEdgeIt(dfs).valid()) {
   1.254 +	    //std::cout << "OutEdgeIt: " << dfs; 
   1.255 +	    //std::cout << " aNode: " << F.aNode(dfs); 
   1.256 +	    //std::cout << " bNode: " << F.bNode(dfs) << " ";
   1.257 +	  
   1.258 +	    typename MutableGraph::NodeIt v=F.aNode(dfs);
   1.259 +	    typename MutableGraph::NodeIt w=F.bNode(dfs);
   1.260 +	    pred.set(w, dfs);
   1.261 +	    if (pred.get(v).valid()) {
   1.262 +	      free.set(w, std::min(free.get(v), free_on_edge.get(dfs)));
   1.263 +	    } else {
   1.264 +	      free.set(w, free_on_edge.get(dfs)/*original_edge.get(dfs).free()*/); 
   1.265 +	    }
   1.266 +	    if (w==tF) { 
   1.267 +	      //std::cout << "AUGMENTATION"<<std::endl;
   1.268 +	      __augment=true; 
   1.269 +	      _augment=true;
   1.270 +	      break; 
   1.271 +	    }
   1.272 +	  } else { 
   1.273 +	    //std::cout << "OutEdgeIt: " << "invalid"; 
   1.274 +	    //std::cout << " aNode: " << dfs.aNode(); 
   1.275 +	    //std::cout << " bNode: " << "invalid" << " ";
   1.276 +	  }
   1.277 +	  if (dfs.isBNodeNewlyReached()) { 
   1.278 +	    //std::cout << "bNodeIsNewlyReached ";
   1.279 +	  } else { 
   1.280 +	    //std::cout << "bNodeIsNotNewlyReached ";
   1.281 +	    if (typename MutableGraph::OutEdgeIt(dfs).valid()) {
   1.282 +	      //std::cout << "DELETE ";
   1.283 +	      F.erase(typename MutableGraph::OutEdgeIt(dfs));
   1.284 +	    }
   1.285 +	  } 
   1.286 +	  //if (dfs.isANodeExamined()) { 
   1.287 +	    //std::cout << "aNodeIsExamined ";
   1.288 +	    //} else { 
   1.289 +	    //std::cout << "aNodeIsNotExamined ";
   1.290 +	    //} 
   1.291 +	  //std::cout<<std::endl;
   1.292  	}
   1.293 +
   1.294 +	if (__augment) {
   1.295 +	  typename MutableGraph::NodeIt n=tF;
   1.296 +	  Number augment_value=free.get(tF);
   1.297 +	  while (pred.get(n).valid()) { 
   1.298 +	    typename MutableGraph::EdgeIt e=pred.get(n);
   1.299 +	    original_edge.get(e).augment(augment_value); 
   1.300 +	    n=F.tail(e);
   1.301 +	    F.erase(e);
   1.302 +	  }
   1.303 +	}
   1.304 +
   1.305        }
   1.306 -      double post_time_copy=currTime();
   1.307 -      std::cout << "copy time: " << post_time_copy-pre_time_copy << " sec"<< std::endl; 
   1.308 -
   1.309 -      return 0;
   1.310 +            
   1.311 +      return _augment;
   1.312      }
   1.313      void run() {
   1.314        //int num_of_augmentations=0;
   1.315 -      while (augment()) { 
   1.316 +      while (augmentOnShortestPath()) { 
   1.317 +	//while (augmentOnBlockingFlow<MutableGraph>()) { 
   1.318 +	//std::cout << ++num_of_augmentations << " ";
   1.319 +	//std::cout<<std::endl;
   1.320 +      } 
   1.321 +    }
   1.322 +    template<typename MutableGraph> void run() {
   1.323 +      //int num_of_augmentations=0;
   1.324 +      //while (augmentOnShortestPath()) { 
   1.325 +	while (augmentOnBlockingFlow<MutableGraph>()) { 
   1.326  	//std::cout << ++num_of_augmentations << " ";
   1.327  	//std::cout<<std::endl;
   1.328        } 
   1.329 @@ -599,7 +767,7 @@
   1.330        typename AugGraph::NodeMap<AugEdgeIt> pred(res_graph); 
   1.331        //filled up with invalid iterators
   1.332        
   1.333 -      typename AugGraph::NodeMap<int> free(res_graph);
   1.334 +      typename AugGraph::NodeMap<Number> free(res_graph);
   1.335  	
   1.336        //searching for augmenting path
   1.337        while ( !res_bfs.finished() ) {