1 // -*- c++ -*- |
1 // -*- c++ -*- |
2 #include <iostream> |
2 #include <iostream> |
3 #include <fstream> |
3 #include <fstream> |
4 |
4 |
5 #include <sage_graph.h> |
|
6 #include <lemon/smart_graph.h> |
5 #include <lemon/smart_graph.h> |
|
6 #include <lemon/list_graph.h> |
7 #include <lemon/dimacs.h> |
7 #include <lemon/dimacs.h> |
8 #include <lemon/time_measure.h> |
8 #include <lemon/time_measure.h> |
9 //#include <graph_wrapper.h> |
9 //#include <graph_wrapper.h> |
10 #include <lemon/max_flow.h> |
10 #include <lemon/preflow.h> |
11 #include <augmenting_flow.h> |
11 #include <augmenting_flow.h> |
12 //#include <preflow_res.h> |
12 //#include <preflow_res.h> |
13 #include <for_each_macros.h> |
|
14 #include <lp_solver_wrapper.h> |
13 #include <lp_solver_wrapper.h> |
15 |
14 |
16 using namespace lemon; |
15 using namespace lemon; |
17 |
16 |
18 // Use a DIMACS max flow file as stdin. |
17 // Use a DIMACS max flow file as stdin. |
45 Graph::EdgeMap<int> cap(g); |
43 Graph::EdgeMap<int> cap(g); |
46 //readDimacsMaxFlow(std::cin, g, s, t, cap); |
44 //readDimacsMaxFlow(std::cin, g, s, t, cap); |
47 readDimacs(std::cin, g, cap, s, t); |
45 readDimacs(std::cin, g, cap, s, t); |
48 Timer ts; |
46 Timer ts; |
49 Graph::EdgeMap<int> flow(g); //0 flow |
47 Graph::EdgeMap<int> flow(g); //0 flow |
50 MaxFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > |
48 Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > |
51 max_flow_test(g, s, t, cap, flow); |
49 max_flow_test(g, s, t, cap, flow); |
52 AugmentingFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > |
50 AugmentingFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> > |
53 augmenting_flow_test(g, s, t, cap, flow); |
51 augmenting_flow_test(g, s, t, cap, flow); |
54 |
52 |
55 Graph::NodeMap<bool> cut(g); |
53 Graph::NodeMap<bool> cut(g); |
58 std::cout << "preflow ..." << std::endl; |
56 std::cout << "preflow ..." << std::endl; |
59 ts.reset(); |
57 ts.reset(); |
60 max_flow_test.run(); |
58 max_flow_test.run(); |
61 std::cout << "elapsed time: " << ts << std::endl; |
59 std::cout << "elapsed time: " << ts << std::endl; |
62 std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl; |
60 std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl; |
63 max_flow_test.actMinCut(cut); |
61 max_flow_test.minCut(cut); |
64 |
62 |
65 FOR_EACH_LOC(Graph::EdgeIt, e, g) { |
63 for (Graph::EdgeIt e(g); e!=INVALID; ++e) { |
66 if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
64 if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
67 std::cout << "Slackness does not hold!" << std::endl; |
65 std::cout << "Slackness does not hold!" << std::endl; |
68 if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0) |
66 if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0) |
69 std::cout << "Slackness does not hold!" << std::endl; |
67 std::cout << "Slackness does not hold!" << std::endl; |
70 } |
68 } |
72 |
70 |
73 // { |
71 // { |
74 // std::cout << "preflow ..." << std::endl; |
72 // std::cout << "preflow ..." << std::endl; |
75 // FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0); |
73 // FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0); |
76 // ts.reset(); |
74 // ts.reset(); |
77 // max_flow_test.preflow(MaxFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >::GEN_FLOW); |
75 // max_flow_test.preflow(Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >::GEN_FLOW); |
78 // std::cout << "elapsed time: " << ts << std::endl; |
76 // std::cout << "elapsed time: " << ts << std::endl; |
79 // std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl; |
77 // std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl; |
80 |
78 |
81 // FOR_EACH_LOC(Graph::EdgeIt, e, g) { |
79 // FOR_EACH_LOC(Graph::EdgeIt, e, g) { |
82 // if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
80 // if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
95 // std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl; |
93 // std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl; |
96 // } |
94 // } |
97 |
95 |
98 { |
96 { |
99 std::cout << "physical blocking flow augmentation ..." << std::endl; |
97 std::cout << "physical blocking flow augmentation ..." << std::endl; |
100 FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0); |
98 for (Graph::EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0); |
101 ts.reset(); |
99 ts.reset(); |
102 int i=0; |
100 int i=0; |
103 while (augmenting_flow_test.augmentOnBlockingFlow<MutableGraph>()) { ++i; } |
101 while (augmenting_flow_test.augmentOnBlockingFlow<MutableGraph>()) { ++i; } |
104 std::cout << "elapsed time: " << ts << std::endl; |
102 std::cout << "elapsed time: " << ts << std::endl; |
105 std::cout << "number of augmentation phases: " << i << std::endl; |
103 std::cout << "number of augmentation phases: " << i << std::endl; |
106 std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl; |
104 std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl; |
107 |
105 |
108 FOR_EACH_LOC(Graph::EdgeIt, e, g) { |
106 for (Graph::EdgeIt e(g); e!=INVALID; ++e) { |
109 if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
107 if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
110 std::cout << "Slackness does not hold!" << std::endl; |
108 std::cout << "Slackness does not hold!" << std::endl; |
111 if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0) |
109 if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0) |
112 std::cout << "Slackness does not hold!" << std::endl; |
110 std::cout << "Slackness does not hold!" << std::endl; |
113 } |
111 } |
124 // std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl; |
122 // std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl; |
125 // } |
123 // } |
126 |
124 |
127 { |
125 { |
128 std::cout << "on-the-fly blocking flow augmentation ..." << std::endl; |
126 std::cout << "on-the-fly blocking flow augmentation ..." << std::endl; |
129 FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0); |
127 for (Graph::EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0); |
130 ts.reset(); |
128 ts.reset(); |
131 int i=0; |
129 int i=0; |
132 while (augmenting_flow_test.augmentOnBlockingFlow2()) { ++i; } |
130 while (augmenting_flow_test.augmentOnBlockingFlow2()) { ++i; } |
133 std::cout << "elapsed time: " << ts << std::endl; |
131 std::cout << "elapsed time: " << ts << std::endl; |
134 std::cout << "number of augmentation phases: " << i << std::endl; |
132 std::cout << "number of augmentation phases: " << i << std::endl; |
135 std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl; |
133 std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl; |
136 |
134 |
137 FOR_EACH_LOC(Graph::EdgeIt, e, g) { |
135 for (Graph::EdgeIt e(g); e!=INVALID; ++e) { |
138 if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
136 if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e]) |
139 std::cout << "Slackness does not hold!" << std::endl; |
137 std::cout << "Slackness does not hold!" << std::endl; |
140 if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0) |
138 if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0) |
141 std::cout << "Slackness does not hold!" << std::endl; |
139 std::cout << "Slackness does not hold!" << std::endl; |
142 } |
140 } |
184 typedef LPSolverWrapper::ColIt ColIt; |
182 typedef LPSolverWrapper::ColIt ColIt; |
185 typedef LPSolverWrapper::RowIt RowIt; |
183 typedef LPSolverWrapper::RowIt RowIt; |
186 typedef Graph::EdgeMap<ColIt> EdgeIndexMap; |
184 typedef Graph::EdgeMap<ColIt> EdgeIndexMap; |
187 EdgeIndexMap edge_index_map(g); |
185 EdgeIndexMap edge_index_map(g); |
188 PrimalMap<Graph::Edge, EdgeIndexMap> lp_flow(lp, edge_index_map); |
186 PrimalMap<Graph::Edge, EdgeIndexMap> lp_flow(lp, edge_index_map); |
189 Graph::EdgeIt e; |
187 for (Graph::EdgeIt e(g); e!=INVALID; ++e) { |
190 for (g.first(e); g.valid(e); g.next(e)) { |
|
191 ColIt col_it=lp.addCol(); |
188 ColIt col_it=lp.addCol(); |
192 edge_index_map.set(e, col_it); |
189 edge_index_map.set(e, col_it); |
193 lp.setColBounds(col_it, LPX_DB, 0.0, cap[e]); |
190 lp.setColBounds(col_it, LPX_DB, 0.0, cap[e]); |
194 } |
191 } |
195 Graph::NodeIt n; |
192 for (Graph::NodeIt n(g); n!=INVALID; ++n) { |
196 for (g.first(n); g.valid(n); g.next(n)) { |
|
197 if (n!=s) { |
193 if (n!=s) { |
198 //hurokelek miatt |
194 //hurokelek miatt |
199 Graph::EdgeMap<int> coeffs(g, 0); |
195 Graph::EdgeMap<int> coeffs(g, 0); |
200 { |
196 for (Graph::InEdgeIt e(g, n); e!=INVALID; ++e) |
201 Graph::InEdgeIt e; |
197 coeffs.set(e, coeffs[e]+1); |
202 for (g.first(e, n); g.valid(e); g.next(e)) coeffs.set(e, coeffs[e]+1); |
198 for (Graph::OutEdgeIt e(g, n); e!=INVALID; ++e) |
203 } |
199 coeffs.set(e, coeffs[e]-1); |
204 { |
|
205 Graph::OutEdgeIt e; |
|
206 for (g.first(e, n); g.valid(e); g.next(e)) coeffs.set(e, coeffs[e]-1); |
|
207 } |
|
208 if (n==t) { |
200 if (n==t) { |
209 Graph::EdgeIt e; |
201 for (Graph::EdgeIt e(g); e!=INVALID; ++e) { |
210 //std::vector< std::pair<ColIt, double> > row; |
|
211 for (g.first(e); g.valid(e); g.next(e)) { |
|
212 if (coeffs[e]!=0) |
202 if (coeffs[e]!=0) |
213 lp.setObjCoef(edge_index_map[e], coeffs[e]); |
203 lp.setObjCoef(edge_index_map[e], coeffs[e]); |
214 } |
204 } |
215 } else { |
205 } else { |
216 RowIt row_it=lp.addRow(); |
206 RowIt row_it=lp.addRow(); |
217 Graph::EdgeIt e; |
|
218 std::vector< std::pair<ColIt, double> > row; |
207 std::vector< std::pair<ColIt, double> > row; |
219 for (g.first(e); g.valid(e); g.next(e)) { |
208 for (Graph::EdgeIt e(g); e!=INVALID; ++e) { |
220 if (coeffs[e]!=0) |
209 if (coeffs[e]!=0) |
221 row.push_back(std::make_pair(edge_index_map[e], coeffs[e])); |
210 row.push_back(std::make_pair(edge_index_map[e], coeffs[e])); |
222 } |
211 } |
223 lp.setRowCoeffs(row_it, row.begin(), row.end()); |
212 lp.setRowCoeffs(row_it, row.begin(), row.end()); |
224 lp.setRowBounds(row_it, LPX_FX, 0.0, 0.0); |
213 lp.setRowBounds(row_it, LPX_FX, 0.0, 0.0); |