marci@764
|
1 |
// -*- c++ -*-
|
marci@764
|
2 |
#include <iostream>
|
marci@764
|
3 |
#include <fstream>
|
marci@764
|
4 |
|
alpar@921
|
5 |
#include <lemon/smart_graph.h>
|
marci@1014
|
6 |
#include <lemon/list_graph.h>
|
alpar@921
|
7 |
#include <lemon/dimacs.h>
|
alpar@921
|
8 |
#include <lemon/time_measure.h>
|
marci@764
|
9 |
//#include <graph_wrapper.h>
|
marci@1014
|
10 |
#include <lemon/preflow.h>
|
marci@764
|
11 |
#include <augmenting_flow.h>
|
marci@764
|
12 |
//#include <preflow_res.h>
|
marci@764
|
13 |
#include <lp_solver_wrapper.h>
|
marci@764
|
14 |
|
alpar@921
|
15 |
using namespace lemon;
|
marci@764
|
16 |
|
marci@764
|
17 |
// Use a DIMACS max flow file as stdin.
|
marci@764
|
18 |
// max_flow_demo < dimacs_max_flow_file
|
marci@764
|
19 |
|
marci@764
|
20 |
template<typename Edge, typename EdgeIndexMap>
|
marci@764
|
21 |
class PrimalMap {
|
marci@764
|
22 |
protected:
|
marci@764
|
23 |
LPSolverWrapper* lp;
|
marci@764
|
24 |
EdgeIndexMap* edge_index_map;
|
marci@764
|
25 |
public:
|
marci@764
|
26 |
PrimalMap(LPSolverWrapper& _lp, EdgeIndexMap& _edge_index_map) :
|
marci@764
|
27 |
lp(&_lp), edge_index_map(&_edge_index_map) { }
|
marci@764
|
28 |
double operator[](Edge e) const {
|
marci@764
|
29 |
return lp->getPrimal((*edge_index_map)[e]);
|
marci@764
|
30 |
}
|
marci@764
|
31 |
};
|
marci@764
|
32 |
|
marci@764
|
33 |
int main(int, char **) {
|
marci@764
|
34 |
|
marci@1014
|
35 |
typedef ListGraph MutableGraph;
|
marci@1014
|
36 |
typedef SmartGraph Graph;
|
marci@764
|
37 |
typedef Graph::Node Node;
|
marci@764
|
38 |
typedef Graph::EdgeIt EdgeIt;
|
marci@764
|
39 |
|
marci@764
|
40 |
Graph g;
|
marci@764
|
41 |
|
marci@764
|
42 |
Node s, t;
|
marci@764
|
43 |
Graph::EdgeMap<int> cap(g);
|
marci@764
|
44 |
//readDimacsMaxFlow(std::cin, g, s, t, cap);
|
marci@764
|
45 |
readDimacs(std::cin, g, cap, s, t);
|
marci@764
|
46 |
Timer ts;
|
marci@764
|
47 |
Graph::EdgeMap<int> flow(g); //0 flow
|
marci@1014
|
48 |
Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >
|
marci@764
|
49 |
max_flow_test(g, s, t, cap, flow);
|
marci@764
|
50 |
AugmentingFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >
|
marci@764
|
51 |
augmenting_flow_test(g, s, t, cap, flow);
|
marci@764
|
52 |
|
marci@764
|
53 |
Graph::NodeMap<bool> cut(g);
|
marci@764
|
54 |
|
marci@764
|
55 |
{
|
marci@764
|
56 |
std::cout << "preflow ..." << std::endl;
|
marci@764
|
57 |
ts.reset();
|
marci@764
|
58 |
max_flow_test.run();
|
marci@764
|
59 |
std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
60 |
std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
|
marci@1014
|
61 |
max_flow_test.minCut(cut);
|
marci@764
|
62 |
|
marci@1014
|
63 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) {
|
alpar@986
|
64 |
if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
|
marci@764
|
65 |
std::cout << "Slackness does not hold!" << std::endl;
|
alpar@986
|
66 |
if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
|
marci@764
|
67 |
std::cout << "Slackness does not hold!" << std::endl;
|
marci@764
|
68 |
}
|
marci@764
|
69 |
}
|
marci@764
|
70 |
|
marci@764
|
71 |
// {
|
marci@764
|
72 |
// std::cout << "preflow ..." << std::endl;
|
marci@764
|
73 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
|
marci@764
|
74 |
// ts.reset();
|
marci@1014
|
75 |
// max_flow_test.preflow(Preflow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >::GEN_FLOW);
|
marci@764
|
76 |
// std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
77 |
// std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
|
marci@764
|
78 |
|
marci@764
|
79 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
|
alpar@986
|
80 |
// if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
|
marci@764
|
81 |
// std::cout << "Slackness does not hold!" << std::endl;
|
alpar@986
|
82 |
// if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
|
marci@764
|
83 |
// std::cout << "Slackness does not hold!" << std::endl;
|
marci@764
|
84 |
// }
|
marci@764
|
85 |
// }
|
marci@764
|
86 |
|
marci@764
|
87 |
// {
|
marci@764
|
88 |
// std::cout << "wrapped preflow ..." << std::endl;
|
marci@764
|
89 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
|
marci@764
|
90 |
// ts.reset();
|
marci@764
|
91 |
// pre_flow_res.run();
|
marci@764
|
92 |
// std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
93 |
// std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl;
|
marci@764
|
94 |
// }
|
marci@764
|
95 |
|
marci@764
|
96 |
{
|
marci@764
|
97 |
std::cout << "physical blocking flow augmentation ..." << std::endl;
|
marci@1014
|
98 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0);
|
marci@764
|
99 |
ts.reset();
|
marci@764
|
100 |
int i=0;
|
marci@764
|
101 |
while (augmenting_flow_test.augmentOnBlockingFlow<MutableGraph>()) { ++i; }
|
marci@764
|
102 |
std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
103 |
std::cout << "number of augmentation phases: " << i << std::endl;
|
marci@764
|
104 |
std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
|
marci@764
|
105 |
|
marci@1014
|
106 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) {
|
alpar@986
|
107 |
if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
|
marci@764
|
108 |
std::cout << "Slackness does not hold!" << std::endl;
|
alpar@986
|
109 |
if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
|
marci@764
|
110 |
std::cout << "Slackness does not hold!" << std::endl;
|
marci@764
|
111 |
}
|
marci@764
|
112 |
}
|
marci@764
|
113 |
|
marci@764
|
114 |
// {
|
marci@764
|
115 |
// std::cout << "faster physical blocking flow augmentation ..." << std::endl;
|
marci@764
|
116 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
|
marci@764
|
117 |
// ts.reset();
|
marci@764
|
118 |
// int i=0;
|
marci@764
|
119 |
// while (max_flow_test.augmentOnBlockingFlow1<MutableGraph>()) { ++i; }
|
marci@764
|
120 |
// std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
121 |
// std::cout << "number of augmentation phases: " << i << std::endl;
|
marci@764
|
122 |
// std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
|
marci@764
|
123 |
// }
|
marci@764
|
124 |
|
marci@764
|
125 |
{
|
marci@764
|
126 |
std::cout << "on-the-fly blocking flow augmentation ..." << std::endl;
|
marci@1014
|
127 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0);
|
marci@764
|
128 |
ts.reset();
|
marci@764
|
129 |
int i=0;
|
marci@764
|
130 |
while (augmenting_flow_test.augmentOnBlockingFlow2()) { ++i; }
|
marci@764
|
131 |
std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
132 |
std::cout << "number of augmentation phases: " << i << std::endl;
|
marci@764
|
133 |
std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
|
marci@764
|
134 |
|
marci@1014
|
135 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) {
|
alpar@986
|
136 |
if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
|
marci@764
|
137 |
std::cout << "Slackness does not hold!" << std::endl;
|
alpar@986
|
138 |
if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
|
marci@764
|
139 |
std::cout << "Slackness does not hold!" << std::endl;
|
marci@764
|
140 |
}
|
marci@764
|
141 |
}
|
marci@764
|
142 |
|
marci@764
|
143 |
// {
|
marci@764
|
144 |
// std::cout << "on-the-fly shortest path augmentation ..." << std::endl;
|
marci@764
|
145 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
|
marci@764
|
146 |
// ts.reset();
|
marci@764
|
147 |
// int i=0;
|
marci@764
|
148 |
// while (augmenting_flow_test.augmentOnShortestPath()) { ++i; }
|
marci@764
|
149 |
// std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
150 |
// std::cout << "number of augmentation phases: " << i << std::endl;
|
marci@764
|
151 |
// std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
|
marci@764
|
152 |
|
marci@764
|
153 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
|
alpar@986
|
154 |
// if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
|
marci@764
|
155 |
// std::cout << "Slackness does not hold!" << std::endl;
|
alpar@986
|
156 |
// if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
|
marci@764
|
157 |
// std::cout << "Slackness does not hold!" << std::endl;
|
marci@764
|
158 |
// }
|
marci@764
|
159 |
// }
|
marci@764
|
160 |
|
marci@764
|
161 |
// {
|
marci@764
|
162 |
// std::cout << "on-the-fly shortest path augmentation ..." << std::endl;
|
marci@764
|
163 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
|
marci@764
|
164 |
// ts.reset();
|
marci@764
|
165 |
// int i=0;
|
marci@764
|
166 |
// while (augmenting_flow_test.augmentOnShortestPath2()) { ++i; }
|
marci@764
|
167 |
// std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
168 |
// std::cout << "number of augmentation phases: " << i << std::endl;
|
marci@764
|
169 |
// std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
|
marci@764
|
170 |
|
marci@764
|
171 |
// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
|
alpar@986
|
172 |
// if (cut[g.source(e)] && !cut[g.target(e)] && !flow[e]==cap[e])
|
marci@764
|
173 |
// std::cout << "Slackness does not hold!" << std::endl;
|
alpar@986
|
174 |
// if (!cut[g.source(e)] && cut[g.target(e)] && flow[e]>0)
|
marci@764
|
175 |
// std::cout << "Slackness does not hold!" << std::endl;
|
marci@764
|
176 |
// }
|
marci@764
|
177 |
// }
|
marci@764
|
178 |
|
marci@764
|
179 |
ts.reset();
|
marci@764
|
180 |
LPSolverWrapper lp;
|
marci@764
|
181 |
lp.setMaximize();
|
marci@764
|
182 |
typedef LPSolverWrapper::ColIt ColIt;
|
marci@764
|
183 |
typedef LPSolverWrapper::RowIt RowIt;
|
marci@764
|
184 |
typedef Graph::EdgeMap<ColIt> EdgeIndexMap;
|
marci@764
|
185 |
EdgeIndexMap edge_index_map(g);
|
marci@764
|
186 |
PrimalMap<Graph::Edge, EdgeIndexMap> lp_flow(lp, edge_index_map);
|
marci@1014
|
187 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) {
|
marci@764
|
188 |
ColIt col_it=lp.addCol();
|
marci@764
|
189 |
edge_index_map.set(e, col_it);
|
marci@764
|
190 |
lp.setColBounds(col_it, LPX_DB, 0.0, cap[e]);
|
marci@764
|
191 |
}
|
marci@1014
|
192 |
for (Graph::NodeIt n(g); n!=INVALID; ++n) {
|
marci@764
|
193 |
if (n!=s) {
|
marci@764
|
194 |
//hurokelek miatt
|
marci@764
|
195 |
Graph::EdgeMap<int> coeffs(g, 0);
|
marci@1014
|
196 |
for (Graph::InEdgeIt e(g, n); e!=INVALID; ++e)
|
marci@1014
|
197 |
coeffs.set(e, coeffs[e]+1);
|
marci@1014
|
198 |
for (Graph::OutEdgeIt e(g, n); e!=INVALID; ++e)
|
marci@1014
|
199 |
coeffs.set(e, coeffs[e]-1);
|
marci@764
|
200 |
if (n==t) {
|
marci@1014
|
201 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) {
|
marci@764
|
202 |
if (coeffs[e]!=0)
|
marci@764
|
203 |
lp.setObjCoef(edge_index_map[e], coeffs[e]);
|
marci@764
|
204 |
}
|
marci@764
|
205 |
} else {
|
marci@764
|
206 |
RowIt row_it=lp.addRow();
|
marci@764
|
207 |
std::vector< std::pair<ColIt, double> > row;
|
marci@1014
|
208 |
for (Graph::EdgeIt e(g); e!=INVALID; ++e) {
|
marci@764
|
209 |
if (coeffs[e]!=0)
|
marci@764
|
210 |
row.push_back(std::make_pair(edge_index_map[e], coeffs[e]));
|
marci@764
|
211 |
}
|
marci@764
|
212 |
lp.setRowCoeffs(row_it, row.begin(), row.end());
|
marci@764
|
213 |
lp.setRowBounds(row_it, LPX_FX, 0.0, 0.0);
|
marci@764
|
214 |
}
|
marci@764
|
215 |
}
|
marci@764
|
216 |
}
|
marci@764
|
217 |
lp.solveSimplex();
|
marci@764
|
218 |
std::cout << "flow value: "<< lp.getObjVal() << std::endl;
|
marci@764
|
219 |
std::cout << "elapsed time: " << ts << std::endl;
|
marci@764
|
220 |
|
marci@764
|
221 |
return 0;
|
marci@764
|
222 |
}
|