1.1 --- a/src/work/athos/mincostflow.h	Fri May 21 12:40:39 2004 +0000
     1.2 +++ b/src/work/athos/mincostflow.h	Mon May 24 10:43:44 2004 +0000
     1.3 @@ -10,7 +10,9 @@
     1.4  #include <hugo/graph_wrapper.h>
     1.5  #include <hugo/maps.h>
     1.6  #include <vector>
     1.7 +#include <list>
     1.8  #include <for_each_macros.h>
     1.9 +#include <hugo/union_find.h>
    1.10  
    1.11  namespace hugo {
    1.12  
    1.13 @@ -118,7 +120,7 @@
    1.14        //total_length = 0;
    1.15  
    1.16        typedef typename Graph::template NodeMap<int> HeapMap;
    1.17 -      typedef Heap<Node, SupplyDemand, typename Graph::template NodeMap<int>,
    1.18 +      typedef Heap< Node, SupplyDemand, typename Graph::template NodeMap<int>,
    1.19  	std::greater<SupplyDemand> > 	HeapType;
    1.20  
    1.21        //A heap for the excess nodes
    1.22 @@ -129,14 +131,25 @@
    1.23        HeapMap deficit_nodes_map(G,-1);
    1.24        HeapType deficit_nodes(deficit_nodes_map);
    1.25  
    1.26 +      //A container to store nonabundant arcs
    1.27 +      list<Edge> nonabundant_arcs;
    1.28        
    1.29        FOR_EACH_LOC(typename Graph::EdgeIt, e, G){
    1.30  	flow.set(e,0);
    1.31 +	nonabundant_arcs.push_back(e);
    1.32        }
    1.33  
    1.34        //Initial value for delta
    1.35        SupplyDemand delta = 0;
    1.36  
    1.37 +      typedef UnionFindEnum<Node, Graph::template NodeMap> UFE;
    1.38 +
    1.39 +      //A union-find structure to store the abundant components
    1.40 +      UFE::MapType abund_comp_map(G);
    1.41 +      UFE abundant_components(abund_comp_map);
    1.42 +
    1.43 +
    1.44 +
    1.45        FOR_EACH_LOC(typename Graph::NodeIt, n, G){
    1.46         	excess_deficit.set(n,supply_demand[n]);
    1.47  	//A supply node
    1.48 @@ -153,6 +166,8 @@
    1.49  	}
    1.50  	//Initialize the copy of the Dijkstra potential to zero
    1.51  	potential.set(n,0);
    1.52 +	//Every single point is an abundant component initially 
    1.53 +	abundant_components.insert(n);
    1.54        }
    1.55  
    1.56        //It'll be allright as an initial value, though this value 
    1.57 @@ -171,11 +186,38 @@
    1.58  
    1.59        while (max_excess > 0){
    1.60  
    1.61 +	//Reset delta if still too big
    1.62 +	if (8*number_of_nodes*max_excess <= delta){
    1.63 +	  delta = max_excess;
    1.64 +	  
    1.65 +	}
    1.66 +
    1.67  	/*
    1.68  	 * Beginning of the delta scaling phase 
    1.69  	*/
    1.70 -	
    1.71  	//Merge and stuff
    1.72 +	{
    1.73 +	  SupplyDemand buf=8*number_of_nodes*delta;
    1.74 +	  list<Edge>::iterator i = nonabundant_arcs.begin();
    1.75 +	  while ( i != nonabundant_arcs.end() ){
    1.76 +	    if (flow[i]>=buf){
    1.77 +	      Node a = abundant_components.find(res_graph.head(i));
    1.78 +	      Node b = abundant_components.find(res_graph.tail(i));
    1.79 +	      //Merge
    1.80 +	      if (a != b){
    1.81 +		//Find path and augment
    1.82 +		//!!!
    1.83 +		//Find path and augment
    1.84 +		abundant_components.join(a,b);
    1.85 +	      }
    1.86 +	      //What happens to i?
    1.87 +	      nonabundant_arcs.erase(i);
    1.88 +	    }
    1.89 +	    else
    1.90 +	      ++i;
    1.91 +	  }
    1.92 +	}
    1.93 +
    1.94  
    1.95  	Node s = excess_nodes.top(); 
    1.96  	SupplyDemand max_excess = excess_nodes[s];
    1.97 @@ -261,11 +303,7 @@
    1.98  	  }
    1.99  	}
   1.100  	*/
   1.101 -	//Reset delta if still too big
   1.102 -	if (8*number_of_nodes*max_excess <= delta){
   1.103 -	  delta = max_excess;
   1.104 -	  
   1.105 -	}
   1.106 +
   1.107  	  
   1.108        }//while(max_excess > 0)
   1.109