#ifndef EDMONDS_KARP_HH #define EDMONDS_KARP_HH #include #include #include #include //#include namespace hugo { template class ResGraph { public: typedef typename Graph::NodeIt NodeIt; typedef typename Graph::EachNodeIt EachNodeIt; private: typedef typename Graph::SymEdgeIt OldSymEdgeIt; const Graph& G; FlowMap& flow; const CapacityMap& capacity; public: ResGraph(const Graph& _G, FlowMap& _flow, const CapacityMap& _capacity) : G(_G), flow(_flow), capacity(_capacity) { } class EdgeIt; class OutEdgeIt; friend class EdgeIt; friend class OutEdgeIt; class EdgeIt { friend class ResGraph; protected: const ResGraph* resG; OldSymEdgeIt sym; public: EdgeIt() { } //EdgeIt(const EdgeIt& e) : resG(e.resG), sym(e.sym) { } Number free() const { if (resG->G.aNode(sym)==resG->G.tail(sym)) { return (resG->capacity.get(sym)-resG->flow.get(sym)); } else { return (resG->flow.get(sym)); } } bool valid() const { return sym.valid(); } void augment(Number a) const { if (resG->G.aNode(sym)==resG->G.tail(sym)) { resG->flow.set(sym, resG->flow.get(sym)+a); //resG->flow[sym]+=a; } else { resG->flow.set(sym, resG->flow.get(sym)-a); //resG->flow[sym]-=a; } } }; class OutEdgeIt : public EdgeIt { friend class ResGraph; public: OutEdgeIt() { } //OutEdgeIt(const OutEdgeIt& e) { resG=e.resG; sym=e.sym; } private: OutEdgeIt(const ResGraph& _resG, NodeIt v) { resG=&_resG; sym=resG->G.template first(v); while( sym.valid() && !(free()>0) ) { ++sym; } } public: OutEdgeIt& operator++() { ++sym; while( sym.valid() && !(free()>0) ) { ++sym; } return *this; } }; void getFirst(OutEdgeIt& e, NodeIt v) const { e=OutEdgeIt(*this, v); } void getFirst(EachNodeIt& v) const { G.getFirst(v); } template< typename It > It first() const { It e; getFirst(e); return e; } template< typename It > It first(NodeIt v) const { It e; getFirst(e, v); return e; } NodeIt tail(EdgeIt e) const { return G.aNode(e.sym); } NodeIt head(EdgeIt e) const { return G.bNode(e.sym); } NodeIt aNode(OutEdgeIt e) const { return G.aNode(e.sym); } NodeIt bNode(OutEdgeIt e) const { return G.bNode(e.sym); } int id(NodeIt v) const { return G.id(v); } template class NodeMap { typename Graph::NodeMap node_map; public: NodeMap(const ResGraph& _G) : node_map(_G.G) { } NodeMap(const ResGraph& _G, S a) : node_map(_G.G, a) { } void set(NodeIt nit, S a) { node_map.set(nit, a); } S get(NodeIt nit) const { return node_map.get(nit); } S& operator[](NodeIt nit) { return node_map[nit]; } const S& operator[](NodeIt nit) const { return node_map[nit]; } }; }; template class ResGraph2 { public: typedef typename Graph::NodeIt NodeIt; typedef typename Graph::EachNodeIt EachNodeIt; private: //typedef typename Graph::SymEdgeIt OldSymEdgeIt; typedef typename Graph::OutEdgeIt OldOutEdgeIt; typedef typename Graph::InEdgeIt OldInEdgeIt; const Graph& G; FlowMap& flow; const CapacityMap& capacity; public: ResGraph2(const Graph& _G, FlowMap& _flow, const CapacityMap& _capacity) : G(_G), flow(_flow), capacity(_capacity) { } class EdgeIt; class OutEdgeIt; friend class EdgeIt; friend class OutEdgeIt; class EdgeIt { friend class ResGraph2; protected: const ResGraph2* resG; //OldSymEdgeIt sym; OldOutEdgeIt out; OldInEdgeIt in; bool out_or_in; //true, iff out public: EdgeIt() : out_or_in(true) { } Number free() const { if (out_or_in) { return (resG->capacity.get(out)-resG->flow.get(out)); } else { return (resG->flow.get(in)); } } bool valid() const { return out_or_in && out.valid() || in.valid(); } void augment(Number a) const { if (out_or_in) { resG->flow.set(out, resG->flow.get(out)+a); } else { resG->flow.set(in, resG->flow.get(in)-a); } } }; class OutEdgeIt : public EdgeIt { friend class ResGraph2; public: OutEdgeIt() { } private: OutEdgeIt(const ResGraph2& _resG, NodeIt v) { resG=&_resG; out=resG->G.template first(v); while( out.valid() && !(free()>0) ) { ++out; } if (!out.valid()) { out_or_in=0; in=resG->G.template first(v); while( in.valid() && !(free()>0) ) { ++in; } } } public: OutEdgeIt& operator++() { if (out_or_in) { NodeIt v=resG->G.aNode(out); ++out; while( out.valid() && !(free()>0) ) { ++out; } if (!out.valid()) { out_or_in=0; in=resG->G.template first(v); while( in.valid() && !(free()>0) ) { ++in; } } } else { ++in; while( in.valid() && !(free()>0) ) { ++in; } } return *this; } }; void getFirst(OutEdgeIt& e, NodeIt v) const { e=OutEdgeIt(*this, v); } void getFirst(EachNodeIt& v) const { G.getFirst(v); } template< typename It > It first() const { It e; getFirst(e); return e; } template< typename It > It first(NodeIt v) const { It e; getFirst(e, v); return e; } NodeIt tail(EdgeIt e) const { return ((e.out_or_in) ? G.aNode(e.out) : G.aNode(e.in)); } NodeIt head(EdgeIt e) const { return ((e.out_or_in) ? G.bNode(e.out) : G.bNode(e.in)); } NodeIt aNode(OutEdgeIt e) const { return ((e.out_or_in) ? G.aNode(e.out) : G.aNode(e.in)); } NodeIt bNode(OutEdgeIt e) const { return ((e.out_or_in) ? G.bNode(e.out) : G.bNode(e.in)); } int id(NodeIt v) const { return G.id(v); } template class NodeMap { typename Graph::NodeMap node_map; public: NodeMap(const ResGraph2& _G) : node_map(_G.G) { } NodeMap(const ResGraph2& _G, S a) : node_map(_G.G, a) { } void set(NodeIt nit, S a) { node_map.set(nit, a); } S get(NodeIt nit) const { return node_map.get(nit); } }; }; template class MaxFlow { public: typedef typename Graph::NodeIt NodeIt; typedef typename Graph::EdgeIt EdgeIt; typedef typename Graph::EachEdgeIt EachEdgeIt; typedef typename Graph::OutEdgeIt OutEdgeIt; typedef typename Graph::InEdgeIt InEdgeIt; private: const Graph* G; NodeIt s; NodeIt t; FlowMap* flow; const CapacityMap* capacity; typedef ResGraphWrapper AugGraph; typedef typename AugGraph::OutEdgeIt AugOutEdgeIt; typedef typename AugGraph::EdgeIt AugEdgeIt; //AugGraph res_graph; //typedef typename AugGraph::NodeMap ReachedMap; //typename AugGraph::NodeMap pred; //typename AugGraph::NodeMap free; public: MaxFlow(const Graph& _G, NodeIt _s, NodeIt _t, FlowMap& _flow, const CapacityMap& _capacity) : G(&_G), s(_s), t(_t), flow(&_flow), capacity(&_capacity) //, //res_graph(G, flow, capacity), pred(res_graph), free(res_graph) { } bool augmentOnShortestPath() { AugGraph res_graph(*G, *flow, *capacity); bool _augment=false; typedef typename AugGraph::NodeMap ReachedMap; BfsIterator5< AugGraph, /*AugOutEdgeIt,*/ ReachedMap > res_bfs(res_graph); res_bfs.pushAndSetReached(s); typename AugGraph::NodeMap pred(res_graph); //filled up with invalid iterators //pred.set(s, AugEdgeIt()); typename AugGraph::NodeMap free(res_graph); //searching for augmenting path while ( !res_bfs.finished() ) { AugOutEdgeIt e=/*AugOutEdgeIt*/(res_bfs); if (res_graph.valid(e) && res_bfs.isBNodeNewlyReached()) { NodeIt v=res_graph.tail(e); NodeIt w=res_graph.head(e); pred.set(w, e); if (res_graph.valid(pred.get(v))) { free.set(w, std::min(free.get(v), res_graph.free(e))); } else { free.set(w, res_graph.free(e)); } if (res_graph.head(e)==t) { _augment=true; break; } } ++res_bfs; } //end of searching augmenting path if (_augment) { NodeIt n=t; Number augment_value=free.get(t); while (res_graph.valid(pred.get(n))) { AugEdgeIt e=pred.get(n); res_graph.augment(e, augment_value); //e.augment(augment_value); n=res_graph.tail(e); } } return _augment; } template bool augmentOnBlockingFlow() { bool _augment=false; AugGraph res_graph(*G, *flow, *capacity); typedef typename AugGraph::NodeMap ReachedMap; BfsIterator4< AugGraph, AugOutEdgeIt, ReachedMap > bfs(res_graph); bfs.pushAndSetReached(s); typename AugGraph::NodeMap dist(res_graph); //filled up with 0's while ( !bfs.finished() ) { AugOutEdgeIt e=/*AugOutEdgeIt*/(bfs); if (res_graph.valid(e) && bfs.isBNodeNewlyReached()) { dist.set(res_graph.head(e), dist.get(res_graph.tail(e))+1); } ++bfs; } //computing distances from s in the residual graph MutableGraph F; typename AugGraph::NodeMap res_graph_to_F(res_graph); for(typename AugGraph::EachNodeIt n=res_graph.template first(); res_graph.valid(n); res_graph.next(n)) { res_graph_to_F.set(n, F.addNode()); } typename MutableGraph::NodeIt sF=res_graph_to_F.get(s); typename MutableGraph::NodeIt tF=res_graph_to_F.get(t); typename MutableGraph::EdgeMap original_edge(F); typename MutableGraph::EdgeMap residual_capacity(F); //Making F to the graph containing the edges of the residual graph //which are in some shortest paths for(typename AugGraph::EachEdgeIt e=res_graph.template first(); res_graph.valid(e); res_graph.next(e)) { if (dist.get(res_graph.head(e))==dist.get(res_graph.tail(e))+1) { typename MutableGraph::EdgeIt f=F.addEdge(res_graph_to_F.get(res_graph.tail(e)), res_graph_to_F.get(res_graph.head(e))); original_edge.update(); original_edge.set(f, e); residual_capacity.update(); residual_capacity.set(f, res_graph.free(e)); } } bool __augment=true; while (__augment) { __augment=false; //computing blocking flow with dfs typedef typename MutableGraph::NodeMap BlockingReachedMap; DfsIterator4< MutableGraph, typename MutableGraph::OutEdgeIt, BlockingReachedMap > dfs(F); typename MutableGraph::NodeMap pred(F); //invalid iterators typename MutableGraph::NodeMap free(F); dfs.pushAndSetReached(sF); while (!dfs.finished()) { ++dfs; if (F.valid(typename MutableGraph::OutEdgeIt(dfs))) { if (dfs.isBNodeNewlyReached()) { // std::cout << "OutEdgeIt: " << dfs; // std::cout << " aNode: " << F.aNode(dfs); // std::cout << " bNode: " << F.bNode(dfs) << " "; typename MutableGraph::NodeIt v=F.aNode(dfs); typename MutableGraph::NodeIt w=F.bNode(dfs); pred.set(w, dfs); if (F.valid(pred.get(v))) { free.set(w, std::min(free.get(v), residual_capacity.get(dfs))); } else { free.set(w, residual_capacity.get(dfs)); } if (w==tF) { //std::cout << "AUGMENTATION"< EAugGraph; typedef FilterGraphWrapper< ErasingResGraphWrapper > EAugGraph; typedef typename EAugGraph::OutEdgeIt EAugOutEdgeIt; typedef typename EAugGraph::EdgeIt EAugEdgeIt; EAugGraph res_graph(*G, *flow, *capacity); //std::cout << "meg jo1" << std::endl; //typedef typename EAugGraph::NodeMap ReachedMap; BfsIterator4< ErasingResGraphWrapper, ErasingResGraphWrapper::OutEdgeIt, ErasingResGraphWrapper::NodeMap > bfs(res_graph); //std::cout << "meg jo2" << std::endl; bfs.pushAndSetReached(s); //std::cout << "meg jo2.5" << std::endl; //typename EAugGraph::NodeMap dist(res_graph); //filled up with 0's typename ErasingResGraphWrapper:: NodeMap& dist=res_graph.dist; //std::cout << "meg jo2.6" << std::endl; while ( !bfs.finished() ) { ErasingResGraphWrapper::OutEdgeIt e=bfs; // EAugOutEdgeIt e=/*AugOutEdgeIt*/(bfs); //if (res_graph.valid(e)) { // std::cout<<"a:"<(); res_graph.valid(n); res_graph.next(n)) { // std::cout << "dist: " << dist.get(n) << std::endl; // } bool __augment=true; while (__augment) { // std::cout << "new iteration"<< std::endl; __augment=false; //computing blocking flow with dfs typedef typename EAugGraph::NodeMap BlockingReachedMap; DfsIterator4< EAugGraph, EAugOutEdgeIt, BlockingReachedMap > dfs(res_graph); typename EAugGraph::NodeMap pred(res_graph); //invalid iterators typename EAugGraph::NodeMap free(res_graph); dfs.pushAndSetReached(s); while (!dfs.finished()) { ++dfs; if (res_graph.valid(EAugOutEdgeIt(dfs))) { if (dfs.isBNodeNewlyReached()) { // std::cout << "OutEdgeIt: " << dfs; // std::cout << " aNode: " << res_graph.aNode(dfs); // std::cout << " res cap: " << EAugOutEdgeIt(dfs).free(); // std::cout << " bNode: " << res_graph.bNode(dfs) << " "; typename EAugGraph::NodeIt v=res_graph.aNode(dfs); typename EAugGraph::NodeIt w=res_graph.bNode(dfs); pred.set(w, EAugOutEdgeIt(dfs)); //std::cout << EAugOutEdgeIt(dfs).free() << std::endl; if (res_graph.valid(pred.get(v))) { free.set(w, std::min(free.get(v), res_graph.free(/*EAugOutEdgeIt*/(dfs)))); } else { free.set(w, res_graph.free(/*EAugOutEdgeIt*/(dfs))); } if (w==t) { // std::cout << "t is reached, AUGMENTATION"<> "; res_graph.erase(dfs); } } } if (__augment) { typename EAugGraph::NodeIt n=t; Number augment_value=free.get(t); // std::cout << "av:" << augment_value << std::endl; while (res_graph.valid(pred.get(n))) { EAugEdgeIt e=pred.get(n); res_graph.augment(e, augment_value); //e.augment(augment_value); n=res_graph.tail(e); if (res_graph.free(e)==0) res_graph.erase(e); } } } return _augment; } void run() { //int num_of_augmentations=0; while (augmentOnShortestPath()) { //while (augmentOnBlockingFlow()) { //std::cout << ++num_of_augmentations << " "; //std::cout< void run() { //int num_of_augmentations=0; //while (augmentOnShortestPath()) { while (augmentOnBlockingFlow()) { //std::cout << ++num_of_augmentations << " "; //std::cout<getFirst(e, s); G->valid(e); G->next(e)) { a+=flow->get(e); } return a; } }; // template // class MaxFlow2 { // public: // typedef typename Graph::NodeIt NodeIt; // typedef typename Graph::EdgeIt EdgeIt; // typedef typename Graph::EachEdgeIt EachEdgeIt; // typedef typename Graph::OutEdgeIt OutEdgeIt; // typedef typename Graph::InEdgeIt InEdgeIt; // private: // const Graph& G; // std::list& S; // std::list& T; // FlowMap& flow; // const CapacityMap& capacity; // typedef ResGraphWrapper AugGraph; // typedef typename AugGraph::OutEdgeIt AugOutEdgeIt; // typedef typename AugGraph::EdgeIt AugEdgeIt; // typename Graph::NodeMap SMap; // typename Graph::NodeMap TMap; // public: // MaxFlow2(const Graph& _G, std::list& _S, std::list& _T, FlowMap& _flow, const CapacityMap& _capacity) : G(_G), S(_S), T(_T), flow(_flow), capacity(_capacity), SMap(_G), TMap(_G) { // for(typename std::list::const_iterator i=S.begin(); // i!=S.end(); ++i) { // SMap.set(*i, true); // } // for (typename std::list::const_iterator i=T.begin(); // i!=T.end(); ++i) { // TMap.set(*i, true); // } // } // bool augment() { // AugGraph res_graph(G, flow, capacity); // bool _augment=false; // NodeIt reached_t_node; // typedef typename AugGraph::NodeMap ReachedMap; // BfsIterator4< AugGraph, AugOutEdgeIt, ReachedMap > res_bfs(res_graph); // for(typename std::list::const_iterator i=S.begin(); // i!=S.end(); ++i) { // res_bfs.pushAndSetReached(*i); // } // //res_bfs.pushAndSetReached(s); // typename AugGraph::NodeMap pred(res_graph); // //filled up with invalid iterators // typename AugGraph::NodeMap free(res_graph); // //searching for augmenting path // while ( !res_bfs.finished() ) { // AugOutEdgeIt e=/*AugOutEdgeIt*/(res_bfs); // if (e.valid() && res_bfs.isBNodeNewlyReached()) { // NodeIt v=res_graph.tail(e); // NodeIt w=res_graph.head(e); // pred.set(w, e); // if (pred.get(v).valid()) { // free.set(w, std::min(free.get(v), e.free())); // } else { // free.set(w, e.free()); // } // if (TMap.get(res_graph.head(e))) { // _augment=true; // reached_t_node=res_graph.head(e); // break; // } // } // ++res_bfs; // } //end of searching augmenting path // if (_augment) { // NodeIt n=reached_t_node; // Number augment_value=free.get(reached_t_node); // while (pred.get(n).valid()) { // AugEdgeIt e=pred.get(n); // e.augment(augment_value); // n=res_graph.tail(e); // } // } // return _augment; // } // void run() { // while (augment()) { } // } // Number flowValue() { // Number a=0; // for(typename std::list::const_iterator i=S.begin(); // i!=S.end(); ++i) { // for(OutEdgeIt e=G.template first(*i); e.valid(); ++e) { // a+=flow.get(e); // } // for(InEdgeIt e=G.template first(*i); e.valid(); ++e) { // a-=flow.get(e); // } // } // return a; // } // }; } // namespace hugo #endif //EDMONDS_KARP_HH