bug correction
authormarci
Thu, 16 Sep 2004 10:59:52 +0000
changeset 862732f2acb7239
parent 861 021e513a2d83
child 863 d27bbe17b0b8
bug correction
src/work/marci/augmenting_flow.h
     1.1 --- a/src/work/marci/augmenting_flow.h	Thu Sep 16 10:59:30 2004 +0000
     1.2 +++ b/src/work/marci/augmenting_flow.h	Thu Sep 16 10:59:52 2004 +0000
     1.3 @@ -3,15 +3,14 @@
     1.4  #define HUGO_AUGMENTING_FLOW_H
     1.5  
     1.6  #include <vector>
     1.7 -#include <queue>
     1.8 -#include <stack>
     1.9 +//#include <queue>
    1.10 +//#include <stack>
    1.11  #include <iostream>
    1.12  
    1.13  #include <hugo/graph_wrapper.h>
    1.14  #include <bfs_dfs.h>
    1.15  #include <hugo/invalid.h>
    1.16  #include <hugo/maps.h>
    1.17 -//#include <for_each_macros.h>
    1.18  
    1.19  /// \file
    1.20  /// \brief Maximum flow algorithms.
    1.21 @@ -19,849 +18,68 @@
    1.22  
    1.23  namespace hugo {
    1.24  
    1.25 +  /// \brief A map for filtering the edge-set to those edges 
    1.26 +  /// which are tight w.r.t. some node_potential map and 
    1.27 +  /// edge_distance map.
    1.28 +  ///
    1.29 +  /// A node-map node_potential is said to be a potential w.r.t. 
    1.30 +  /// an edge-map edge_distance 
    1.31 +  /// if and only if for each edge e, node_potential[g.head(e)] 
    1.32 +  /// <= edge_distance[e]+node_potential[g.tail(e)] 
    1.33 +  /// (or the reverse inequality holds for each edge).
    1.34 +  /// An edge is said to be tight if this inequality holds with equality, 
    1.35 +  /// and the map returns true exactly for those edges.
    1.36 +  /// To avoid rounding errors, it is recommended to use this class with exact 
    1.37 +  /// types, e.g. with int.
    1.38 +  template<typename Graph, 
    1.39 +	   typename NodePotentialMap, typename EdgeDistanceMap>
    1.40 +  class TightEdgeFilterMap {
    1.41 +  protected:
    1.42 +    const Graph* g;
    1.43 +    NodePotentialMap* node_potential;
    1.44 +    EdgeDistanceMap* edge_distance;
    1.45 +  public:
    1.46 +    TightEdgeFilterMap(Graph& _g, NodePotentialMap& _node_potential, 
    1.47 +		       EdgeDistanceMap& _edge_distance) : 
    1.48 +      g(&_g), node_potential(&_node_potential), 
    1.49 +      edge_distance(&_edge_distance) { }
    1.50 +//     void set(const typename Graph::Node& n, int a) {
    1.51 +//       pot->set(n, a);
    1.52 +//     }
    1.53 +//     int operator[](const typename Graph::Node& n) const { 
    1.54 +//       return (*node_potential)[n]; 
    1.55 +//     }
    1.56 +    bool operator[](const typename Graph::Edge& e) const {
    1.57 +      return ((*node_potential)[g->head(e)] == 
    1.58 +	      (*edge_distance)[e]+(*node_potential)[g->tail(e)]);
    1.59 +    }
    1.60 +  };
    1.61 +
    1.62    /// \addtogroup galgs
    1.63    /// @{                                                                                                                                        
    1.64 -  ///Maximum flow algorithms class.
    1.65 +  /// Class for augmenting path flow algorithms.
    1.66  
    1.67 -  ///This class provides various algorithms for finding a flow of
    1.68 -  ///maximum value in a directed graph. The \e source node, the \e
    1.69 -  ///target node, the \e capacity of the edges and the \e starting \e
    1.70 -  ///flow value of the edges should be passed to the algorithm through the
    1.71 -  ///constructor. It is possible to change these quantities using the
    1.72 -  ///functions \ref resetSource, \ref resetTarget, \ref resetCap and
    1.73 -  ///\ref resetFlow. Before any subsequent runs of any algorithm of
    1.74 -  ///the class \ref resetFlow should be called. 
    1.75 +  /// This class provides various algorithms for finding a flow of
    1.76 +  /// maximum value in a directed graph. The \e source node, the \e
    1.77 +  /// target node, the \e capacity of the edges and the \e starting \e
    1.78 +  /// flow value of the edges should be passed to the algorithm through the
    1.79 +  /// constructor. 
    1.80 +//   /// It is possible to change these quantities using the
    1.81 +//   /// functions \ref resetSource, \ref resetTarget, \ref resetCap and
    1.82 +//   /// \ref resetFlow. Before any subsequent runs of any algorithm of
    1.83 +//   /// the class \ref resetFlow should be called. 
    1.84  
    1.85 -  ///After running an algorithm of the class, the actual flow value 
    1.86 -  ///can be obtained by calling \ref flowValue(). The minimum
    1.87 -  ///value cut can be written into a \c node map of \c bools by
    1.88 -  ///calling \ref minCut. (\ref minMinCut and \ref maxMinCut writes
    1.89 -  ///the inclusionwise minimum and maximum of the minimum value
    1.90 -  ///cuts, resp.)                                                                                                                               
    1.91 +  /// After running an algorithm of the class, the actual flow value 
    1.92 +  /// can be obtained by calling \ref flowValue(). The minimum
    1.93 +  /// value cut can be written into a \c node map of \c bools by
    1.94 +  /// calling \ref minCut. (\ref minMinCut and \ref maxMinCut writes
    1.95 +  /// the inclusionwise minimum and maximum of the minimum value
    1.96 +  /// cuts, resp.)                                                                                                                               
    1.97    ///\param Graph The directed graph type the algorithm runs on.
    1.98    ///\param Num The number type of the capacities and the flow values.
    1.99    ///\param CapMap The capacity map type.
   1.100    ///\param FlowMap The flow map type.                                                                                                           
   1.101 -  ///\author Marton Makai, Jacint Szabo 
   1.102 -//   template <typename Graph, typename Num,
   1.103 -// 	    typename CapMap=typename Graph::template EdgeMap<Num>,
   1.104 -//             typename FlowMap=typename Graph::template EdgeMap<Num> >
   1.105 -//   class MaxFlow {
   1.106 -//   protected:
   1.107 -//     typedef typename Graph::Node Node;
   1.108 -//     typedef typename Graph::NodeIt NodeIt;
   1.109 -//     typedef typename Graph::EdgeIt EdgeIt;
   1.110 -//     typedef typename Graph::OutEdgeIt OutEdgeIt;
   1.111 -//     typedef typename Graph::InEdgeIt InEdgeIt;
   1.112 -
   1.113 -//     typedef typename std::vector<std::stack<Node> > VecStack;
   1.114 -//     typedef typename Graph::template NodeMap<Node> NNMap;
   1.115 -//     typedef typename std::vector<Node> VecNode;
   1.116 -
   1.117 -//     const Graph* g;
   1.118 -//     Node s;
   1.119 -//     Node t;
   1.120 -//     const CapMap* capacity;
   1.121 -//     FlowMap* flow;
   1.122 -//     int n;      //the number of nodes of G
   1.123 -//     typedef ResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW;   
   1.124 -//     //typedef ExpResGraphWrapper<const Graph, Num, CapMap, FlowMap> ResGW;
   1.125 -//     typedef typename ResGW::OutEdgeIt ResGWOutEdgeIt;
   1.126 -//     typedef typename ResGW::Edge ResGWEdge;
   1.127 -//     //typedef typename ResGW::template NodeMap<bool> ReachedMap;
   1.128 -//     typedef typename Graph::template NodeMap<int> ReachedMap;
   1.129 -
   1.130 -
   1.131 -//     //level works as a bool map in augmenting path algorithms and is
   1.132 -//     //used by bfs for storing reached information.  In preflow, it
   1.133 -//     //shows the levels of nodes.     
   1.134 -//     ReachedMap level;
   1.135 -
   1.136 -//     //excess is needed only in preflow
   1.137 -//     typename Graph::template NodeMap<Num> excess;
   1.138 -
   1.139 -//     //fixme    
   1.140 -// //   protected:
   1.141 -//     //     MaxFlow() { }
   1.142 -//     //     void set(const Graph& _G, Node _s, Node _t, const CapMap& _capacity,
   1.143 -//     // 	     FlowMap& _flow)
   1.144 -//     //       {
   1.145 -//     // 	g=&_G;
   1.146 -//     // 	s=_s;
   1.147 -//     // 	t=_t;
   1.148 -//     // 	capacity=&_capacity;
   1.149 -//     // 	flow=&_flow;
   1.150 -//     // 	n=_G.nodeNum;
   1.151 -//     // 	level.set (_G); //kellene vmi ilyesmi fv
   1.152 -//     // 	excess(_G,0); //itt is
   1.153 -//     //       }
   1.154 -
   1.155 -//     // constants used for heuristics
   1.156 -//     static const int H0=20;
   1.157 -//     static const int H1=1;
   1.158 -
   1.159 -//   public:
   1.160 -
   1.161 -//     ///Indicates the property of the starting flow.
   1.162 -
   1.163 -//     ///Indicates the property of the starting flow. The meanings are as follows:
   1.164 -//     ///- \c ZERO_FLOW: constant zero flow
   1.165 -//     ///- \c GEN_FLOW: any flow, i.e. the sum of the in-flows equals to
   1.166 -//     ///the sum of the out-flows in every node except the \e source and
   1.167 -//     ///the \e target.
   1.168 -//     ///- \c PRE_FLOW: any preflow, i.e. the sum of the in-flows is at 
   1.169 -//     ///least the sum of the out-flows in every node except the \e source.
   1.170 -//     ///- \c NO_FLOW: indicates an unspecified edge map. \ref flow will be 
   1.171 -//     ///set to the constant zero flow in the beginning of the algorithm in this case.
   1.172 -//     enum FlowEnum{
   1.173 -//       ZERO_FLOW,
   1.174 -//       GEN_FLOW,
   1.175 -//       PRE_FLOW,
   1.176 -//       NO_FLOW
   1.177 -//     };
   1.178 -
   1.179 -//     enum StatusEnum {
   1.180 -//       AFTER_NOTHING,
   1.181 -//       AFTER_AUGMENTING,
   1.182 -//       AFTER_FAST_AUGMENTING, 
   1.183 -//       AFTER_PRE_FLOW_PHASE_1,      
   1.184 -//       AFTER_PRE_FLOW_PHASE_2
   1.185 -//     };
   1.186 -
   1.187 -//     /// Don not needle this flag only if necessary.
   1.188 -//     StatusEnum status;
   1.189 -// //     int number_of_augmentations;
   1.190 -
   1.191 -
   1.192 -// //     template<typename IntMap>
   1.193 -// //     class TrickyReachedMap {
   1.194 -// //     protected:
   1.195 -// //       IntMap* map;
   1.196 -// //       int* number_of_augmentations;
   1.197 -// //     public:
   1.198 -// //       TrickyReachedMap(IntMap& _map, int& _number_of_augmentations) : 
   1.199 -// // 	map(&_map), number_of_augmentations(&_number_of_augmentations) { }
   1.200 -// //       void set(const Node& n, bool b) {
   1.201 -// // 	if (b)
   1.202 -// // 	  map->set(n, *number_of_augmentations);
   1.203 -// // 	else 
   1.204 -// // 	  map->set(n, *number_of_augmentations-1);
   1.205 -// //       }
   1.206 -// //       bool operator[](const Node& n) const { 
   1.207 -// // 	return (*map)[n]==*number_of_augmentations; 
   1.208 -// //       }
   1.209 -// //     };
   1.210 -    
   1.211 -//     MaxFlow(const Graph& _G, Node _s, Node _t, const CapMap& _capacity,
   1.212 -// 	    FlowMap& _flow) :
   1.213 -//       g(&_G), s(_s), t(_t), capacity(&_capacity),
   1.214 -//       flow(&_flow), n(_G.nodeNum()), level(_G), excess(_G,0), 
   1.215 -//       status(AFTER_NOTHING) { }
   1.216 -
   1.217 -//     ///Runs a maximum flow algorithm.
   1.218 -
   1.219 -//     ///Runs a preflow algorithm, which is the fastest maximum flow
   1.220 -//     ///algorithm up-to-date. The default for \c fe is ZERO_FLOW.
   1.221 -//     ///\pre The starting flow must be
   1.222 -//     /// - a constant zero flow if \c fe is \c ZERO_FLOW,
   1.223 -//     /// - an arbitary flow if \c fe is \c GEN_FLOW,
   1.224 -//     /// - an arbitary preflow if \c fe is \c PRE_FLOW,
   1.225 -//     /// - any map if \c fe is NO_FLOW.
   1.226 -//     void run(FlowEnum fe=ZERO_FLOW) {
   1.227 -//       preflow(fe);
   1.228 -//     }
   1.229 -
   1.230 -                                                                              
   1.231 -//     ///Runs a preflow algorithm.  
   1.232 -
   1.233 -//     ///Runs a preflow algorithm. The preflow algorithms provide the
   1.234 -//     ///fastest way to compute a maximum flow in a directed graph.
   1.235 -//     ///\pre The starting flow must be
   1.236 -//     /// - a constant zero flow if \c fe is \c ZERO_FLOW,
   1.237 -//     /// - an arbitary flow if \c fe is \c GEN_FLOW,
   1.238 -//     /// - an arbitary preflow if \c fe is \c PRE_FLOW,
   1.239 -//     /// - any map if \c fe is NO_FLOW.
   1.240 -//     void preflow(FlowEnum fe) {
   1.241 -//       preflowPhase1(fe);
   1.242 -//       preflowPhase2();
   1.243 -//     }
   1.244 -//     // Heuristics:
   1.245 -//     //   2 phase
   1.246 -//     //   gap
   1.247 -//     //   list 'level_list' on the nodes on level i implemented by hand
   1.248 -//     //   stack 'active' on the active nodes on level i                                                                                    
   1.249 -//     //   runs heuristic 'highest label' for H1*n relabels
   1.250 -//     //   runs heuristic 'bound decrease' for H0*n relabels, starts with 'highest label'
   1.251 -//     //   Parameters H0 and H1 are initialized to 20 and 1.
   1.252 -
   1.253 -//     ///Runs the first phase of the preflow algorithm.
   1.254 -
   1.255 -//     ///The preflow algorithm consists of two phases, this method runs the
   1.256 -//     ///first phase. After the first phase the maximum flow value and a
   1.257 -//     ///minimum value cut can already be computed, though a maximum flow
   1.258 -//     ///is net yet obtained. So after calling this method \ref flowValue
   1.259 -//     ///and \ref actMinCut gives proper results.
   1.260 -//     ///\warning: \ref minCut, \ref minMinCut and \ref maxMinCut do not
   1.261 -//     ///give minimum value cuts unless calling \ref preflowPhase2.
   1.262 -//     ///\pre The starting flow must be
   1.263 -//     /// - a constant zero flow if \c fe is \c ZERO_FLOW,
   1.264 -//     /// - an arbitary flow if \c fe is \c GEN_FLOW,
   1.265 -//     /// - an arbitary preflow if \c fe is \c PRE_FLOW,
   1.266 -//     /// - any map if \c fe is NO_FLOW.
   1.267 -//     void preflowPhase1(FlowEnum fe);
   1.268 -
   1.269 -//     ///Runs the second phase of the preflow algorithm.
   1.270 -
   1.271 -//     ///The preflow algorithm consists of two phases, this method runs
   1.272 -//     ///the second phase. After calling \ref preflowPhase1 and then
   1.273 -//     ///\ref preflowPhase2 the methods \ref flowValue, \ref minCut,
   1.274 -//     ///\ref minMinCut and \ref maxMinCut give proper results.
   1.275 -//     ///\pre \ref preflowPhase1 must be called before.
   1.276 -//     void preflowPhase2();
   1.277 -
   1.278 -//     /// Returns the maximum value of a flow.
   1.279 -
   1.280 -//     /// Returns the maximum value of a flow, by counting the 
   1.281 -//     /// over-flow of the target node \ref t.
   1.282 -//     /// It can be called already after running \ref preflowPhase1.
   1.283 -//     Num flowValue() const {
   1.284 -//       Num a=0;
   1.285 -//       FOR_EACH_INC_LOC(InEdgeIt, e, *g, t) a+=(*flow)[e];
   1.286 -//       FOR_EACH_INC_LOC(OutEdgeIt, e, *g, t) a-=(*flow)[e];
   1.287 -//       return a;
   1.288 -//       //marci figyu: excess[t] epp ezt adja preflow 1. fazisa utan   
   1.289 -//     }
   1.290 -
   1.291 -//     ///Returns a minimum value cut after calling \ref preflowPhase1.
   1.292 -
   1.293 -//     ///After the first phase of the preflow algorithm the maximum flow
   1.294 -//     ///value and a minimum value cut can already be computed. This
   1.295 -//     ///method can be called after running \ref preflowPhase1 for
   1.296 -//     ///obtaining a minimum value cut.
   1.297 -//     /// \warning Gives proper result only right after calling \ref
   1.298 -//     /// preflowPhase1.
   1.299 -//     /// \todo We have to make some status variable which shows the
   1.300 -//     /// actual state
   1.301 -//     /// of the class. This enables us to determine which methods are valid
   1.302 -//     /// for MinCut computation
   1.303 -//     template<typename _CutMap>
   1.304 -//     void actMinCut(_CutMap& M) const {
   1.305 -//       NodeIt v;
   1.306 -//       switch (status) {
   1.307 -//       case AFTER_PRE_FLOW_PHASE_1:
   1.308 -// 	for(g->first(v); g->valid(v); g->next(v)) {
   1.309 -// 	  if (level[v] < n) {
   1.310 -// 	    M.set(v, false);
   1.311 -// 	  } else {
   1.312 -// 	    M.set(v, true);
   1.313 -// 	  }
   1.314 -// 	}
   1.315 -// 	break;
   1.316 -//       case AFTER_PRE_FLOW_PHASE_2:
   1.317 -//       case AFTER_NOTHING:
   1.318 -//       case AFTER_AUGMENTING:
   1.319 -//       case AFTER_FAST_AUGMENTING:
   1.320 -// 	minMinCut(M);
   1.321 -// 	break;
   1.322 -// //       case AFTER_AUGMENTING:
   1.323 -// // 	for(g->first(v); g->valid(v); g->next(v)) {
   1.324 -// // 	  if (level[v]) {
   1.325 -// // 	    M.set(v, true);
   1.326 -// // 	  } else {
   1.327 -// // 	    M.set(v, false);
   1.328 -// // 	  }
   1.329 -// // 	}
   1.330 -// // 	break;
   1.331 -// //       case AFTER_FAST_AUGMENTING:
   1.332 -// // 	for(g->first(v); g->valid(v); g->next(v)) {
   1.333 -// // 	  if (level[v]==number_of_augmentations) {
   1.334 -// // 	    M.set(v, true);
   1.335 -// // 	  } else {
   1.336 -// // 	    M.set(v, false);
   1.337 -// // 	  }
   1.338 -// // 	}
   1.339 -// // 	break;
   1.340 -//       }
   1.341 -//     }
   1.342 -
   1.343 -//     ///Returns the inclusionwise minimum of the minimum value cuts.
   1.344 -
   1.345 -//     ///Sets \c M to the characteristic vector of the minimum value cut
   1.346 -//     ///which is inclusionwise minimum. It is computed by processing
   1.347 -//     ///a bfs from the source node \c s in the residual graph.
   1.348 -//     ///\pre M should be a node map of bools initialized to false.
   1.349 -//     ///\pre \c flow must be a maximum flow.
   1.350 -//     template<typename _CutMap>
   1.351 -//     void minMinCut(_CutMap& M) const {
   1.352 -//       std::queue<Node> queue;
   1.353 -
   1.354 -//       M.set(s,true);
   1.355 -//       queue.push(s);
   1.356 -
   1.357 -//       while (!queue.empty()) {
   1.358 -//         Node w=queue.front();
   1.359 -// 	queue.pop();
   1.360 -
   1.361 -// 	OutEdgeIt e;
   1.362 -// 	for(g->first(e,w) ; g->valid(e); g->next(e)) {
   1.363 -// 	  Node v=g->head(e);
   1.364 -// 	  if (!M[v] && (*flow)[e] < (*capacity)[e] ) {
   1.365 -// 	    queue.push(v);
   1.366 -// 	    M.set(v, true);
   1.367 -// 	  }
   1.368 -// 	}
   1.369 -
   1.370 -// 	InEdgeIt f;
   1.371 -// 	for(g->first(f,w) ; g->valid(f); g->next(f)) {
   1.372 -// 	  Node v=g->tail(f);
   1.373 -// 	  if (!M[v] && (*flow)[f] > 0 ) {
   1.374 -// 	    queue.push(v);
   1.375 -// 	    M.set(v, true);
   1.376 -// 	  }
   1.377 -// 	}
   1.378 -//       }
   1.379 -//     }
   1.380 -
   1.381 -//     ///Returns the inclusionwise maximum of the minimum value cuts.
   1.382 -
   1.383 -//     ///Sets \c M to the characteristic vector of the minimum value cut
   1.384 -//     ///which is inclusionwise maximum. It is computed by processing a
   1.385 -//     ///backward bfs from the target node \c t in the residual graph.
   1.386 -//     ///\pre M should be a node map of bools initialized to false.
   1.387 -//     ///\pre \c flow must be a maximum flow. 
   1.388 -//     template<typename _CutMap>
   1.389 -//     void maxMinCut(_CutMap& M) const {
   1.390 -
   1.391 -//       NodeIt v;
   1.392 -//       for(g->first(v) ; g->valid(v); g->next(v)) {
   1.393 -// 	M.set(v, true);
   1.394 -//       }
   1.395 -
   1.396 -//       std::queue<Node> queue;
   1.397 -
   1.398 -//       M.set(t,false);
   1.399 -//       queue.push(t);
   1.400 -
   1.401 -//       while (!queue.empty()) {
   1.402 -//         Node w=queue.front();
   1.403 -// 	queue.pop();
   1.404 -
   1.405 -// 	InEdgeIt e;
   1.406 -// 	for(g->first(e,w) ; g->valid(e); g->next(e)) {
   1.407 -// 	  Node v=g->tail(e);
   1.408 -// 	  if (M[v] && (*flow)[e] < (*capacity)[e] ) {
   1.409 -// 	    queue.push(v);
   1.410 -// 	    M.set(v, false);
   1.411 -// 	  }
   1.412 -// 	}
   1.413 -
   1.414 -// 	OutEdgeIt f;
   1.415 -// 	for(g->first(f,w) ; g->valid(f); g->next(f)) {
   1.416 -// 	  Node v=g->head(f);
   1.417 -// 	  if (M[v] && (*flow)[f] > 0 ) {
   1.418 -// 	    queue.push(v);
   1.419 -// 	    M.set(v, false);
   1.420 -// 	  }
   1.421 -// 	}
   1.422 -//       }
   1.423 -//     }
   1.424 -
   1.425 -//     ///Returns a minimum value cut.
   1.426 -
   1.427 -//     ///Sets \c M to the characteristic vector of a minimum value cut.
   1.428 -//     ///\pre M should be a node map of bools initialized to false.
   1.429 -//     ///\pre \c flow must be a maximum flow.    
   1.430 -//     template<typename CutMap>
   1.431 -//     void minCut(CutMap& M) const { minMinCut(M); }
   1.432 -
   1.433 -//     ///Resets the source node to \c _s.
   1.434 -
   1.435 -//     ///Resets the source node to \c _s.
   1.436 -//     /// 
   1.437 -//     void resetSource(Node _s) { s=_s; status=AFTER_NOTHING; }
   1.438 -
   1.439 -//     ///Resets the target node to \c _t.
   1.440 -
   1.441 -//     ///Resets the target node to \c _t.
   1.442 -//     ///
   1.443 -//     void resetTarget(Node _t) { t=_t; status=AFTER_NOTHING; }
   1.444 -
   1.445 -//     /// Resets the edge map of the capacities to _cap.
   1.446 -
   1.447 -//     /// Resets the edge map of the capacities to _cap.
   1.448 -//     /// 
   1.449 -//     void resetCap(const CapMap& _cap) { capacity=&_cap; status=AFTER_NOTHING; }
   1.450 -
   1.451 -//     /// Resets the edge map of the flows to _flow.
   1.452 -
   1.453 -//     /// Resets the edge map of the flows to _flow.
   1.454 -//     /// 
   1.455 -//     void resetFlow(FlowMap& _flow) { flow=&_flow; status=AFTER_NOTHING; }
   1.456 -
   1.457 -
   1.458 -//   private:
   1.459 -
   1.460 -//     int push(Node w, VecStack& active) {
   1.461 -
   1.462 -//       int lev=level[w];
   1.463 -//       Num exc=excess[w];
   1.464 -//       int newlevel=n;       //bound on the next level of w
   1.465 -
   1.466 -//       OutEdgeIt e;
   1.467 -//       for(g->first(e,w); g->valid(e); g->next(e)) {
   1.468 -
   1.469 -// 	if ( (*flow)[e] >= (*capacity)[e] ) continue;
   1.470 -// 	Node v=g->head(e);
   1.471 -
   1.472 -// 	if( lev > level[v] ) { //Push is allowed now
   1.473 -
   1.474 -// 	  if ( excess[v]<=0 && v!=t && v!=s ) {
   1.475 -// 	    int lev_v=level[v];
   1.476 -// 	    active[lev_v].push(v);
   1.477 -// 	  }
   1.478 -
   1.479 -// 	  Num cap=(*capacity)[e];
   1.480 -// 	  Num flo=(*flow)[e];
   1.481 -// 	  Num remcap=cap-flo;
   1.482 -
   1.483 -// 	  if ( remcap >= exc ) { //A nonsaturating push.
   1.484 -
   1.485 -// 	    flow->set(e, flo+exc);
   1.486 -// 	    excess.set(v, excess[v]+exc);
   1.487 -// 	    exc=0;
   1.488 -// 	    break;
   1.489 -
   1.490 -// 	  } else { //A saturating push.
   1.491 -// 	    flow->set(e, cap);
   1.492 -// 	    excess.set(v, excess[v]+remcap);
   1.493 -// 	    exc-=remcap;
   1.494 -// 	  }
   1.495 -// 	} else if ( newlevel > level[v] ) newlevel = level[v];
   1.496 -//       } //for out edges wv
   1.497 -
   1.498 -//       if ( exc > 0 ) {
   1.499 -// 	InEdgeIt e;
   1.500 -// 	for(g->first(e,w); g->valid(e); g->next(e)) {
   1.501 -
   1.502 -// 	  if( (*flow)[e] <= 0 ) continue;
   1.503 -// 	  Node v=g->tail(e);
   1.504 -
   1.505 -// 	  if( lev > level[v] ) { //Push is allowed now
   1.506 -
   1.507 -// 	    if ( excess[v]<=0 && v!=t && v!=s ) {
   1.508 -// 	      int lev_v=level[v];
   1.509 -// 	      active[lev_v].push(v);
   1.510 -// 	    }
   1.511 -
   1.512 -// 	    Num flo=(*flow)[e];
   1.513 -
   1.514 -// 	    if ( flo >= exc ) { //A nonsaturating push.
   1.515 -
   1.516 -// 	      flow->set(e, flo-exc);
   1.517 -// 	      excess.set(v, excess[v]+exc);
   1.518 -// 	      exc=0;
   1.519 -// 	      break;
   1.520 -// 	    } else {  //A saturating push.
   1.521 -
   1.522 -// 	      excess.set(v, excess[v]+flo);
   1.523 -// 	      exc-=flo;
   1.524 -// 	      flow->set(e,0);
   1.525 -// 	    }
   1.526 -// 	  } else if ( newlevel > level[v] ) newlevel = level[v];
   1.527 -// 	} //for in edges vw
   1.528 -
   1.529 -//       } // if w still has excess after the out edge for cycle
   1.530 -
   1.531 -//       excess.set(w, exc);
   1.532 -
   1.533 -//       return newlevel;
   1.534 -//     }
   1.535 -
   1.536 -
   1.537 -//     void preflowPreproc(FlowEnum fe, VecStack& active,
   1.538 -// 			VecNode& level_list, NNMap& left, NNMap& right)
   1.539 -//     {
   1.540 -//       std::queue<Node> bfs_queue;
   1.541 -
   1.542 -//       switch (fe) {
   1.543 -//       case NO_FLOW:   //flow is already set to const zero in this case
   1.544 -//       case ZERO_FLOW:
   1.545 -// 	{
   1.546 -// 	  //Reverse_bfs from t, to find the starting level.
   1.547 -// 	  level.set(t,0);
   1.548 -// 	  bfs_queue.push(t);
   1.549 -
   1.550 -// 	  while (!bfs_queue.empty()) {
   1.551 -
   1.552 -// 	    Node v=bfs_queue.front();
   1.553 -// 	    bfs_queue.pop();
   1.554 -// 	    int l=level[v]+1;
   1.555 -
   1.556 -// 	    InEdgeIt e;
   1.557 -// 	    for(g->first(e,v); g->valid(e); g->next(e)) {
   1.558 -// 	      Node w=g->tail(e);
   1.559 -// 	      if ( level[w] == n && w != s ) {
   1.560 -// 		bfs_queue.push(w);
   1.561 -// 		Node first=level_list[l];
   1.562 -// 		if ( g->valid(first) ) left.set(first,w);
   1.563 -// 		right.set(w,first);
   1.564 -// 		level_list[l]=w;
   1.565 -// 		level.set(w, l);
   1.566 -// 	      }
   1.567 -// 	    }
   1.568 -// 	  }
   1.569 -
   1.570 -// 	  //the starting flow
   1.571 -// 	  OutEdgeIt e;
   1.572 -// 	  for(g->first(e,s); g->valid(e); g->next(e))
   1.573 -// 	    {
   1.574 -// 	      Num c=(*capacity)[e];
   1.575 -// 	      if ( c <= 0 ) continue;
   1.576 -// 	      Node w=g->head(e);
   1.577 -// 	      if ( level[w] < n ) {
   1.578 -// 		if ( excess[w] <= 0 && w!=t ) active[level[w]].push(w);
   1.579 -// 		flow->set(e, c);
   1.580 -// 		excess.set(w, excess[w]+c);
   1.581 -// 	      }
   1.582 -// 	    }
   1.583 -// 	  break;
   1.584 -// 	}
   1.585 -
   1.586 -//       case GEN_FLOW:
   1.587 -//       case PRE_FLOW:
   1.588 -// 	{
   1.589 -// 	  //Reverse_bfs from t in the residual graph,
   1.590 -// 	  //to find the starting level.
   1.591 -// 	  level.set(t,0);
   1.592 -// 	  bfs_queue.push(t);
   1.593 -
   1.594 -// 	  while (!bfs_queue.empty()) {
   1.595 -
   1.596 -// 	    Node v=bfs_queue.front();
   1.597 -// 	    bfs_queue.pop();
   1.598 -// 	    int l=level[v]+1;
   1.599 -
   1.600 -// 	    InEdgeIt e;
   1.601 -// 	    for(g->first(e,v); g->valid(e); g->next(e)) {
   1.602 -// 	      if ( (*capacity)[e] <= (*flow)[e] ) continue;
   1.603 -// 	      Node w=g->tail(e);
   1.604 -// 	      if ( level[w] == n && w != s ) {
   1.605 -// 		bfs_queue.push(w);
   1.606 -// 		Node first=level_list[l];
   1.607 -// 		if ( g->valid(first) ) left.set(first,w);
   1.608 -// 		right.set(w,first);
   1.609 -// 		level_list[l]=w;
   1.610 -// 		level.set(w, l);
   1.611 -// 	      }
   1.612 -// 	    }
   1.613 -
   1.614 -// 	    OutEdgeIt f;
   1.615 -// 	    for(g->first(f,v); g->valid(f); g->next(f)) {
   1.616 -// 	      if ( 0 >= (*flow)[f] ) continue;
   1.617 -// 	      Node w=g->head(f);
   1.618 -// 	      if ( level[w] == n && w != s ) {
   1.619 -// 		bfs_queue.push(w);
   1.620 -// 		Node first=level_list[l];
   1.621 -// 		if ( g->valid(first) ) left.set(first,w);
   1.622 -// 		right.set(w,first);
   1.623 -// 		level_list[l]=w;
   1.624 -// 		level.set(w, l);
   1.625 -// 	      }
   1.626 -// 	    }
   1.627 -// 	  }
   1.628 -
   1.629 -
   1.630 -// 	  //the starting flow
   1.631 -// 	  OutEdgeIt e;
   1.632 -// 	  for(g->first(e,s); g->valid(e); g->next(e))
   1.633 -// 	    {
   1.634 -// 	      Num rem=(*capacity)[e]-(*flow)[e];
   1.635 -// 	      if ( rem <= 0 ) continue;
   1.636 -// 	      Node w=g->head(e);
   1.637 -// 	      if ( level[w] < n ) {
   1.638 -// 		if ( excess[w] <= 0 && w!=t ) active[level[w]].push(w);
   1.639 -// 		flow->set(e, (*capacity)[e]);
   1.640 -// 		excess.set(w, excess[w]+rem);
   1.641 -// 	      }
   1.642 -// 	    }
   1.643 -
   1.644 -// 	  InEdgeIt f;
   1.645 -// 	  for(g->first(f,s); g->valid(f); g->next(f))
   1.646 -// 	    {
   1.647 -// 	      if ( (*flow)[f] <= 0 ) continue;
   1.648 -// 	      Node w=g->tail(f);
   1.649 -// 	      if ( level[w] < n ) {
   1.650 -// 		if ( excess[w] <= 0 && w!=t ) active[level[w]].push(w);
   1.651 -// 		excess.set(w, excess[w]+(*flow)[f]);
   1.652 -// 		flow->set(f, 0);
   1.653 -// 	      }
   1.654 -// 	    }
   1.655 -// 	  break;
   1.656 -// 	} //case PRE_FLOW
   1.657 -//       }
   1.658 -//     } //preflowPreproc
   1.659 -
   1.660 -
   1.661 -
   1.662 -//     void relabel(Node w, int newlevel, VecStack& active,
   1.663 -// 		 VecNode& level_list, NNMap& left,
   1.664 -// 		 NNMap& right, int& b, int& k, bool what_heur )
   1.665 -//     {
   1.666 -
   1.667 -//       //FIXME jacint: ez mitol num
   1.668 -// //      Num lev=level[w];
   1.669 -//       int lev=level[w];
   1.670 -
   1.671 -//       Node right_n=right[w];
   1.672 -//       Node left_n=left[w];
   1.673 -
   1.674 -//       //unlacing starts
   1.675 -//       if ( g->valid(right_n) ) {
   1.676 -// 	if ( g->valid(left_n) ) {
   1.677 -// 	  right.set(left_n, right_n);
   1.678 -// 	  left.set(right_n, left_n);
   1.679 -// 	} else {
   1.680 -// 	  level_list[lev]=right_n;
   1.681 -// 	  left.set(right_n, INVALID);
   1.682 -// 	}
   1.683 -//       } else {
   1.684 -// 	if ( g->valid(left_n) ) {
   1.685 -// 	  right.set(left_n, INVALID);
   1.686 -// 	} else {
   1.687 -// 	  level_list[lev]=INVALID;
   1.688 -// 	}
   1.689 -//       }
   1.690 -//       //unlacing ends
   1.691 -
   1.692 -//       if ( !g->valid(level_list[lev]) ) {
   1.693 -
   1.694 -// 	//gapping starts
   1.695 -// 	for (int i=lev; i!=k ; ) {
   1.696 -// 	  Node v=level_list[++i];
   1.697 -// 	  while ( g->valid(v) ) {
   1.698 -// 	    level.set(v,n);
   1.699 -// 	    v=right[v];
   1.700 -// 	  }
   1.701 -// 	  level_list[i]=INVALID;
   1.702 -// 	  if ( !what_heur ) {
   1.703 -// 	    while ( !active[i].empty() ) {
   1.704 -// 	      active[i].pop();    //FIXME: ezt szebben kene
   1.705 -// 	    }
   1.706 -// 	  }
   1.707 -// 	}
   1.708 -
   1.709 -// 	level.set(w,n);
   1.710 -// 	b=lev-1;
   1.711 -// 	k=b;
   1.712 -// 	//gapping ends
   1.713 -
   1.714 -//       } else {
   1.715 -
   1.716 -// 	if ( newlevel == n ) level.set(w,n);
   1.717 -// 	else {
   1.718 -// 	  level.set(w,++newlevel);
   1.719 -// 	  active[newlevel].push(w);
   1.720 -// 	  if ( what_heur ) b=newlevel;
   1.721 -// 	  if ( k < newlevel ) ++k;      //now k=newlevel
   1.722 -// 	  Node first=level_list[newlevel];
   1.723 -// 	  if ( g->valid(first) ) left.set(first,w);
   1.724 -// 	  right.set(w,first);
   1.725 -// 	  left.set(w,INVALID);
   1.726 -// 	  level_list[newlevel]=w;
   1.727 -// 	}
   1.728 -//       }
   1.729 -
   1.730 -//     } //relabel
   1.731 -
   1.732 -//   };
   1.733 -
   1.734 -
   1.735 -
   1.736 -//   template <typename Graph, typename Num, typename CapMap, typename FlowMap>
   1.737 -//   void MaxFlow<Graph, Num, CapMap, FlowMap>::preflowPhase1(FlowEnum fe)
   1.738 -//   {
   1.739 -
   1.740 -//     int heur0=(int)(H0*n);  //time while running 'bound decrease'
   1.741 -//     int heur1=(int)(H1*n);  //time while running 'highest label'
   1.742 -//     int heur=heur1;         //starting time interval (#of relabels)
   1.743 -//     int numrelabel=0;
   1.744 -
   1.745 -//     bool what_heur=1;
   1.746 -//     //It is 0 in case 'bound decrease' and 1 in case 'highest label'
   1.747 -
   1.748 -//     bool end=false;
   1.749 -//     //Needed for 'bound decrease', true means no active nodes are above bound
   1.750 -//     //b.
   1.751 -
   1.752 -//     int k=n-2;  //bound on the highest level under n containing a node
   1.753 -//     int b=k;    //bound on the highest level under n of an active node
   1.754 -
   1.755 -//     VecStack active(n);
   1.756 -
   1.757 -//     NNMap left(*g, INVALID);
   1.758 -//     NNMap right(*g, INVALID);
   1.759 -//     VecNode level_list(n,INVALID);
   1.760 -//     //List of the nodes in level i<n, set to n.
   1.761 -
   1.762 -//     NodeIt v;
   1.763 -//     for(g->first(v); g->valid(v); g->next(v)) level.set(v,n);
   1.764 -//     //setting each node to level n
   1.765 -
   1.766 -//     if ( fe == NO_FLOW ) {
   1.767 -//       EdgeIt e;
   1.768 -//       for(g->first(e); g->valid(e); g->next(e)) flow->set(e,0);
   1.769 -//     }
   1.770 -
   1.771 -//     switch (fe) { //computing the excess
   1.772 -//     case PRE_FLOW:
   1.773 -//       {
   1.774 -// 	NodeIt v;
   1.775 -// 	for(g->first(v); g->valid(v); g->next(v)) {
   1.776 -// 	  Num exc=0;
   1.777 -
   1.778 -// 	  InEdgeIt e;
   1.779 -// 	  for(g->first(e,v); g->valid(e); g->next(e)) exc+=(*flow)[e];
   1.780 -// 	  OutEdgeIt f;
   1.781 -// 	  for(g->first(f,v); g->valid(f); g->next(f)) exc-=(*flow)[f];
   1.782 -
   1.783 -// 	  excess.set(v,exc);
   1.784 -
   1.785 -// 	  //putting the active nodes into the stack
   1.786 -// 	  int lev=level[v];
   1.787 -// 	  if ( exc > 0 && lev < n && v != t ) active[lev].push(v);
   1.788 -// 	}
   1.789 -// 	break;
   1.790 -//       }
   1.791 -//     case GEN_FLOW:
   1.792 -//       {
   1.793 -// 	NodeIt v;
   1.794 -// 	for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0);
   1.795 -
   1.796 -// 	Num exc=0;
   1.797 -// 	InEdgeIt e;
   1.798 -// 	for(g->first(e,t); g->valid(e); g->next(e)) exc+=(*flow)[e];
   1.799 -// 	OutEdgeIt f;
   1.800 -// 	for(g->first(f,t); g->valid(f); g->next(f)) exc-=(*flow)[f];
   1.801 -// 	excess.set(t,exc);
   1.802 -// 	break;
   1.803 -//       }
   1.804 -//     case ZERO_FLOW:
   1.805 -//     case NO_FLOW:
   1.806 -//       {
   1.807 -// 	NodeIt v;
   1.808 -//         for(g->first(v); g->valid(v); g->next(v)) excess.set(v,0);
   1.809 -// 	break;
   1.810 -//       }
   1.811 -//     }
   1.812 -
   1.813 -//     preflowPreproc(fe, active, level_list, left, right);
   1.814 -//     //End of preprocessing
   1.815 -
   1.816 -
   1.817 -//     //Push/relabel on the highest level active nodes.
   1.818 -//     while ( true ) {
   1.819 -//       if ( b == 0 ) {
   1.820 -// 	if ( !what_heur && !end && k > 0 ) {
   1.821 -// 	  b=k;
   1.822 -// 	  end=true;
   1.823 -// 	} else break;
   1.824 -//       }
   1.825 -
   1.826 -//       if ( active[b].empty() ) --b;
   1.827 -//       else {
   1.828 -// 	end=false;
   1.829 -// 	Node w=active[b].top();
   1.830 -// 	active[b].pop();
   1.831 -// 	int newlevel=push(w,active);
   1.832 -// 	if ( excess[w] > 0 ) relabel(w, newlevel, active, level_list,
   1.833 -// 				     left, right, b, k, what_heur);
   1.834 -
   1.835 -// 	++numrelabel;
   1.836 -// 	if ( numrelabel >= heur ) {
   1.837 -// 	  numrelabel=0;
   1.838 -// 	  if ( what_heur ) {
   1.839 -// 	    what_heur=0;
   1.840 -// 	    heur=heur0;
   1.841 -// 	    end=false;
   1.842 -// 	  } else {
   1.843 -// 	    what_heur=1;
   1.844 -// 	    heur=heur1;
   1.845 -// 	    b=k;
   1.846 -// 	  }
   1.847 -// 	}
   1.848 -//       }
   1.849 -//     }
   1.850 -
   1.851 -//     status=AFTER_PRE_FLOW_PHASE_1;
   1.852 -//   }
   1.853 -
   1.854 -
   1.855 -
   1.856 -//   template <typename Graph, typename Num, typename CapMap, typename FlowMap>
   1.857 -//   void MaxFlow<Graph, Num, CapMap, FlowMap>::preflowPhase2()
   1.858 -//   {
   1.859 -
   1.860 -//     int k=n-2;  //bound on the highest level under n containing a node
   1.861 -//     int b=k;    //bound on the highest level under n of an active node
   1.862 -
   1.863 -//     VecStack active(n);
   1.864 -//     level.set(s,0);
   1.865 -//     std::queue<Node> bfs_queue;
   1.866 -//     bfs_queue.push(s);
   1.867 -
   1.868 -//     while (!bfs_queue.empty()) {
   1.869 -
   1.870 -//       Node v=bfs_queue.front();
   1.871 -//       bfs_queue.pop();
   1.872 -//       int l=level[v]+1;
   1.873 -
   1.874 -//       InEdgeIt e;
   1.875 -//       for(g->first(e,v); g->valid(e); g->next(e)) {
   1.876 -// 	if ( (*capacity)[e] <= (*flow)[e] ) continue;
   1.877 -// 	Node u=g->tail(e);
   1.878 -// 	if ( level[u] >= n ) {
   1.879 -// 	  bfs_queue.push(u);
   1.880 -// 	  level.set(u, l);
   1.881 -// 	  if ( excess[u] > 0 ) active[l].push(u);
   1.882 -// 	}
   1.883 -//       }
   1.884 -
   1.885 -//       OutEdgeIt f;
   1.886 -//       for(g->first(f,v); g->valid(f); g->next(f)) {
   1.887 -// 	if ( 0 >= (*flow)[f] ) continue;
   1.888 -// 	Node u=g->head(f);
   1.889 -// 	if ( level[u] >= n ) {
   1.890 -// 	  bfs_queue.push(u);
   1.891 -// 	  level.set(u, l);
   1.892 -// 	  if ( excess[u] > 0 ) active[l].push(u);
   1.893 -// 	}
   1.894 -//       }
   1.895 -//     }
   1.896 -//     b=n-2;
   1.897 -
   1.898 -//     while ( true ) {
   1.899 -
   1.900 -//       if ( b == 0 ) break;
   1.901 -
   1.902 -//       if ( active[b].empty() ) --b;
   1.903 -//       else {
   1.904 -// 	Node w=active[b].top();
   1.905 -// 	active[b].pop();
   1.906 -// 	int newlevel=push(w,active);
   1.907 -
   1.908 -// 	//relabel
   1.909 -// 	if ( excess[w] > 0 ) {
   1.910 -// 	  level.set(w,++newlevel);
   1.911 -// 	  active[newlevel].push(w);
   1.912 -// 	  b=newlevel;
   1.913 -// 	}
   1.914 -//       }  // if stack[b] is nonempty
   1.915 -//     } // while(true)
   1.916 -
   1.917 -//     status=AFTER_PRE_FLOW_PHASE_2;
   1.918 -//   }
   1.919 -
   1.920 -
   1.921 +  ///\author Marton Makai
   1.922    template <typename Graph, typename Num,
   1.923  	    typename CapMap=typename Graph::template EdgeMap<Num>,
   1.924              typename FlowMap=typename Graph::template EdgeMap<Num> >
   1.925 @@ -873,10 +91,6 @@
   1.926      typedef typename Graph::OutEdgeIt OutEdgeIt;
   1.927      typedef typename Graph::InEdgeIt InEdgeIt;
   1.928  
   1.929 -//    typedef typename std::vector<std::stack<Node> > VecStack;
   1.930 -//    typedef typename Graph::template NodeMap<Node> NNMap;
   1.931 -//    typedef typename std::vector<Node> VecNode;
   1.932 -
   1.933      const Graph* g;
   1.934      Node s;
   1.935      Node t;
   1.936 @@ -890,37 +104,12 @@
   1.937      //typedef typename ResGW::template NodeMap<bool> ReachedMap;
   1.938      typedef typename Graph::template NodeMap<int> ReachedMap;
   1.939  
   1.940 -
   1.941      //level works as a bool map in augmenting path algorithms and is
   1.942      //used by bfs for storing reached information.  In preflow, it
   1.943      //shows the levels of nodes.     
   1.944      ReachedMap level;
   1.945  
   1.946 -    //excess is needed only in preflow
   1.947 -//    typename Graph::template NodeMap<Num> excess;
   1.948 -
   1.949 -    //fixme    
   1.950 -//   protected:
   1.951 -    //     MaxFlow() { }
   1.952 -    //     void set(const Graph& _G, Node _s, Node _t, const CapMap& _capacity,
   1.953 -    // 	     FlowMap& _flow)
   1.954 -    //       {
   1.955 -    // 	g=&_G;
   1.956 -    // 	s=_s;
   1.957 -    // 	t=_t;
   1.958 -    // 	capacity=&_capacity;
   1.959 -    // 	flow=&_flow;
   1.960 -    // 	n=_G.nodeNum;
   1.961 -    // 	level.set (_G); //kellene vmi ilyesmi fv
   1.962 -    // 	excess(_G,0); //itt is
   1.963 -    //       }
   1.964 -
   1.965 -    // constants used for heuristics
   1.966 -//    static const int H0=20;
   1.967 -//    static const int H1=1;
   1.968 -
   1.969    public:
   1.970 -
   1.971      ///Indicates the property of the starting flow.
   1.972  
   1.973      ///Indicates the property of the starting flow. The meanings are as follows:
   1.974 @@ -1088,28 +277,6 @@
   1.975        //marci figyu: excess[t] epp ezt adja preflow 1. fazisa utan   
   1.976      }
   1.977  
   1.978 -    template<typename MapGraphWrapper>
   1.979 -    class DistanceMap {
   1.980 -    protected:
   1.981 -      const MapGraphWrapper* g;
   1.982 -      typename MapGraphWrapper::template NodeMap<int> dist;
   1.983 -    public:
   1.984 -      DistanceMap(MapGraphWrapper& _g) : g(&_g), dist(*g, g->nodeNum()) { }
   1.985 -      void set(const typename MapGraphWrapper::Node& n, int a) {
   1.986 -	dist.set(n, a);
   1.987 -      }
   1.988 -      int operator[](const typename MapGraphWrapper::Node& n) const { 
   1.989 -	return dist[n]; 
   1.990 -      }
   1.991 -      //       int get(const typename MapGraphWrapper::Node& n) const {
   1.992 -      // 	return dist[n]; }
   1.993 -      //       bool get(const typename MapGraphWrapper::Edge& e) const {
   1.994 -      // 	return (dist.get(g->tail(e))<dist.get(g->head(e))); }
   1.995 -      bool operator[](const typename MapGraphWrapper::Edge& e) const {
   1.996 -	return (dist[g->tail(e)]<dist[g->head(e)]);
   1.997 -      }
   1.998 -    };
   1.999 -
  1.1000    };
  1.1001  
  1.1002  
  1.1003 @@ -1244,9 +411,8 @@
  1.1004        res_graph_to_F(res_graph);
  1.1005      {
  1.1006        typename ResGW::NodeIt n;
  1.1007 -      for(res_graph.first(n); n!=INVALID; ++n) {
  1.1008 +      for(res_graph.first(n); n!=INVALID; ++n) 
  1.1009  	res_graph_to_F.set(n, F.addNode());
  1.1010 -      }
  1.1011      }
  1.1012  
  1.1013      typename MG::Node sF=res_graph_to_F[s];
  1.1014 @@ -1336,7 +502,8 @@
  1.1015      return _augment;
  1.1016    }
  1.1017  
  1.1018 -
  1.1019 +  /// Blocking flow augmentation without constructing the layered 
  1.1020 +  /// graph physically in which the blocking flow is computed.
  1.1021    template <typename Graph, typename Num, typename CapMap, typename FlowMap>
  1.1022    bool AugmentingFlow<Graph, Num, CapMap, FlowMap>::augmentOnBlockingFlow2()
  1.1023    {
  1.1024 @@ -1344,37 +511,41 @@
  1.1025  
  1.1026      ResGW res_graph(*g, *capacity, *flow);
  1.1027  
  1.1028 -    //ReachedMap level(res_graph);
  1.1029 +    //Potential map, for distances from s
  1.1030 +    typename ResGW::template NodeMap<int> potential(res_graph, 0);
  1.1031 +    typedef ConstMap<typename ResGW::Edge, int> Const1Map; 
  1.1032 +    Const1Map const_1_map(1);
  1.1033 +    TightEdgeFilterMap<ResGW, typename ResGW::template NodeMap<int>,
  1.1034 +      Const1Map> tight_edge_filter(res_graph, potential, const_1_map);
  1.1035 +
  1.1036      for (typename Graph::NodeIt n(*g); n!=INVALID; ++n) level.set(n, 0);
  1.1037      BfsIterator<ResGW, ReachedMap> bfs(res_graph, level);
  1.1038 +    bfs.pushAndSetReached(s);
  1.1039  
  1.1040 -    bfs.pushAndSetReached(s);
  1.1041 -    DistanceMap<ResGW> dist(res_graph);
  1.1042 +    //computing distances from s in the residual graph
  1.1043      while ( !bfs.finished() ) {
  1.1044        ResGWEdge e=bfs;
  1.1045 -      if (e!=INVALID && bfs.isBNodeNewlyReached()) {
  1.1046 -	dist.set(res_graph.head(e), dist[res_graph.tail(e)]+1);
  1.1047 -      }
  1.1048 +      if (e!=INVALID && bfs.isBNodeNewlyReached())
  1.1049 +	potential.set(res_graph.head(e), potential[res_graph.tail(e)]+1);
  1.1050        ++bfs;
  1.1051 -    } //computing distances from s in the residual graph
  1.1052 +    } 
  1.1053  
  1.1054 -    //Subgraph containing the edges on some shortest paths
  1.1055 +    //Subgraph containing the edges on some shortest paths 
  1.1056 +    //(i.e. tight edges)
  1.1057      ConstMap<typename ResGW::Node, bool> true_map(true);
  1.1058      typedef SubGraphWrapper<ResGW, ConstMap<typename ResGW::Node, bool>,
  1.1059 -      DistanceMap<ResGW> > FilterResGW;
  1.1060 -    FilterResGW filter_res_graph(res_graph, true_map, dist);
  1.1061 +      TightEdgeFilterMap<ResGW, typename ResGW::template NodeMap<int>, 
  1.1062 +      Const1Map> > FilterResGW;
  1.1063 +    FilterResGW filter_res_graph(res_graph, true_map, tight_edge_filter);
  1.1064  
  1.1065      //Subgraph, which is able to delete edges which are already
  1.1066      //met by the dfs
  1.1067      typename FilterResGW::template NodeMap<typename FilterResGW::Edge>
  1.1068        first_out_edges(filter_res_graph);
  1.1069 -    typename FilterResGW::NodeIt v;
  1.1070 -    for(filter_res_graph.first(v); v!=INVALID; ++v)
  1.1071 -      {
  1.1072 -  	typename FilterResGW::OutEdgeIt e;
  1.1073 -  	filter_res_graph.first(e, v);
  1.1074 -  	first_out_edges.set(v, e);
  1.1075 -      }
  1.1076 +    for (typename FilterResGW::NodeIt v(filter_res_graph); v!=INVALID; ++v)
  1.1077 +      first_out_edges.set
  1.1078 +	(v, typename FilterResGW::OutEdgeIt(filter_res_graph, v));
  1.1079 +
  1.1080      typedef ErasingFirstGraphWrapper<FilterResGW, typename FilterResGW::
  1.1081        template NodeMap<typename FilterResGW::Edge> > ErasingResGW;
  1.1082      ErasingResGW erasing_res_graph(filter_res_graph, first_out_edges);
  1.1083 @@ -1407,47 +578,37 @@
  1.1084  	
  1.1085        while (!dfs.finished()) {
  1.1086  	++dfs;
  1.1087 -	if (typename ErasingResGW::Edge(dfs)!=INVALID)
  1.1088 - 	  {
  1.1089 -  	    if (dfs.isBNodeNewlyReached()) {
  1.1090 +	if (typename ErasingResGW::Edge(dfs)!=INVALID) {
  1.1091 +	  if (dfs.isBNodeNewlyReached()) {
  1.1092 +	    
  1.1093 +	    typename ErasingResGW::Node v=erasing_res_graph.tail(dfs);
  1.1094 +	    typename ErasingResGW::Node w=erasing_res_graph.head(dfs);
  1.1095  
  1.1096 - 	      typename ErasingResGW::Node v=erasing_res_graph.tail(dfs);
  1.1097 - 	      typename ErasingResGW::Node w=erasing_res_graph.head(dfs);
  1.1098 +	    pred.set(w, typename ErasingResGW::Edge(dfs));
  1.1099 +	    if (pred[v]!=INVALID) {
  1.1100 +	      free1.set
  1.1101 +		(w, std::min(free1[v], res_graph.resCap
  1.1102 +			     (typename ErasingResGW::Edge(dfs))));
  1.1103 +	    } else {
  1.1104 +	      free1.set
  1.1105 +		(w, res_graph.resCap
  1.1106 +		 (typename ErasingResGW::Edge(dfs)));
  1.1107 +	    }
  1.1108  
  1.1109 - 	      pred.set(w, typename ErasingResGW::Edge(dfs));
  1.1110 - 	      if (pred[v]!=INVALID) {
  1.1111 - 		free1.set
  1.1112 -		  (w, std::min(free1[v], res_graph.resCap
  1.1113 -			       (typename ErasingResGW::Edge(dfs))));
  1.1114 - 	      } else {
  1.1115 - 		free1.set
  1.1116 -		  (w, res_graph.resCap
  1.1117 -		   (typename ErasingResGW::Edge(dfs)));
  1.1118 - 	      }
  1.1119 -
  1.1120 - 	      if (w==t) {
  1.1121 - 		__augment=true;
  1.1122 - 		_augment=true;
  1.1123 - 		break;
  1.1124 - 	      }
  1.1125 - 	    } else {
  1.1126 - 	      erasing_res_graph.erase(dfs);
  1.1127 +	    if (w==t) {
  1.1128 +	      __augment=true;
  1.1129 +	      _augment=true;
  1.1130 +	      break;
  1.1131  	    }
  1.1132 +	  } else {
  1.1133 +	    erasing_res_graph.erase(dfs);
  1.1134  	  }
  1.1135 +	}
  1.1136        }
  1.1137  
  1.1138        if (__augment) {
  1.1139  	typename ErasingResGW::Node
  1.1140  	  n=typename FilterResGW::Node(typename ResGW::Node(t));
  1.1141 -	// 	  typename ResGW::NodeMap<Num> a(res_graph);
  1.1142 -	// 	  typename ResGW::Node b;
  1.1143 -	// 	  Num j=a[b];
  1.1144 -	// 	  typename FilterResGW::NodeMap<Num> a1(filter_res_graph);
  1.1145 -	// 	  typename FilterResGW::Node b1;
  1.1146 -	// 	  Num j1=a1[b1];
  1.1147 -	// 	  typename ErasingResGW::NodeMap<Num> a2(erasing_res_graph);
  1.1148 -	// 	  typename ErasingResGW::Node b2;
  1.1149 -	// 	  Num j2=a2[b2];
  1.1150  	Num augment_value=free1[n];
  1.1151  	while (pred[n]!=INVALID) {
  1.1152  	  typename ErasingResGW::Edge e=pred[n];
  1.1153 @@ -1470,5 +631,3 @@
  1.1154  #endif //HUGO_AUGMENTING_FLOW_H
  1.1155  
  1.1156  
  1.1157 -
  1.1158 -