marci@762: // -*- C++ -*-
alpar@921: #ifndef LEMON_AUGMENTING_FLOW_H
alpar@921: #define LEMON_AUGMENTING_FLOW_H
marci@762: 
marci@762: #include <vector>
marci@762: #include <iostream>
marci@762: 
alpar@921: #include <lemon/graph_wrapper.h>
marci@944: //#include <bfs_dfs.h>
marci@944: #include <bfs_mm.h>
alpar@921: #include <lemon/invalid.h>
alpar@921: #include <lemon/maps.h>
marci@944: #include <demo/tight_edge_filter_map.h>
marci@762: 
marci@762: /// \file
marci@762: /// \brief Maximum flow algorithms.
marci@762: /// \ingroup galgs
marci@762: 
alpar@921: namespace lemon {
marci@944:   using lemon::marci::BfsIterator;
marci@944:   using lemon::marci::DfsIterator;
marci@762: 
marci@762:   /// \addtogroup galgs
marci@762:   /// @{                                                                                                                                        
marci@862:   /// Class for augmenting path flow algorithms.
marci@762: 
marci@862:   /// This class provides various algorithms for finding a flow of
marci@862:   /// maximum value in a directed graph. The \e source node, the \e
marci@862:   /// target node, the \e capacity of the edges and the \e starting \e
marci@862:   /// flow value of the edges should be passed to the algorithm through the
marci@862:   /// constructor. 
marci@862: //   /// It is possible to change these quantities using the
marci@862: //   /// functions \ref resetSource, \ref resetTarget, \ref resetCap and
marci@862: //   /// \ref resetFlow. Before any subsequent runs of any algorithm of
marci@862: //   /// the class \ref resetFlow should be called. 
marci@762: 
marci@862:   /// After running an algorithm of the class, the actual flow value 
marci@862:   /// can be obtained by calling \ref flowValue(). The minimum
marci@862:   /// value cut can be written into a \c node map of \c bools by
marci@862:   /// calling \ref minCut. (\ref minMinCut and \ref maxMinCut writes
marci@862:   /// the inclusionwise minimum and maximum of the minimum value
marci@862:   /// cuts, resp.)                                                                                                                               
marci@762:   ///\param Graph The directed graph type the algorithm runs on.
marci@762:   ///\param Num The number type of the capacities and the flow values.
marci@762:   ///\param CapMap The capacity map type.
marci@762:   ///\param FlowMap The flow map type.                                                                                                           
marci@862:   ///\author Marton Makai
marci@762:   template <typename Graph, typename Num,
marci@762: 	    typename CapMap=typename Graph::template EdgeMap<Num>,
marci@762:             typename FlowMap=typename Graph::template EdgeMap<Num> >
marci@762:   class AugmentingFlow {
marci@762:   protected:
marci@762:     typedef typename Graph::Node Node;
marci@762:     typedef typename Graph::NodeIt NodeIt;
marci@762:     typedef typename Graph::EdgeIt EdgeIt;
marci@762:     typedef typename Graph::OutEdgeIt OutEdgeIt;
marci@762:     typedef typename Graph::InEdgeIt InEdgeIt;
marci@762: 
marci@762:     const Graph* g;
marci@762:     Node s;
marci@762:     Node t;
marci@762:     const CapMap* capacity;
marci@762:     FlowMap* flow;
marci@762: //    int n;      //the number of nodes of G
marci@762:     typedef ResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW;   
marci@762:     //typedef ExpResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW;
marci@762:     typedef typename ResGW::OutEdgeIt ResGWOutEdgeIt;
marci@762:     typedef typename ResGW::Edge ResGWEdge;
marci@762:     //typedef typename ResGW::template NodeMap<bool> ReachedMap;
marci@762:     typedef typename Graph::template NodeMap<int> ReachedMap;
marci@762: 
marci@762:     //level works as a bool map in augmenting path algorithms and is
marci@762:     //used by bfs for storing reached information.  In preflow, it
marci@762:     //shows the levels of nodes.     
marci@762:     ReachedMap level;
marci@762: 
marci@762:   public:
marci@762:     ///Indicates the property of the starting flow.
marci@762: 
marci@762:     ///Indicates the property of the starting flow. The meanings are as follows:
marci@762:     ///- \c ZERO_FLOW: constant zero flow
marci@762:     ///- \c GEN_FLOW: any flow, i.e. the sum of the in-flows equals to
marci@762:     ///the sum of the out-flows in every node except the \e source and
marci@762:     ///the \e target.
marci@762:     ///- \c PRE_FLOW: any preflow, i.e. the sum of the in-flows is at 
marci@762:     ///least the sum of the out-flows in every node except the \e source.
marci@762:     ///- \c NO_FLOW: indicates an unspecified edge map. \ref flow will be 
marci@762:     ///set to the constant zero flow in the beginning of the algorithm in this case.
marci@762:     enum FlowEnum{
marci@762:       ZERO_FLOW,
marci@762:       GEN_FLOW,
marci@762:       PRE_FLOW,
marci@762:       NO_FLOW
marci@762:     };
marci@762: 
marci@762:     enum StatusEnum {
marci@762:       AFTER_NOTHING,
marci@762:       AFTER_AUGMENTING,
marci@762:       AFTER_FAST_AUGMENTING, 
marci@762:       AFTER_PRE_FLOW_PHASE_1,      
marci@762:       AFTER_PRE_FLOW_PHASE_2
marci@762:     };
marci@762: 
marci@762:     /// Don not needle this flag only if necessary.
marci@762:     StatusEnum status;
marci@762:     int number_of_augmentations;
marci@762: 
marci@762: 
marci@762:     template<typename IntMap>
marci@762:     class TrickyReachedMap {
marci@762:     protected:
marci@762:       IntMap* map;
marci@762:       int* number_of_augmentations;
marci@762:     public:
alpar@987:       typedef Node Key;
alpar@987:       typedef bool Value;
marci@762:       TrickyReachedMap(IntMap& _map, int& _number_of_augmentations) : 
marci@762: 	map(&_map), number_of_augmentations(&_number_of_augmentations) { }
marci@762:       void set(const Node& n, bool b) {
marci@762: 	if (b)
marci@762: 	  map->set(n, *number_of_augmentations);
marci@762: 	else 
marci@762: 	  map->set(n, *number_of_augmentations-1);
marci@762:       }
marci@762:       bool operator[](const Node& n) const { 
marci@762: 	return (*map)[n]==*number_of_augmentations; 
marci@762:       }
marci@762:     };
marci@762:     
marci@762:     AugmentingFlow(const Graph& _G, Node _s, Node _t, const CapMap& _capacity,
marci@762: 		   FlowMap& _flow) :
marci@762:       g(&_G), s(_s), t(_t), capacity(&_capacity),
marci@762:       flow(&_flow), //n(_G.nodeNum()), 
marci@762:       level(_G), //excess(_G,0), 
marci@762:       status(AFTER_NOTHING), number_of_augmentations(0) { }
marci@762: 
marci@762:     /// Starting from a flow, this method searches for an augmenting path
marci@762:     /// according to the Edmonds-Karp algorithm
marci@762:     /// and augments the flow on if any.
marci@762:     /// The return value shows if the augmentation was succesful.
marci@762:     bool augmentOnShortestPath();
marci@762:     bool augmentOnShortestPath2();
marci@762: 
marci@762:     /// Starting from a flow, this method searches for an augmenting blocking
marci@762:     /// flow according to Dinits' algorithm and augments the flow on if any.
marci@762:     /// The blocking flow is computed in a physically constructed
marci@762:     /// residual graph of type \c Mutablegraph.
marci@762:     /// The return value show sif the augmentation was succesful.
marci@762:     template<typename MutableGraph> bool augmentOnBlockingFlow();
marci@762: 
marci@762:     /// The same as \c augmentOnBlockingFlow<MutableGraph> but the
marci@762:     /// residual graph is not constructed physically.
marci@762:     /// The return value shows if the augmentation was succesful.
marci@762:     bool augmentOnBlockingFlow2();
marci@762: 
marci@762:     template<typename _CutMap>
marci@762:     void actMinCut(_CutMap& M) const {
marci@762:       NodeIt v;
marci@762:       switch (status) {
marci@762: 	case AFTER_PRE_FLOW_PHASE_1:
marci@762: //	std::cout << "AFTER_PRE_FLOW_PHASE_1" << std::endl;
marci@762: // 	for(g->first(v); g->valid(v); g->next(v)) {
marci@762: // 	  if (level[v] < n) {
marci@762: // 	    M.set(v, false);
marci@762: // 	  } else {
marci@762: // 	    M.set(v, true);
marci@762: // 	  }
marci@762: // 	}
marci@762: 	break;
marci@762:       case AFTER_PRE_FLOW_PHASE_2:
marci@762: //	std::cout << "AFTER_PRE_FLOW_PHASE_2" << std::endl;
marci@762: 	break;
marci@762:       case AFTER_NOTHING:
marci@762: //	std::cout << "AFTER_NOTHING" << std::endl;
marci@762: 	minMinCut(M);
marci@762: 	break;
marci@762:       case AFTER_AUGMENTING:
marci@762: //	std::cout << "AFTER_AUGMENTING" << std::endl;
marci@775: 	for(g->first(v); v!=INVALID; ++v) {
marci@762: 	  if (level[v]) {
marci@762: 	    M.set(v, true);
marci@762: 	  } else {
marci@762: 	    M.set(v, false);
marci@762: 	  }
marci@762: 	}
marci@762: 	break;
marci@762:       case AFTER_FAST_AUGMENTING:
marci@762: //	std::cout << "AFTER_FAST_AUGMENTING" << std::endl;
marci@775: 	for(g->first(v); v!=INVALID; ++v) {
marci@762: 	  if (level[v]==number_of_augmentations) {
marci@762: 	    M.set(v, true);
marci@762: 	  } else {
marci@762: 	    M.set(v, false);
marci@762: 	  }
marci@762: 	}
marci@762: 	break;
marci@762:       }
marci@762:     }
marci@762: 
marci@762:     template<typename _CutMap>
marci@762:     void minMinCut(_CutMap& M) const {
marci@762:       std::queue<Node> queue;
marci@762: 
marci@762:       M.set(s,true);
marci@762:       queue.push(s);
marci@762: 
marci@762:       while (!queue.empty()) {
marci@762:         Node w=queue.front();
marci@762: 	queue.pop();
marci@762: 
marci@762: 	OutEdgeIt e;
marci@775: 	for(g->first(e,w) ; e!=INVALID; ++e) {
alpar@986: 	  Node v=g->target(e);
marci@762: 	  if (!M[v] && (*flow)[e] < (*capacity)[e] ) {
marci@762: 	    queue.push(v);
marci@762: 	    M.set(v, true);
marci@762: 	  }
marci@762: 	}
marci@762: 
marci@762: 	InEdgeIt f;
marci@775: 	for(g->first(f,w) ; f!=INVALID; ++f) {
alpar@986: 	  Node v=g->source(f);
marci@762: 	  if (!M[v] && (*flow)[f] > 0 ) {
marci@762: 	    queue.push(v);
marci@762: 	    M.set(v, true);
marci@762: 	  }
marci@762: 	}
marci@762:       }
marci@762:     }
marci@762: 
marci@762:     template<typename _CutMap>
marci@762:     void minMinCut2(_CutMap& M) const {
marci@762:       ResGW res_graph(*g, *capacity, *flow);
marci@762:       BfsIterator<ResGW, _CutMap> bfs(res_graph, M);
marci@762:       bfs.pushAndSetReached(s);
marci@762:       while (!bfs.finished()) ++bfs;
marci@762:     }
marci@762: 
marci@762:     Num flowValue() const {
marci@762:       Num a=0;
marci@777:       for (InEdgeIt e(*g, t); e!=INVALID; ++e) a+=(*flow)[e];
marci@777:       for (OutEdgeIt e(*g, t); e!=INVALID; ++e) a-=(*flow)[e];
marci@762:       return a;
marci@762:       //marci figyu: excess[t] epp ezt adja preflow 1. fazisa utan   
marci@762:     }
marci@762: 
marci@762:   };
marci@762: 
marci@762: 
marci@762: 
marci@762:   template <typename Graph, typename Num, typename CapMap, typename FlowMap>
marci@762:   bool AugmentingFlow<Graph, Num, CapMap, FlowMap>::augmentOnShortestPath()
marci@762:   {
marci@762:     ResGW res_graph(*g, *capacity, *flow);
marci@888:     typename ResGW::ResCap res_cap(res_graph);
marci@888: 
marci@762:     bool _augment=false;
marci@762: 
marci@762:     //ReachedMap level(res_graph);
marci@777:     for (typename Graph::NodeIt n(*g); n!=INVALID; ++n) level.set(n, 0);
marci@762:     BfsIterator<ResGW, ReachedMap> bfs(res_graph, level);
marci@762:     bfs.pushAndSetReached(s);
marci@762: 
marci@762:     typename ResGW::template NodeMap<ResGWEdge> pred(res_graph);
marci@762:     pred.set(s, INVALID);
marci@762: 
marci@762:     typename ResGW::template NodeMap<Num> free(res_graph);
marci@762: 
marci@762:     //searching for augmenting path
marci@762:     while ( !bfs.finished() ) {
marci@777:       ResGWEdge e=bfs;
marci@775:       if (e!=INVALID && bfs.isBNodeNewlyReached()) {
alpar@986: 	Node v=res_graph.source(e);
alpar@986: 	Node w=res_graph.target(e);
marci@762: 	pred.set(w, e);
marci@775: 	if (pred[v]!=INVALID) {
marci@888: 	  free.set(w, std::min(free[v], res_cap[e]));
marci@762: 	} else {
marci@888: 	  free.set(w, res_cap[e]);
marci@762: 	}
alpar@986: 	if (res_graph.target(e)==t) { _augment=true; break; }
marci@762:       }
marci@762: 
marci@762:       ++bfs;
marci@762:     } //end of searching augmenting path
marci@762: 
marci@762:     if (_augment) {
marci@762:       Node n=t;
marci@762:       Num augment_value=free[t];
marci@775:       while (pred[n]!=INVALID) {
marci@762: 	ResGWEdge e=pred[n];
marci@762: 	res_graph.augment(e, augment_value);
alpar@986: 	n=res_graph.source(e);
marci@762:       }
marci@762:     }
marci@762: 
marci@762:     status=AFTER_AUGMENTING;
marci@762:     return _augment;
marci@762:   }
marci@762: 
marci@762:   template <typename Graph, typename Num, typename CapMap, typename FlowMap>
marci@762:   bool AugmentingFlow<Graph, Num, CapMap, FlowMap>::augmentOnShortestPath2()
marci@762:   {
marci@762:     ResGW res_graph(*g, *capacity, *flow);
marci@888:     typename ResGW::ResCap res_cap(res_graph);
marci@888: 
marci@762:     bool _augment=false;
marci@762: 
marci@762:     if (status!=AFTER_FAST_AUGMENTING) {
marci@777:       for (typename Graph::NodeIt n(*g); n!=INVALID; ++n) level.set(n, 0); 
marci@762:       number_of_augmentations=1;
marci@762:     } else {
marci@762:       ++number_of_augmentations;
marci@762:     }
marci@762:     TrickyReachedMap<ReachedMap> 
marci@762:       tricky_reached_map(level, number_of_augmentations);
marci@762:     //ReachedMap level(res_graph);
marci@762: //    FOR_EACH_LOC(typename Graph::NodeIt, e, *g) level.set(e, 0);
marci@762:     BfsIterator<ResGW, TrickyReachedMap<ReachedMap> > 
marci@762:       bfs(res_graph, tricky_reached_map);
marci@762:     bfs.pushAndSetReached(s);
marci@762: 
marci@762:     typename ResGW::template NodeMap<ResGWEdge> pred(res_graph);
marci@762:     pred.set(s, INVALID);
marci@762: 
marci@762:     typename ResGW::template NodeMap<Num> free(res_graph);
marci@762: 
marci@762:     //searching for augmenting path
marci@762:     while ( !bfs.finished() ) {
marci@777:       ResGWEdge e=bfs;
marci@775:       if (e!=INVALID && bfs.isBNodeNewlyReached()) {
alpar@986: 	Node v=res_graph.source(e);
alpar@986: 	Node w=res_graph.target(e);
marci@762: 	pred.set(w, e);
marci@775: 	if (pred[v]!=INVALID) {
marci@888: 	  free.set(w, std::min(free[v], res_cap[e]));
marci@762: 	} else {
marci@888: 	  free.set(w, res_cap[e]);
marci@762: 	}
alpar@986: 	if (res_graph.target(e)==t) { _augment=true; break; }
marci@762:       }
marci@762: 
marci@762:       ++bfs;
marci@762:     } //end of searching augmenting path
marci@762: 
marci@762:     if (_augment) {
marci@762:       Node n=t;
marci@762:       Num augment_value=free[t];
marci@775:       while (pred[n]!=INVALID) {
marci@762: 	ResGWEdge e=pred[n];
marci@762: 	res_graph.augment(e, augment_value);
alpar@986: 	n=res_graph.source(e);
marci@762:       }
marci@762:     }
marci@762: 
marci@762:     status=AFTER_FAST_AUGMENTING;
marci@762:     return _augment;
marci@762:   }
marci@762: 
marci@762: 
marci@762:   template <typename Graph, typename Num, typename CapMap, typename FlowMap>
marci@762:   template<typename MutableGraph>
marci@762:   bool AugmentingFlow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow()
marci@762:   {
marci@762:     typedef MutableGraph MG;
marci@762:     bool _augment=false;
marci@762: 
marci@762:     ResGW res_graph(*g, *capacity, *flow);
marci@888:     typename ResGW::ResCap res_cap(res_graph);
marci@762: 
marci@762:     //bfs for distances on the residual graph
marci@762:     //ReachedMap level(res_graph);
marci@777:     for (typename Graph::NodeIt n(*g); n!=INVALID; ++n) level.set(n, 0);
marci@762:     BfsIterator<ResGW, ReachedMap> bfs(res_graph, level);
marci@762:     bfs.pushAndSetReached(s);
marci@762:     typename ResGW::template NodeMap<int>
marci@762:       dist(res_graph); //filled up with 0's
marci@762: 
marci@762:     //F will contain the physical copy of the residual graph
marci@762:     //with the set of edges which are on shortest paths
marci@762:     MG F;
marci@762:     typename ResGW::template NodeMap<typename MG::Node>
marci@762:       res_graph_to_F(res_graph);
marci@762:     {
marci@970:       for(typename ResGW::NodeIt n(res_graph); n!=INVALID; ++n) 
marci@762: 	res_graph_to_F.set(n, F.addNode());
marci@762:     }
marci@762: 
marci@762:     typename MG::Node sF=res_graph_to_F[s];
marci@762:     typename MG::Node tF=res_graph_to_F[t];
marci@762:     typename MG::template EdgeMap<ResGWEdge> original_edge(F);
marci@762:     typename MG::template EdgeMap<Num> residual_capacity(F);
marci@762: 
marci@762:     while ( !bfs.finished() ) {
marci@777:       ResGWEdge e=bfs;
marci@775:       if (e!=INVALID) {
marci@762: 	if (bfs.isBNodeNewlyReached()) {
alpar@986: 	  dist.set(res_graph.target(e), dist[res_graph.source(e)]+1);
alpar@986: 	  typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.source(e)],
alpar@986: 					res_graph_to_F[res_graph.target(e)]);
marci@854: 	  //original_edge.update();
marci@762: 	  original_edge.set(f, e);
marci@854: 	  //residual_capacity.update();
marci@888: 	  residual_capacity.set(f, res_cap[e]);
marci@762: 	} else {
alpar@986: 	  if (dist[res_graph.target(e)]==(dist[res_graph.source(e)]+1)) {
alpar@986: 	    typename MG::Edge f=F.addEdge(res_graph_to_F[res_graph.source(e)],
alpar@986: 					  res_graph_to_F[res_graph.target(e)]);
marci@854: 	    //original_edge.update();
marci@762: 	    original_edge.set(f, e);
marci@854: 	    //residual_capacity.update();
marci@888: 	    residual_capacity.set(f, res_cap[e]);
marci@762: 	  }
marci@762: 	}
marci@762:       }
marci@762:       ++bfs;
marci@762:     } //computing distances from s in the residual graph
marci@762: 
marci@762:     bool __augment=true;
marci@762: 
marci@762:     while (__augment) {
marci@762:       __augment=false;
marci@762:       //computing blocking flow with dfs
marci@762:       DfsIterator< MG, typename MG::template NodeMap<bool> > dfs(F);
marci@762:       typename MG::template NodeMap<typename MG::Edge> pred(F);
marci@762:       pred.set(sF, INVALID);
marci@762:       //invalid iterators for sources
marci@762: 
marci@762:       typename MG::template NodeMap<Num> free(F);
marci@762: 
marci@762:       dfs.pushAndSetReached(sF);
marci@762:       while (!dfs.finished()) {
marci@762: 	++dfs;
marci@854: 	if (typename MG::Edge(dfs)!=INVALID) {
marci@762: 	  if (dfs.isBNodeNewlyReached()) {
alpar@986: 	    typename MG::Node v=F.source(dfs);
alpar@986: 	    typename MG::Node w=F.target(dfs);
marci@762: 	    pred.set(w, dfs);
marci@775: 	    if (pred[v]!=INVALID) {
marci@762: 	      free.set(w, std::min(free[v], residual_capacity[dfs]));
marci@762: 	    } else {
marci@762: 	      free.set(w, residual_capacity[dfs]);
marci@762: 	    }
marci@762: 	    if (w==tF) {
marci@762: 	      __augment=true;
marci@762: 	      _augment=true;
marci@762: 	      break;
marci@762: 	    }
marci@762: 
marci@762: 	  } else {
marci@854: 	    F.erase(typename MG::Edge(dfs));
marci@762: 	  }
marci@762: 	}
marci@762:       }
marci@762: 
marci@762:       if (__augment) {
marci@762: 	typename MG::Node n=tF;
marci@762: 	Num augment_value=free[tF];
marci@775: 	while (pred[n]!=INVALID) {
marci@762: 	  typename MG::Edge e=pred[n];
marci@762: 	  res_graph.augment(original_edge[e], augment_value);
alpar@986: 	  n=F.source(e);
marci@762: 	  if (residual_capacity[e]==augment_value)
marci@762: 	    F.erase(e);
marci@762: 	  else
marci@762: 	    residual_capacity.set(e, residual_capacity[e]-augment_value);
marci@762: 	}
marci@762:       }
marci@762: 
marci@762:     }
marci@762: 
marci@762:     status=AFTER_AUGMENTING;
marci@762:     return _augment;
marci@762:   }
marci@762: 
marci@862:   /// Blocking flow augmentation without constructing the layered 
marci@862:   /// graph physically in which the blocking flow is computed.
marci@762:   template <typename Graph, typename Num, typename CapMap, typename FlowMap>
marci@762:   bool AugmentingFlow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow2()
marci@762:   {
marci@762:     bool _augment=false;
marci@762: 
marci@762:     ResGW res_graph(*g, *capacity, *flow);
marci@888:     typename ResGW::ResCap res_cap(res_graph);
marci@762: 
marci@862:     //Potential map, for distances from s
marci@862:     typename ResGW::template NodeMap<int> potential(res_graph, 0);
marci@862:     typedef ConstMap<typename ResGW::Edge, int> Const1Map; 
marci@862:     Const1Map const_1_map(1);
marci@862:     TightEdgeFilterMap<ResGW, typename ResGW::template NodeMap<int>,
marci@862:       Const1Map> tight_edge_filter(res_graph, potential, const_1_map);
marci@862: 
marci@777:     for (typename Graph::NodeIt n(*g); n!=INVALID; ++n) level.set(n, 0);
marci@762:     BfsIterator<ResGW, ReachedMap> bfs(res_graph, level);
marci@862:     bfs.pushAndSetReached(s);
marci@762: 
marci@862:     //computing distances from s in the residual graph
marci@762:     while ( !bfs.finished() ) {
marci@777:       ResGWEdge e=bfs;
marci@862:       if (e!=INVALID && bfs.isBNodeNewlyReached())
alpar@986: 	potential.set(res_graph.target(e), potential[res_graph.source(e)]+1);
marci@762:       ++bfs;
marci@862:     } 
marci@762: 
marci@862:     //Subgraph containing the edges on some shortest paths 
marci@862:     //(i.e. tight edges)
marci@762:     ConstMap<typename ResGW::Node, bool> true_map(true);
marci@762:     typedef SubGraphWrapper<ResGW, ConstMap<typename ResGW::Node, bool>,
marci@862:       TightEdgeFilterMap<ResGW, typename ResGW::template NodeMap<int>, 
marci@862:       Const1Map> > FilterResGW;
marci@862:     FilterResGW filter_res_graph(res_graph, true_map, tight_edge_filter);
marci@762: 
marci@762:     //Subgraph, which is able to delete edges which are already
marci@762:     //met by the dfs
marci@777:     typename FilterResGW::template NodeMap<typename FilterResGW::Edge>
marci@762:       first_out_edges(filter_res_graph);
marci@862:     for (typename FilterResGW::NodeIt v(filter_res_graph); v!=INVALID; ++v)
marci@862:       first_out_edges.set
marci@862: 	(v, typename FilterResGW::OutEdgeIt(filter_res_graph, v));
marci@862: 
marci@762:     typedef ErasingFirstGraphWrapper<FilterResGW, typename FilterResGW::
marci@777:       template NodeMap<typename FilterResGW::Edge> > ErasingResGW;
marci@762:     ErasingResGW erasing_res_graph(filter_res_graph, first_out_edges);
marci@762: 
marci@762:     bool __augment=true;
marci@762: 
marci@762:     while (__augment) {
marci@762: 
marci@762:       __augment=false;
marci@762:       //computing blocking flow with dfs
marci@762:       DfsIterator< ErasingResGW,
marci@762: 	typename ErasingResGW::template NodeMap<bool> >
marci@762: 	dfs(erasing_res_graph);
marci@762:       typename ErasingResGW::
marci@777: 	template NodeMap<typename ErasingResGW::Edge> pred(erasing_res_graph);
marci@762:       pred.set(s, INVALID);
marci@762:       //invalid iterators for sources
marci@762: 
marci@762:       typename ErasingResGW::template NodeMap<Num>
marci@762: 	free1(erasing_res_graph);
marci@762: 
marci@762:       dfs.pushAndSetReached
alpar@921: 	/// \bug lemon 0.2
marci@762: 	(typename ErasingResGW::Node
marci@762: 	 (typename FilterResGW::Node
marci@762: 	  (typename ResGW::Node(s)
marci@762: 	   )
marci@762: 	  )
marci@762: 	 );
marci@777: 	
marci@762:       while (!dfs.finished()) {
marci@762: 	++dfs;
marci@862: 	if (typename ErasingResGW::Edge(dfs)!=INVALID) {
marci@862: 	  if (dfs.isBNodeNewlyReached()) {
marci@862: 	    
alpar@986: 	    typename ErasingResGW::Node v=erasing_res_graph.source(dfs);
alpar@986: 	    typename ErasingResGW::Node w=erasing_res_graph.target(dfs);
marci@762: 
marci@862: 	    pred.set(w, typename ErasingResGW::Edge(dfs));
marci@862: 	    if (pred[v]!=INVALID) {
marci@862: 	      free1.set
marci@888: 		(w, std::min(free1[v], res_cap
marci@888: 			     [typename ErasingResGW::Edge(dfs)]));
marci@862: 	    } else {
marci@862: 	      free1.set
marci@888: 		(w, res_cap
marci@888: 		 [typename ErasingResGW::Edge(dfs)]);
marci@862: 	    }
marci@762: 
marci@862: 	    if (w==t) {
marci@862: 	      __augment=true;
marci@862: 	      _augment=true;
marci@862: 	      break;
marci@762: 	    }
marci@862: 	  } else {
marci@862: 	    erasing_res_graph.erase(dfs);
marci@762: 	  }
marci@862: 	}
marci@762:       }
marci@762: 
marci@762:       if (__augment) {
marci@762: 	typename ErasingResGW::Node
marci@762: 	  n=typename FilterResGW::Node(typename ResGW::Node(t));
marci@762: 	Num augment_value=free1[n];
marci@777: 	while (pred[n]!=INVALID) {
marci@777: 	  typename ErasingResGW::Edge e=pred[n];
marci@762: 	  res_graph.augment(e, augment_value);
alpar@986: 	  n=erasing_res_graph.source(e);
marci@888: 	  if (res_cap[e]==0)
marci@762: 	    erasing_res_graph.erase(e);
marci@762: 	}
marci@762:       }
marci@762: 
marci@762:     } //while (__augment)
marci@762: 
marci@762:     status=AFTER_AUGMENTING;
marci@762:     return _augment;
marci@762:   }
marci@762: 
marci@762: 
alpar@921: } //namespace lemon
marci@762: 
alpar@921: #endif //LEMON_AUGMENTING_FLOW_H
marci@762: 
marci@762: