COIN-OR::LEMON - Graph Library

source: lemon-0.x/test/min_cost_flow_test.cc

Last change on this file was 2621:814ba94d9989, checked in by Peter Kovacs, 11 years ago

Bug fix in min_cost_flow_test.cc

File size: 15.3 KB
Line 
1/* -*- C++ -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library
4 *
5 * Copyright (C) 2003-2008
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 *
9 * Permission to use, modify and distribute this software is granted
10 * provided that this copyright notice appears in all copies. For
11 * precise terms see the accompanying LICENSE file.
12 *
13 * This software is provided "AS IS" with no warranty of any kind,
14 * express or implied, and with no claim as to its suitability for any
15 * purpose.
16 *
17 */
18
19#include <iostream>
20#include <fstream>
21
22#include <lemon/list_graph.h>
23#include <lemon/smart_graph.h>
24#include <lemon/graph_reader.h>
25#include <lemon/dimacs.h>
26#include <lemon/time_measure.h>
27
28#include <lemon/cycle_canceling.h>
29#include <lemon/capacity_scaling.h>
30#include <lemon/cost_scaling.h>
31#include <lemon/network_simplex.h>
32
33#include <lemon/min_cost_flow.h>
34#include <lemon/min_cost_max_flow.h>
35
36#include "test_tools.h"
37
38using namespace lemon;
39
40// Checks the feasibility of a flow
41template < typename Graph, typename LowerMap, typename CapacityMap,
42           typename SupplyMap, typename FlowMap >
43bool checkFlow( const Graph& gr,
44                const LowerMap& lower, const CapacityMap& upper,
45                const SupplyMap& supply, const FlowMap& flow )
46{
47  GRAPH_TYPEDEFS(typename Graph);
48  for (EdgeIt e(gr); e != INVALID; ++e)
49    if (flow[e] < lower[e] || flow[e] > upper[e]) return false;
50
51  for (NodeIt n(gr); n != INVALID; ++n) {
52    typename SupplyMap::Value sum = 0;
53    for (OutEdgeIt e(gr, n); e != INVALID; ++e)
54      sum += flow[e];
55    for (InEdgeIt e(gr, n); e != INVALID; ++e)
56      sum -= flow[e];
57    if (sum != supply[n]) return false;
58  }
59
60  return true;
61}
62
63// Checks the optimalitiy of a flow
64template < typename Graph, typename LowerMap, typename CapacityMap,
65           typename CostMap, typename FlowMap, typename PotentialMap >
66bool checkOptimality( const Graph& gr, const LowerMap& lower,
67                      const CapacityMap& upper, const CostMap& cost,
68                      const FlowMap& flow, const PotentialMap& pi )
69{
70  GRAPH_TYPEDEFS(typename Graph);
71  // Checking the Complementary Slackness optimality condition
72  bool opt = true;
73  for (EdgeIt e(gr); e != INVALID; ++e) {
74    typename CostMap::Value red_cost =
75      cost[e] + pi[gr.source(e)] - pi[gr.target(e)];
76    opt = red_cost == 0 ||
77          (red_cost > 0 && flow[e] == lower[e]) ||
78          (red_cost < 0 && flow[e] == upper[e]);
79    if (!opt) break;
80  }
81  return opt;
82}
83
84// Runs a minimum cost flow algorithm and checks the results
85template < typename MinCostFlowImpl, typename Graph,
86           typename LowerMap, typename CapacityMap,
87           typename CostMap, typename SupplyMap >
88void checkMcf( std::string test_id,
89               const MinCostFlowImpl& mcf, const Graph& gr,
90               const LowerMap& lower, const CapacityMap& upper,
91               const CostMap& cost, const SupplyMap& supply,
92               bool mcf_result, bool result,
93               typename CostMap::Value total )
94{
95  check(mcf_result == result, "Wrong result " + test_id);
96  if (result) {
97    check(checkFlow(gr, lower, upper, supply, mcf.flowMap()),
98          "The flow is not feasible " + test_id);
99    check(mcf.totalCost() == total, "The flow is not optimal " + test_id);
100    check(checkOptimality(gr, lower, upper, cost, mcf.flowMap(), mcf.potentialMap()),
101          "Wrong potentials " + test_id);
102  }
103}
104
105int main()
106{
107  // Various tests on a small graph
108  {
109    typedef ListGraph Graph;
110    GRAPH_TYPEDEFS(ListGraph);
111
112    // Reading the test graph
113    Graph gr;
114    Graph::EdgeMap<int> c(gr), l1(gr), l2(gr), u(gr);
115    Graph::NodeMap<int> s1(gr), s2(gr), s3(gr);
116    Node v, w;
117
118    std::string fname;
119    if(getenv("srcdir"))
120      fname = std::string(getenv("srcdir"));
121    else fname = ".";
122    fname += "/test/min_cost_flow_test.lgf";
123
124    std::ifstream input(fname.c_str());
125    check(input, "Input file '" << fname << "' not found");
126    GraphReader<Graph>(input, gr).
127      readEdgeMap("cost", c).
128      readEdgeMap("capacity", u).
129      readEdgeMap("lower1", l1).
130      readEdgeMap("lower2", l2).
131      readNodeMap("supply1", s1).
132      readNodeMap("supply2", s2).
133      readNodeMap("supply3", s3).
134      readNode("source", v).
135      readNode("target", w).
136      run();
137    input.close();
138
139    // Testing CapacityScaling (scaling enabled)
140    {
141      CapacityScaling<Graph> mcf1(gr,u,c,s1);        checkMcf("#A1",mcf1,gr,l1,u,c,s1,mcf1.run(),true,    0);
142      CapacityScaling<Graph> mcf2(gr,u,c,s2);        checkMcf("#A2",mcf2,gr,l1,u,c,s2,mcf2.run(),true, 5240);
143      CapacityScaling<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#A3",mcf3,gr,l1,u,c,s3,mcf3.run(),true, 7620);
144      CapacityScaling<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#A4",mcf4,gr,l2,u,c,s1,mcf4.run(),false,   0);
145      CapacityScaling<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#A5",mcf5,gr,l2,u,c,s2,mcf5.run(),true, 5970);
146      CapacityScaling<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#A6",mcf6,gr,l2,u,c,s3,mcf6.run(),true, 8010);
147    }
148    // Testing CapacityScaling (scaling disabled)
149    {
150      CapacityScaling<Graph> mcf1(gr,u,c,s1);        checkMcf("#B1",mcf1,gr,l1,u,c,s1,mcf1.run(false),true,    0);
151      CapacityScaling<Graph> mcf2(gr,u,c,s2);        checkMcf("#B2",mcf2,gr,l1,u,c,s2,mcf2.run(false),true, 5240);
152      CapacityScaling<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#B3",mcf3,gr,l1,u,c,s3,mcf3.run(false),true, 7620);
153      CapacityScaling<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#B4",mcf4,gr,l2,u,c,s1,mcf4.run(false),false,   0);
154      CapacityScaling<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#B5",mcf5,gr,l2,u,c,s2,mcf5.run(false),true, 5970);
155      CapacityScaling<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#B6",mcf6,gr,l2,u,c,s3,mcf6.run(false),true, 8010);
156    }
157
158    // Testing CostScaling
159    {
160      CostScaling<Graph> mcf1(gr,u,c,s1);        checkMcf("#C1",mcf1,gr,l1,u,c,s1,mcf1.run(),true,    0);
161      CostScaling<Graph> mcf2(gr,u,c,s2);        checkMcf("#C2",mcf2,gr,l1,u,c,s2,mcf2.run(),true, 5240);
162      CostScaling<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#C3",mcf3,gr,l1,u,c,s3,mcf3.run(),true, 7620);
163      CostScaling<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#C4",mcf4,gr,l2,u,c,s1,mcf4.run(),false,   0);
164      CostScaling<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#C5",mcf5,gr,l2,u,c,s2,mcf5.run(),true, 5970);
165      CostScaling<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#C6",mcf6,gr,l2,u,c,s3,mcf6.run(),true, 8010);
166    }
167
168    // Testing NetworkSimplex (with the default pivot rule)
169    {
170      NetworkSimplex<Graph> mcf1(gr,u,c,s1);        checkMcf("#D1",mcf1,gr,l1,u,c,s1,mcf1.run(),true,    0);
171      NetworkSimplex<Graph> mcf2(gr,u,c,s2);        checkMcf("#D2",mcf2,gr,l1,u,c,s2,mcf2.run(),true, 5240);
172      NetworkSimplex<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#D3",mcf3,gr,l1,u,c,s3,mcf3.run(),true, 7620);
173      NetworkSimplex<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#D4",mcf4,gr,l2,u,c,s1,mcf4.run(),false,   0);
174      NetworkSimplex<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#D5",mcf5,gr,l2,u,c,s2,mcf5.run(),true, 5970);
175      NetworkSimplex<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#D6",mcf6,gr,l2,u,c,s3,mcf6.run(),true, 8010);
176    }
177    // Testing NetworkSimplex (with FIRST_ELIGIBLE_PIVOT)
178    {
179      NetworkSimplex<Graph>::PivotRuleEnum pr = NetworkSimplex<Graph>::FIRST_ELIGIBLE_PIVOT;
180      NetworkSimplex<Graph> mcf1(gr,u,c,s1);        checkMcf("#E1",mcf1,gr,l1,u,c,s1,mcf1.run(pr),true,    0);
181      NetworkSimplex<Graph> mcf2(gr,u,c,s2);        checkMcf("#E2",mcf2,gr,l1,u,c,s2,mcf2.run(pr),true, 5240);
182      NetworkSimplex<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#E3",mcf3,gr,l1,u,c,s3,mcf3.run(pr),true, 7620);
183      NetworkSimplex<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#E4",mcf4,gr,l2,u,c,s1,mcf4.run(pr),false,   0);
184      NetworkSimplex<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#E5",mcf5,gr,l2,u,c,s2,mcf5.run(pr),true, 5970);
185      NetworkSimplex<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#E6",mcf6,gr,l2,u,c,s3,mcf6.run(pr),true, 8010);
186    }
187    // Testing NetworkSimplex (with BEST_ELIGIBLE_PIVOT)
188    {
189      NetworkSimplex<Graph>::PivotRuleEnum pr = NetworkSimplex<Graph>::BEST_ELIGIBLE_PIVOT;
190      NetworkSimplex<Graph> mcf1(gr,u,c,s1);        checkMcf("#F1",mcf1,gr,l1,u,c,s1,mcf1.run(pr),true,    0);
191      NetworkSimplex<Graph> mcf2(gr,u,c,s2);        checkMcf("#F2",mcf2,gr,l1,u,c,s2,mcf2.run(pr),true, 5240);
192      NetworkSimplex<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#F3",mcf3,gr,l1,u,c,s3,mcf3.run(pr),true, 7620);
193      NetworkSimplex<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#F4",mcf4,gr,l2,u,c,s1,mcf4.run(pr),false,   0);
194      NetworkSimplex<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#F5",mcf5,gr,l2,u,c,s2,mcf5.run(pr),true, 5970);
195      NetworkSimplex<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#F6",mcf6,gr,l2,u,c,s3,mcf6.run(pr),true, 8010);
196    }
197    // Testing NetworkSimplex (with BLOCK_SEARCH_PIVOT)
198    {
199      NetworkSimplex<Graph>::PivotRuleEnum pr = NetworkSimplex<Graph>::BLOCK_SEARCH_PIVOT;
200      NetworkSimplex<Graph> mcf1(gr,u,c,s1);        checkMcf("#G1",mcf1,gr,l1,u,c,s1,mcf1.run(pr),true,    0);
201      NetworkSimplex<Graph> mcf2(gr,u,c,s2);        checkMcf("#G2",mcf2,gr,l1,u,c,s2,mcf2.run(pr),true, 5240);
202      NetworkSimplex<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#G3",mcf3,gr,l1,u,c,s3,mcf3.run(pr),true, 7620);
203      NetworkSimplex<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#G4",mcf4,gr,l2,u,c,s1,mcf4.run(pr),false,   0);
204      NetworkSimplex<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#G5",mcf5,gr,l2,u,c,s2,mcf5.run(pr),true, 5970);
205      NetworkSimplex<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#G6",mcf6,gr,l2,u,c,s3,mcf6.run(pr),true, 8010);
206    }
207    // Testing NetworkSimplex (with CANDIDATE_LIST_PIVOT)
208    {
209      NetworkSimplex<Graph>::PivotRuleEnum pr = NetworkSimplex<Graph>::CANDIDATE_LIST_PIVOT;
210      NetworkSimplex<Graph> mcf1(gr,u,c,s1);        checkMcf("#I1",mcf1,gr,l1,u,c,s1,mcf1.run(pr),true,    0);
211      NetworkSimplex<Graph> mcf2(gr,u,c,s2);        checkMcf("#I2",mcf2,gr,l1,u,c,s2,mcf2.run(pr),true, 5240);
212      NetworkSimplex<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#I3",mcf3,gr,l1,u,c,s3,mcf3.run(pr),true, 7620);
213      NetworkSimplex<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#I4",mcf4,gr,l2,u,c,s1,mcf4.run(pr),false,   0);
214      NetworkSimplex<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#I5",mcf5,gr,l2,u,c,s2,mcf5.run(pr),true, 5970);
215      NetworkSimplex<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#I6",mcf6,gr,l2,u,c,s3,mcf6.run(pr),true, 8010);
216    }
217    // Testing NetworkSimplex (with ALTERING_LIST_PIVOT)
218    {
219      NetworkSimplex<Graph>::PivotRuleEnum pr = NetworkSimplex<Graph>::ALTERING_LIST_PIVOT;
220      NetworkSimplex<Graph> mcf1(gr,u,c,s1);        checkMcf("#H1",mcf1,gr,l1,u,c,s1,mcf1.run(pr),true,    0);
221      NetworkSimplex<Graph> mcf2(gr,u,c,s2);        checkMcf("#H2",mcf2,gr,l1,u,c,s2,mcf2.run(pr),true, 5240);
222      NetworkSimplex<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#H3",mcf3,gr,l1,u,c,s3,mcf3.run(pr),true, 7620);
223      NetworkSimplex<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#H4",mcf4,gr,l2,u,c,s1,mcf4.run(pr),false,   0);
224      NetworkSimplex<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#H5",mcf5,gr,l2,u,c,s2,mcf5.run(pr),true, 5970);
225      NetworkSimplex<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#H6",mcf6,gr,l2,u,c,s3,mcf6.run(pr),true, 8010);
226    }
227
228    // Testing CycleCanceling (with BellmanFord)
229    {
230      CycleCanceling<Graph> mcf1(gr,u,c,s1);        checkMcf("#J1",mcf1,gr,l1,u,c,s1,mcf1.run(),true,    0);
231      CycleCanceling<Graph> mcf2(gr,u,c,s2);        checkMcf("#J2",mcf2,gr,l1,u,c,s2,mcf2.run(),true, 5240);
232      CycleCanceling<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#J3",mcf3,gr,l1,u,c,s3,mcf3.run(),true, 7620);
233      CycleCanceling<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#J4",mcf4,gr,l2,u,c,s1,mcf4.run(),false,   0);
234      CycleCanceling<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#J5",mcf5,gr,l2,u,c,s2,mcf5.run(),true, 5970);
235      CycleCanceling<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#J6",mcf6,gr,l2,u,c,s3,mcf6.run(),true, 8010);
236    }
237    // Testing CycleCanceling (with MinMeanCycle)
238    {
239      CycleCanceling<Graph> mcf1(gr,u,c,s1);        checkMcf("#K1",mcf1,gr,l1,u,c,s1,mcf1.run(true),true,    0);
240      CycleCanceling<Graph> mcf2(gr,u,c,s2);        checkMcf("#K2",mcf2,gr,l1,u,c,s2,mcf2.run(true),true, 5240);
241      CycleCanceling<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#K3",mcf3,gr,l1,u,c,s3,mcf3.run(true),true, 7620);
242      CycleCanceling<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#K4",mcf4,gr,l2,u,c,s1,mcf4.run(true),false,   0);
243      CycleCanceling<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#K5",mcf5,gr,l2,u,c,s2,mcf5.run(true),true, 5970);
244      CycleCanceling<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#K6",mcf6,gr,l2,u,c,s3,mcf6.run(true),true, 8010);
245    }
246
247    // Testing MinCostFlow
248    {
249      MinCostFlow<Graph> mcf1(gr,u,c,s1);        checkMcf("#L1",mcf1,gr,l1,u,c,s1,mcf1.run(),true,    0);
250      MinCostFlow<Graph> mcf2(gr,u,c,s2);        checkMcf("#L2",mcf2,gr,l1,u,c,s2,mcf2.run(),true, 5240);
251      MinCostFlow<Graph> mcf3(gr,u,c,v,w,27);    checkMcf("#L3",mcf3,gr,l1,u,c,s3,mcf3.run(),true, 7620);
252      MinCostFlow<Graph> mcf4(gr,l2,u,c,s1);     checkMcf("#L4",mcf4,gr,l2,u,c,s1,mcf4.run(),false,   0);
253      MinCostFlow<Graph> mcf5(gr,l2,u,c,s2);     checkMcf("#L5",mcf5,gr,l2,u,c,s2,mcf5.run(),true, 5970);
254      MinCostFlow<Graph> mcf6(gr,l2,u,c,v,w,27); checkMcf("#L6",mcf6,gr,l2,u,c,s3,mcf6.run(),true, 8010);
255    }
256
257    // Testing MinCostMaxFlow
258    {
259      MinCostMaxFlow<Graph> mcmf(gr,u,c,v,w);
260      mcmf.run();
261      checkMcf("#M1",mcmf,gr,l1,u,c,s3,true,true,7620);
262    }
263  }
264
265  // Benchmark test on a DIMACS network
266  {
267    typedef SmartGraph Graph;
268    GRAPH_TYPEDEFS(SmartGraph);
269
270    // Reading the test graph
271    Graph graph;
272    Graph::EdgeMap<int> lower(graph), capacity(graph), cost(graph);
273    Graph::NodeMap<int> supply(graph);
274
275    std::string fname;
276    if(getenv("srcdir"))
277      fname = std::string(getenv("srcdir"));
278    else fname = ".";
279    fname += "/test/min_cost_flow_test.net";
280
281    std::ifstream input(fname.c_str());
282    check(input, "Input file '" << fname << "' not found");
283    readDimacs(input, graph, lower, capacity, cost, supply);
284    input.close();
285
286    // NetworkSimplex
287    {
288      Timer t;
289      NetworkSimplex<Graph> mcf(graph, lower, capacity, cost, supply);
290      bool res = mcf.run();
291      t.stop();
292      checkMcf("#T3", mcf, graph, lower, capacity, cost, supply, res, true, 196587626);
293      std::cout << "NetworkSimplex";
294      std::cout << std::endl << t << std::endl;
295    }
296    // CapacityScaling
297    {
298      Timer t;
299      CapacityScaling<Graph> mcf(graph, lower, capacity, cost, supply);
300      bool res = mcf.run();
301      t.stop();
302      checkMcf("#T1", mcf, graph, lower, capacity, cost, supply, res, true, 196587626);
303      std::cout << "CapacityScaling";
304      std::cout << std::endl << t << std::endl;
305    }
306    // CostScaling
307    {
308      Timer t;
309      CostScaling<Graph> mcf(graph, lower, capacity, cost, supply);
310      bool res = mcf.run();
311      t.stop();
312      checkMcf("#T2", mcf, graph, lower, capacity, cost, supply, res, true, 196587626);
313      std::cout << "CostScaling";
314      std::cout << std::endl << t << std::endl;
315    }
316    // CycleCanceling
317    {
318      Timer t;
319      CycleCanceling<Graph> mcf(graph, lower, capacity, cost, supply);
320      bool res = mcf.run();
321      t.stop();
322      checkMcf("#T4", mcf, graph, lower, capacity, cost, supply, res, true, 196587626);
323      std::cout << "CycleCanceling";
324      std::cout << std::endl << t << std::endl;
325    }
326  }
327
328  return 0;
329}
Note: See TracBrowser for help on using the repository browser.