# HG changeset patch
# User athos
# Date 1084470143 0
# Node ID 933f593824c24740078a3aed562fb762995c916f
# Parent aacabcd724f0d99de5fa3858495e8aa02d2ab384
Started mincostflow.
diff -r aacabcd724f0 -r 933f593824c2 src/work/athos/mincostflow.h
--- a/src/work/athos/mincostflow.h Thu May 13 17:33:40 2004 +0000
+++ b/src/work/athos/mincostflow.h Thu May 13 17:42:23 2004 +0000
@@ -36,13 +36,13 @@
/// where \c M is the value of the maximal flow.
///
///\author Attila Bernath
- template <typename Graph, typename LengthMap, typename SupplyMap>
+ template <typename Graph, typename LengthMap, typename SupplyDemandMap>
class MinCostFlow {
typedef typename LengthMap::ValueType Length;
- typedef typename SupplyMap::ValueType Supply;
+ typedef typename SupplyDemandMap::ValueType SupplyDemand;
typedef typename Graph::Node Node;
typedef typename Graph::NodeIt NodeIt;
@@ -84,7 +84,7 @@
//Input
const Graph& G;
const LengthMap& length;
- const SupplyMap& supply;//supply or demand of nodes
+ const SupplyDemandMap& supply_demand;//supply or demand of nodes
//auxiliary variables
@@ -94,7 +94,7 @@
//To store the potentila (dual variables)
typename Graph::template NodeMap<Length> potential;
//To store excess-deficit values
- SupplyMap excess;
+ SupplyDemandMap excess_deficit;
Length total_length;
@@ -103,39 +103,62 @@
public :
- MinCostFlow(Graph& _G, LengthMap& _length, SupplyMap& _supply) : G(_G),
- length(_length), supply(_supply), flow(_G), potential(_G){ }
+ MinCostFlow(Graph& _G, LengthMap& _length, SupplyDemandMap& _supply_demand) : G(_G),
+ length(_length), supply_demand(_supply_demand), flow(_G), potential(_G){ }
///Runs the algorithm.
///Runs the algorithm.
- ///Returns k if there are at least k edge-disjoint paths from s to t.
- ///Otherwise it returns the number of found edge-disjoint paths from s to t.
+
///\todo May be it does make sense to be able to start with a nonzero
/// feasible primal-dual solution pair as well.
int run() {
//Resetting variables from previous runs
- total_length = 0;
+ //total_length = 0;
+
+ typedef typename Graph::template NodeMap<int> HeapMap;
+ typedef Heap<Node, SupplyDemand, typename Graph::template NodeMap<int>,
+ std::greater<SupplyDemand> > HeapType;
+
+ //A heap for the excess nodes
+ HeapMap excess_nodes_map(G,-1);
+ HeapType excess_nodes(excess_nodes_map);
+
+ //A heap for the deficit nodes
+ HeapMap deficit_nodes_map(G,-1);
+ HeapType deficit_nodes(deficit_nodes_map);
+
FOR_EACH_LOC(typename Graph::EdgeIt, e, G){
flow.set(e,0);
}
//Initial value for delta
- Supply delta = 0;
-
+ SupplyDemand delta = 0;
+
FOR_EACH_LOC(typename Graph::NodeIt, n, G){
- if (delta < abs(supply[e])){
- delta = abs(supply[e]);
+ excess_deficit.set(n,supply_demand[n]);
+ //A supply node
+ if (excess_deficit[n] > 0){
+ excess_nodes.push(n,excess_deficit[n]);
}
- excess.set(n,supply[e]);
+ //A demand node
+ if (excess_deficit[n] < 0){
+ deficit_nodes.push(n, - excess_deficit[n]);
+ }
+ //Finding out starting value of delta
+ if (delta < abs(excess_deficit[n])){
+ delta = abs(excess_deficit[n]);
+ }
//Initialize the copy of the Dijkstra potential to zero
potential.set(n,0);
}
-
+ //It'll be allright as an initial value, though this value
+ //can be the maximum deficit here
+ SupplyDemand max_excess = delta;
//We need a residual graph which is uncapacitated
ResGraphType res_graph(G, flow);
@@ -147,45 +170,102 @@
Dijkstra<ResGraphType, ModLengthMap> dijkstra(res_graph, mod_length);
- int i;
- for (i=0; i<k; ++i){
- dijkstra.run(s);
- if (!dijkstra.reached(t)){
- //There are no k paths from s to t
- break;
- };
+ while (max_excess > 0){
+
- //We have to copy the potential
- FOR_EACH_LOC(typename ResGraphType::NodeIt, n, res_graph){
- potential[n] += dijkstra.distMap()[n];
+ //Merge and stuff
+
+ Node s = excess_nodes.top();
+ SupplyDemand max_excess = excess_nodes[s];
+ Node t = deficit_nodes.top();
+ if (max_excess < dificit_nodes[t]){
+ max_excess = dificit_nodes[t];
+ }
+
+
+ while(max_excess > ){
+
+
+ //s es t valasztasa
+
+ //Dijkstra part
+ dijkstra.run(s);
+
+ /*We know from theory that t can be reached
+ if (!dijkstra.reached(t)){
+ //There are no k paths from s to t
+ break;
+ };
+ */
+
+ //We have to change the potential
+ FOR_EACH_LOC(typename ResGraphType::NodeIt, n, res_graph){
+ potential[n] += dijkstra.distMap()[n];
+ }
+
+
+ //Augmenting on the sortest path
+ Node n=t;
+ ResGraphEdge e;
+ while (n!=s){
+ e = dijkstra.pred(n);
+ n = dijkstra.predNode(n);
+ res_graph.augment(e,delta);
+ /*
+ //Let's update the total length
+ if (res_graph.forward(e))
+ total_length += length[e];
+ else
+ total_length -= length[e];
+ */
+ }
+
+ //Update the excess_nodes heap
+ if (delta >= excess_nodes[s]){
+ if (delta > excess_nodes[s])
+ deficit_nodes.push(s,delta - excess_nodes[s]);
+ excess_nodes.pop();
+
+ }
+ else{
+ excess_nodes[s] -= delta;
+ }
+ //Update the deficit_nodes heap
+ if (delta >= deficit_nodes[t]){
+ if (delta > deficit_nodes[t])
+ excess_nodes.push(t,delta - deficit_nodes[t]);
+ deficit_nodes.pop();
+
+ }
+ else{
+ deficit_nodes[t] -= delta;
+ }
+ //Dijkstra part ends here
}
/*
- {
+ * End of the delta scaling phase
+ */
- typename ResGraphType::NodeIt n;
- for ( res_graph.first(n) ; res_graph.valid(n) ; res_graph.next(n) ) {
- potential[n] += dijkstra.distMap()[n];
+ //Whatever this means
+ delta = delta / 2;
+
+ /*This is not necessary here
+ //Update the max_excess
+ max_excess = 0;
+ FOR_EACH_LOC(typename Graph::NodeIt, n, G){
+ if (max_excess < excess_deficit[n]){
+ max_excess = excess_deficit[n];
}
}
*/
-
- //Augmenting on the sortest path
- Node n=t;
- ResGraphEdge e;
- while (n!=s){
- e = dijkstra.pred(n);
- n = dijkstra.predNode(n);
- res_graph.augment(e,delta);
- //Let's update the total length
- if (res_graph.forward(e))
- total_length += length[e];
- else
- total_length -= length[e];
+ //Reset delta if still too big
+ if (8*number_of_nodes*max_excess <= delta){
+ delta = max_excess;
+
}
-
- }
+ }//while(max_excess > 0)
return i;