|
1 #include <iostream> |
|
2 #include <set> |
|
3 #include <vector> |
|
4 #include <iterator> |
|
5 |
|
6 #include <cmath> |
|
7 #include <cstdlib> |
|
8 |
|
9 #include <lemon/smart_graph.h> |
|
10 #include <lemon/min_cost_arborescence.h> |
|
11 |
|
12 #include <lemon/graph_utils.h> |
|
13 #include <lemon/time_measure.h> |
|
14 |
|
15 #include <lemon/tolerance.h> |
|
16 |
|
17 using namespace lemon; |
|
18 using namespace std; |
|
19 |
|
20 int main(int argc, const char *argv[]) { |
|
21 srand(time(0)); |
|
22 typedef SmartGraph Graph; |
|
23 GRAPH_TYPEDEFS(Graph); |
|
24 |
|
25 typedef Graph::EdgeMap<double> CostMap; |
|
26 |
|
27 const int n = argc > 1 ? atoi(argv[1]) : 100; |
|
28 const int e = argc > 2 ? atoi(argv[2]) : (int)(n * log(n)); |
|
29 |
|
30 Graph graph; |
|
31 CostMap cost(graph); |
|
32 vector<Node> nodes; |
|
33 |
|
34 for (int i = 0; i < n; ++i) { |
|
35 nodes.push_back(graph.addNode()); |
|
36 } |
|
37 |
|
38 for (int i = 0; i < e; ++i) { |
|
39 int s = (int)(n * (double)rand() / (RAND_MAX + 1.0)); |
|
40 int t = (int)(n * (double)rand() / (RAND_MAX + 1.0)); |
|
41 double c = rand() / (1.0 + RAND_MAX) * 100.0 + 20.0; |
|
42 Edge edge = graph.addEdge(nodes[s], nodes[t]); |
|
43 cost[edge] = c; |
|
44 } |
|
45 |
|
46 Node source = nodes[(int)(n * (double)rand() / (RAND_MAX + 1.0))]; |
|
47 |
|
48 MinCostArborescence<Graph, CostMap> mca(graph, cost); |
|
49 mca.run(source); |
|
50 |
|
51 vector<pair<double, set<Node> > > dualSolution(mca.dualSize()); |
|
52 |
|
53 for (int i = 0; i < mca.dualSize(); ++i) { |
|
54 dualSolution[i].first = mca.dualValue(i); |
|
55 for (MinCostArborescence<Graph, CostMap>::DualIt it(mca, i); |
|
56 it != INVALID; ++it) { |
|
57 dualSolution[i].second.insert(it); |
|
58 } |
|
59 } |
|
60 |
|
61 Tolerance<double> tolerance; |
|
62 |
|
63 for (EdgeIt it(graph); it != INVALID; ++it) { |
|
64 if (mca.reached(graph.source(it))) { |
|
65 double sum = 0.0; |
|
66 for (int i = 0; i < (int)dualSolution.size(); ++i) { |
|
67 if (dualSolution[i].second.find(graph.target(it)) |
|
68 != dualSolution[i].second.end() && |
|
69 dualSolution[i].second.find(graph.source(it)) |
|
70 == dualSolution[i].second.end()) { |
|
71 sum += dualSolution[i].first; |
|
72 } |
|
73 } |
|
74 if (mca.arborescence(it)) { |
|
75 LEMON_ASSERT(!tolerance.less(sum, cost[it]), "INVALID DUAL"); |
|
76 } |
|
77 LEMON_ASSERT(!tolerance.less(cost[it], sum), "INVALID DUAL"); |
|
78 } |
|
79 } |
|
80 |
|
81 |
|
82 LEMON_ASSERT(!tolerance.different(mca.dualValue(), mca.arborescenceValue()), |
|
83 "INVALID DUAL"); |
|
84 |
|
85 |
|
86 LEMON_ASSERT(mca.reached(source), "INVALID ARBORESCENCE"); |
|
87 for (EdgeIt it(graph); it != INVALID; ++it) { |
|
88 LEMON_ASSERT(!mca.reached(graph.source(it)) || |
|
89 mca.reached(graph.target(it)), "INVALID ARBORESCENCE"); |
|
90 } |
|
91 |
|
92 for (NodeIt it(graph); it != INVALID; ++it) { |
|
93 if (!mca.reached(it)) continue; |
|
94 int cnt = 0; |
|
95 for (InEdgeIt jt(graph, it); jt != INVALID; ++jt) { |
|
96 if (mca.arborescence(jt)) { |
|
97 ++cnt; |
|
98 } |
|
99 } |
|
100 LEMON_ASSERT((it == source ? cnt == 0 : cnt == 1), "INVALID ARBORESCENCE"); |
|
101 } |
|
102 |
|
103 return 0; |
|
104 } |