Changeset 659:c5984e925384 in lemon0.x
 Timestamp:
 05/25/04 14:31:18 (20 years ago)
 Branch:
 default
 Phase:
 public
 Convert:
 svn:c9d7d8f590d60310b91f818b3a526b0e/lemon/trunk@862
 Location:
 src/work/athos
 Files:

 1 added
 1 edited
Legend:
 Unmodified
 Added
 Removed

src/work/athos/mincostflow.h
r657 r659 60 60 //typedef typename ResGraphType::template NodeMap<Length> NodeMap; 61 61 typedef typename Graph::template NodeMap<Length> NodeMap; 62 const ResGraphType& G;62 const ResGraphType& res_graph; 63 63 // const EdgeIntMap& rev; 64 64 const LengthMap &ol; … … 69 69 70 70 ValueType operator[](typename ResGraphType::Edge e) const { 71 if ( G.forward(e))72 return ol[e](pot[ G.head(e)]pot[G.tail(e)]);71 if (res_graph.forward(e)) 72 return ol[e](pot[res_graph.head(e)]pot[res_graph.tail(e)]); 73 73 else 74 return ol[e](pot[ G.head(e)]pot[G.tail(e)]);74 return ol[e](pot[res_graph.head(e)]pot[res_graph.tail(e)]); 75 75 } 76 76 77 ModLengthMap(const ResGraphType& _ G,77 ModLengthMap(const ResGraphType& _res_graph, 78 78 const LengthMap &o, const NodeMap &p) : 79 G(_G), /*rev(_rev),*/ ol(o), pot(p){};79 res_graph(_res_graph), /*rev(_rev),*/ ol(o), pot(p){}; 80 80 };//ModLengthMap 81 81 … … 84 84 85 85 //Input 86 const Graph& G;86 const Graph& graph; 87 87 const LengthMap& length; 88 88 const SupplyDemandMap& supply_demand;//supply or demand of nodes … … 105 105 106 106 107 MinCostFlow(Graph& _ G, LengthMap& _length, SupplyDemandMap& _supply_demand) : G(_G),108 length(_length), supply_demand(_supply_demand), flow(_ G), potential(_G){ }107 MinCostFlow(Graph& _graph, LengthMap& _length, SupplyDemandMap& _supply_demand) : graph(_graph), 108 length(_length), supply_demand(_supply_demand), flow(_graph), potential(_graph){ } 109 109 110 110 … … 115 115 ///\todo May be it does make sense to be able to start with a nonzero 116 116 /// feasible primaldual solution pair as well. 117 intrun() {117 void run() { 118 118 119 119 //Resetting variables from previous runs … … 125 125 126 126 //A heap for the excess nodes 127 HeapMap excess_nodes_map( G,1);127 HeapMap excess_nodes_map(graph,1); 128 128 HeapType excess_nodes(excess_nodes_map); 129 129 130 130 //A heap for the deficit nodes 131 HeapMap deficit_nodes_map( G,1);131 HeapMap deficit_nodes_map(graph,1); 132 132 HeapType deficit_nodes(deficit_nodes_map); 133 133 134 134 //A container to store nonabundant arcs 135 135 list<Edge> nonabundant_arcs; 136 137 FOR_EACH_LOC(typename Graph::EdgeIt, e, G){ 136 137 138 FOR_EACH_LOC(typename Graph::EdgeIt, e, graph){ 138 139 flow.set(e,0); 139 140 nonabundant_arcs.push_back(e); … … 146 147 147 148 //A unionfind structure to store the abundant components 148 UFE::MapType abund_comp_map( G);149 UFE::MapType abund_comp_map(graph); 149 150 UFE abundant_components(abund_comp_map); 150 151 151 152 152 153 153 FOR_EACH_LOC(typename Graph::NodeIt, n, G){154 FOR_EACH_LOC(typename Graph::NodeIt, n, graph){ 154 155 excess_deficit.set(n,supply_demand[n]); 155 156 //A supply node … … 175 176 SupplyDemand max_excess = delta; 176 177 178 //ConstMap<Edge,SupplyDemand> ConstEdgeMap; 179 177 180 //We need a residual graph which is uncapacitated 178 ResGraphType res_graph(G, flow); 179 180 181 ResGraphType res_graph(graph, flow); 182 183 //An EdgeMap to tell which arcs are abundant 184 template typename Graph::EdgeMap<bool> abundant_arcs(graph); 185 186 //Let's construct the sugraph consisting only of the abundant edges 187 typedef ConstMap< typename Graph::Node, bool > ConstNodeMap; 188 ConstNodeMap const_true_map(true); 189 typedef SubGraphWrapper< Graph, ConstNodeMap, 190 template typename Graph::EdgeMap<bool> > 191 AbundantGraph; 192 AbundantGraph abundant_graph(graph, const_true_map, abundant_arcs ); 193 194 //Let's construct the residual graph for the abundant graph 195 typedef ResGraphWrapper<const AbundantGraph,int,CapacityMap,EdgeIntMap> 196 ResAbGraph; 197 //Again uncapacitated 198 ResAbGraph res_ab_graph(abundant_graph, flow); 199 200 //We need things for the bfs 201 typename ResAbGraph::NodeMap<bool> bfs_reached(res_ab_graph); 202 typename ResAbGraph::NodeMap<typename ResAbGraph::Edge> 203 bfs_pred(res_ab_graph); 204 NullMap<typename ResAbGraph::Node, int> bfs_dist_dummy(res_ab_graph); 205 //We want to run bfses (more) on this graph 'res_ab_graph' 206 Bfs < ResAbGraph , 207 typename ResAbGraph::NodeMap<bool>, 208 typename ResAbGraph::NodeMap<typename ResAbGraph::Edge>, 209 NullMap<typename ResAbGraph::Node, int> > 210 bfs(res_ab_graph, bfs_reached, bfs_pred, bfs_dist_dummy); 181 211 182 212 ModLengthMap mod_length(res_graph, length, potential); … … 206 236 //Merge 207 237 if (a != b){ 208 //Find path and augment209 //!!!210 //Find path and augment211 238 abundant_components.join(a,b); 239 //We want to push the smaller 240 //Which has greater absolut value excess/deficit 241 Node root=(abs(excess_deficit[a])>abs(excess_deficit[b]))?a:b; 242 //Which is the other 243 Node non_root = ( a == root ) ? b : a ; 244 abundant_components.makeRep(root); 245 SupplyDemand qty_to_augment = abs(excess_deficit[non_root]); 246 //Push the positive value 247 if (excess_deficit[non_root] < 0) 248 swap(root, non_root); 249 //If the non_root node has excess/deficit at all 250 if (qty_to_augment>0){ 251 //Find path and augment 252 bfs.run(non_root); 253 //root should be reached 254 255 //Augmenting on the found path 256 Node n=root; 257 ResGraphEdge e; 258 while (n!=non_root){ 259 e = bfs_pred(n); 260 n = res_graph.tail(e); 261 res_graph.augment(e,qty_to_augment); 262 } 263 264 //We know that non_root had positive excess 265 excess_nodes[non_root] = qty_to_augment; 266 //But what about root node 267 //It might have been positive and so became larger 268 if (excess_deficit[root]>0){ 269 excess_nodes[root] += qty_to_augment; 270 } 271 else{ 272 //Or negative but not turned into positive 273 deficit_nodes[root] = qty_to_augment; 274 } 275 276 //Update the excess_deficit map 277 excess_deficit[non_root] = qty_to_augment; 278 excess_deficit[root] += qty_to_augment; 279 280 281 } 212 282 } 213 283 //What happens to i? 214 nonabundant_arcs.erase(i); 284 //Marci and Zsolt says I shouldn't do such things 285 nonabundant_arcs.erase(i++); 286 abundant_arcs[i] = true; 215 287 } 216 288 else … … 223 295 SupplyDemand max_excess = excess_nodes[s]; 224 296 Node t = deficit_nodes.top(); 225 if (max_excess < d ificit_nodes[t]){226 max_excess = d ificit_nodes[t];227 } 228 229 230 while(max_excess > 0){231 297 if (max_excess < deficit_nodes[t]){ 298 max_excess = deficit_nodes[t]; 299 } 300 301 302 while(max_excess > (n1)*delta/n){ 303 232 304 233 305 //s es t valasztasa 234 306 235 307 //Dijkstra part 236 308 dijkstra.run(s); 237 309 238 310 /*We know from theory that t can be reached 239 311 if (!dijkstra.reached(t)){ … … 264 336 */ 265 337 } 338 339 //Update the excess_deficit map 340 excess_deficit[s] = delta; 341 excess_deficit[t] += delta; 342 266 343 267 344 //Update the excess_nodes heap … … 286 363 } 287 364 //Dijkstra part ends here 365 366 //Choose s and t again 367 s = excess_nodes.top(); 368 max_excess = excess_nodes[s]; 369 t = deficit_nodes.top(); 370 if (max_excess < deficit_nodes[t]){ 371 max_excess = deficit_nodes[t]; 372 } 373 288 374 } 289 375 … … 298 384 //Update the max_excess 299 385 max_excess = 0; 300 FOR_EACH_LOC(typename Graph::NodeIt, n, G){386 FOR_EACH_LOC(typename Graph::NodeIt, n, graph){ 301 387 if (max_excess < excess_deficit[n]){ 302 388 max_excess = excess_deficit[n]; … … 337 423 Length mod_pot; 338 424 Length fl_e; 339 FOR_EACH_LOC(typename Graph::EdgeIt, e, G){425 FOR_EACH_LOC(typename Graph::EdgeIt, e, graph){ 340 426 //C^{\Pi}_{i,j} 341 mod_pot = length[e]potential[ G.head(e)]+potential[G.tail(e)];427 mod_pot = length[e]potential[graph.head(e)]+potential[graph.tail(e)]; 342 428 fl_e = flow[e]; 343 429 // std::cout << fl_e << std::endl;
Note: See TracChangeset
for help on using the changeset viewer.