1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/src/work/athos/lp/max_flow_by_lp.cc Tue Mar 22 11:45:47 2005 +0000
1.3 @@ -0,0 +1,186 @@
1.4 +// -*- c++ -*-
1.5 +#include <iostream>
1.6 +#include <fstream>
1.7 +
1.8 +#include <lemon/smart_graph.h>
1.9 +#include <lemon/list_graph.h>
1.10 +#include <lemon/dimacs.h>
1.11 +#include <lemon/time_measure.h>
1.12 +//#include <graph_wrapper.h>
1.13 +#include <lemon/preflow.h>
1.14 +#include <augmenting_flow.h>
1.15 +//#include <preflow_res.h>
1.16 +//#include <lp_solver_wrapper_2.h>
1.17 +#include <min_cost_gen_flow.h>
1.18 +
1.19 +// Use a DIMACS max flow file as stdin.
1.20 +// max_flow_demo < dimacs_max_flow_file
1.21 +
1.22 +using namespace lemon;
1.23 +
1.24 +int main(int, char **) {
1.25 +
1.26 + typedef ListGraph MutableGraph;
1.27 + typedef ListGraph Graph;
1.28 + typedef Graph::Node Node;
1.29 + typedef Graph::Edge Edge;
1.30 + typedef Graph::EdgeIt EdgeIt;
1.31 +
1.32 + Graph g;
1.33 +
1.34 + Node s, t;
1.35 + Graph::EdgeMap<int> cap(g);
1.36 + //readDimacsMaxFlow(std::cin, g, s, t, cap);
1.37 + readDimacs(std::cin, g, cap, s, t);
1.38 + Timer ts;
1.39 + Graph::EdgeMap<int> flow(g); //0 flow
1.40 + Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >
1.41 + max_flow_test(g, s, t, cap, flow);
1.42 + AugmentingFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >
1.43 + augmenting_flow_test(g, s, t, cap, flow);
1.44 +
1.45 + Graph::NodeMap<bool> cut(g);
1.46 +
1.47 + {
1.48 + std::cout << "preflow ..." << std::endl;
1.49 + ts.reset();
1.50 + max_flow_test.run();
1.51 + std::cout << "elapsed time: " << ts << std::endl;
1.52 + std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
1.53 + max_flow_test.minCut(cut);
1.54 +
1.55 + for (EdgeIt e(g); e!=INVALID; ++e) {
1.56 + if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
1.57 + std::cout << "Slackness does not hold!" << std::endl;
1.58 + if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
1.59 + std::cout << "Slackness does not hold!" << std::endl;
1.60 + }
1.61 + }
1.62 +
1.63 +// {
1.64 +// std::cout << "preflow ..." << std::endl;
1.65 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
1.66 +// ts.reset();
1.67 +// max_flow_test.preflow(Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >::GEN_FLOW);
1.68 +// std::cout << "elapsed time: " << ts << std::endl;
1.69 +// std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
1.70 +
1.71 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
1.72 +// if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
1.73 +// std::cout << "Slackness does not hold!" << std::endl;
1.74 +// if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
1.75 +// std::cout << "Slackness does not hold!" << std::endl;
1.76 +// }
1.77 +// }
1.78 +
1.79 +// {
1.80 +// std::cout << "wrapped preflow ..." << std::endl;
1.81 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
1.82 +// ts.reset();
1.83 +// pre_flow_res.run();
1.84 +// std::cout << "elapsed time: " << ts << std::endl;
1.85 +// std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl;
1.86 +// }
1.87 +
1.88 + {
1.89 + std::cout << "physical blocking flow augmentation ..." << std::endl;
1.90 + for (EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0);
1.91 + ts.reset();
1.92 + int i=0;
1.93 + while (augmenting_flow_test.augmentOnBlockingFlow<MutableGraph>()) { ++i; }
1.94 + std::cout << "elapsed time: " << ts << std::endl;
1.95 + std::cout << "number of augmentation phases: " << i << std::endl;
1.96 + std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
1.97 +
1.98 + for (EdgeIt e(g); e!=INVALID; ++e) {
1.99 + if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
1.100 + std::cout << "Slackness does not hold!" << std::endl;
1.101 + if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
1.102 + std::cout << "Slackness does not hold!" << std::endl;
1.103 + }
1.104 + }
1.105 +
1.106 +// {
1.107 +// std::cout << "faster physical blocking flow augmentation ..." << std::endl;
1.108 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
1.109 +// ts.reset();
1.110 +// int i=0;
1.111 +// while (max_flow_test.augmentOnBlockingFlow1<MutableGraph>()) { ++i; }
1.112 +// std::cout << "elapsed time: " << ts << std::endl;
1.113 +// std::cout << "number of augmentation phases: " << i << std::endl;
1.114 +// std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
1.115 +// }
1.116 +
1.117 + {
1.118 + std::cout << "on-the-fly blocking flow augmentation ..." << std::endl;
1.119 + for (EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0);
1.120 + ts.reset();
1.121 + int i=0;
1.122 + while (augmenting_flow_test.augmentOnBlockingFlow2()) { ++i; }
1.123 + std::cout << "elapsed time: " << ts << std::endl;
1.124 + std::cout << "number of augmentation phases: " << i << std::endl;
1.125 + std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
1.126 +
1.127 + for (EdgeIt e(g); e!=INVALID; ++e) {
1.128 + if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
1.129 + std::cout << "Slackness does not hold!" << std::endl;
1.130 + if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
1.131 + std::cout << "Slackness does not hold!" << std::endl;
1.132 + }
1.133 + }
1.134 +
1.135 +// {
1.136 +// std::cout << "on-the-fly shortest path augmentation ..." << std::endl;
1.137 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
1.138 +// ts.reset();
1.139 +// int i=0;
1.140 +// while (augmenting_flow_test.augmentOnShortestPath()) { ++i; }
1.141 +// std::cout << "elapsed time: " << ts << std::endl;
1.142 +// std::cout << "number of augmentation phases: " << i << std::endl;
1.143 +// std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
1.144 +
1.145 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
1.146 +// if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
1.147 +// std::cout << "Slackness does not hold!" << std::endl;
1.148 +// if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
1.149 +// std::cout << "Slackness does not hold!" << std::endl;
1.150 +// }
1.151 +// }
1.152 +
1.153 +// {
1.154 +// std::cout << "on-the-fly shortest path augmentation ..." << std::endl;
1.155 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
1.156 +// ts.reset();
1.157 +// int i=0;
1.158 +// while (augmenting_flow_test.augmentOnShortestPath2()) { ++i; }
1.159 +// std::cout << "elapsed time: " << ts << std::endl;
1.160 +// std::cout << "number of augmentation phases: " << i << std::endl;
1.161 +// std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
1.162 +
1.163 +// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
1.164 +// if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
1.165 +// std::cout << "Slackness does not hold!" << std::endl;
1.166 +// if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
1.167 +// std::cout << "Slackness does not hold!" << std::endl;
1.168 +// }
1.169 +// }
1.170 +
1.171 + ts.reset();
1.172 +
1.173 + Edge e=g.addEdge(t, s);
1.174 + Graph::EdgeMap<int> cost(g, 0);
1.175 + cost.set(e, -1);
1.176 + cap.set(e, 10000);
1.177 + typedef ConstMap<Node, int> Excess;
1.178 + Excess excess(0);
1.179 + typedef ConstMap<Edge, int> LCap;
1.180 + LCap lcap(0);
1.181 +
1.182 + MinCostGenFlow<Graph, int, Excess, LCap>
1.183 + min_cost(g, excess, lcap, cap, flow, cost);
1.184 + min_cost.feasible();
1.185 + min_cost.runByLP();
1.186 +
1.187 + std::cout << "elapsed time: " << ts << std::endl;
1.188 + std::cout << "flow value: "<< flow[e] << std::endl;
1.189 +}