COIN-OR::LEMON - Graph Library

source: lemon-1.2/test/min_cost_flow_test.cc @ 818:bc75ee2ad082

Last change on this file since 818:bc75ee2ad082 was 818:bc75ee2ad082, checked in by Peter Kovacs <kpeter@…>, 10 years ago

Rework the MCF test file to help extending it (#180)

File size: 14.7 KB
Line 
1/* -*- mode: C++; indent-tabs-mode: nil; -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library.
4 *
5 * Copyright (C) 2003-2009
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#include <limits>
22
23#include <lemon/list_graph.h>
24#include <lemon/lgf_reader.h>
25
26#include <lemon/network_simplex.h>
27
28#include <lemon/concepts/digraph.h>
29#include <lemon/concept_check.h>
30
31#include "test_tools.h"
32
33using namespace lemon;
34
35// Test networks
36char test_lgf[] =
37  "@nodes\n"
38  "label  sup1 sup2 sup3 sup4 sup5 sup6\n"
39  "    1    20   27    0   30   20   30\n"
40  "    2    -4    0    0    0   -8   -3\n"
41  "    3     0    0    0    0    0    0\n"
42  "    4     0    0    0    0    0    0\n"
43  "    5     9    0    0    0    6   11\n"
44  "    6    -6    0    0    0   -5   -6\n"
45  "    7     0    0    0    0    0    0\n"
46  "    8     0    0    0    0    0    3\n"
47  "    9     3    0    0    0    0    0\n"
48  "   10    -2    0    0    0   -7   -2\n"
49  "   11     0    0    0    0  -10    0\n"
50  "   12   -20  -27    0  -30  -30  -20\n"
51  "\n"               
52  "@arcs\n"
53  "       cost  cap low1 low2 low3\n"
54  " 1  2    70   11    0    8    8\n"
55  " 1  3   150    3    0    1    0\n"
56  " 1  4    80   15    0    2    2\n"
57  " 2  8    80   12    0    0    0\n"
58  " 3  5   140    5    0    3    1\n"
59  " 4  6    60   10    0    1    0\n"
60  " 4  7    80    2    0    0    0\n"
61  " 4  8   110    3    0    0    0\n"
62  " 5  7    60   14    0    0    0\n"
63  " 5 11   120   12    0    0    0\n"
64  " 6  3     0    3    0    0    0\n"
65  " 6  9   140    4    0    0    0\n"
66  " 6 10    90    8    0    0    0\n"
67  " 7  1    30    5    0    0   -5\n"
68  " 8 12    60   16    0    4    3\n"
69  " 9 12    50    6    0    0    0\n"
70  "10 12    70   13    0    5    2\n"
71  "10  2   100    7    0    0    0\n"
72  "10  7    60   10    0    0   -3\n"
73  "11 10    20   14    0    6  -20\n"
74  "12 11    30   10    0    0  -10\n"
75  "\n"
76  "@attributes\n"
77  "source 1\n"
78  "target 12\n";
79
80char test_neg1_lgf[] =
81  "@nodes\n"
82  "label   sup\n"
83  "    1   100\n"
84  "    2     0\n"
85  "    3     0\n"
86  "    4  -100\n"
87  "    5     0\n"
88  "    6     0\n"
89  "    7     0\n"
90  "@arcs\n"
91  "      cost   low1   low2\n"
92  "1 2    100      0      0\n"
93  "1 3     30      0      0\n"
94  "2 4     20      0      0\n"
95  "3 4     80      0      0\n"
96  "3 2     50      0      0\n"
97  "5 3     10      0      0\n"
98  "5 6     80      0   1000\n"
99  "6 7     30      0  -1000\n"
100  "7 5   -120      0      0\n";
101 
102char test_neg2_lgf[] =
103  "@nodes\n"
104  "label   sup\n"
105  "    1   100\n"
106  "    2  -300\n"
107  "@arcs\n"
108  "      cost\n"
109  "1 2     -1\n";
110
111
112// Test data
113typedef ListDigraph Digraph;
114DIGRAPH_TYPEDEFS(ListDigraph);
115
116Digraph gr;
117Digraph::ArcMap<int> c(gr), l1(gr), l2(gr), l3(gr), u(gr);
118Digraph::NodeMap<int> s1(gr), s2(gr), s3(gr), s4(gr), s5(gr), s6(gr);
119ConstMap<Arc, int> cc(1), cu(std::numeric_limits<int>::max());
120Node v, w;
121
122Digraph neg1_gr;
123Digraph::ArcMap<int> neg1_c(neg1_gr), neg1_l1(neg1_gr), neg1_l2(neg1_gr);
124ConstMap<Arc, int> neg1_u1(std::numeric_limits<int>::max()), neg1_u2(5000);
125Digraph::NodeMap<int> neg1_s(neg1_gr);
126
127Digraph neg2_gr;
128Digraph::ArcMap<int> neg2_c(neg2_gr);
129ConstMap<Arc, int> neg2_l(0), neg2_u(1000);
130Digraph::NodeMap<int> neg2_s(neg2_gr);
131
132
133enum SupplyType {
134  EQ,
135  GEQ,
136  LEQ
137};
138
139
140// Check the interface of an MCF algorithm
141template <typename GR, typename Value, typename Cost>
142class McfClassConcept
143{
144public:
145
146  template <typename MCF>
147  struct Constraints {
148    void constraints() {
149      checkConcept<concepts::Digraph, GR>();
150     
151      const Constraints& me = *this;
152
153      MCF mcf(me.g);
154      const MCF& const_mcf = mcf;
155
156      b = mcf.reset()
157             .lowerMap(me.lower)
158             .upperMap(me.upper)
159             .costMap(me.cost)
160             .supplyMap(me.sup)
161             .stSupply(me.n, me.n, me.k)
162             .run();
163
164      c = const_mcf.totalCost();
165      x = const_mcf.template totalCost<double>();
166      v = const_mcf.flow(me.a);
167      c = const_mcf.potential(me.n);
168      const_mcf.flowMap(fm);
169      const_mcf.potentialMap(pm);
170    }
171
172    typedef typename GR::Node Node;
173    typedef typename GR::Arc Arc;
174    typedef concepts::ReadMap<Node, Value> NM;
175    typedef concepts::ReadMap<Arc, Value> VAM;
176    typedef concepts::ReadMap<Arc, Cost> CAM;
177    typedef concepts::WriteMap<Arc, Value> FlowMap;
178    typedef concepts::WriteMap<Node, Cost> PotMap;
179 
180    GR g;
181    VAM lower;
182    VAM upper;
183    CAM cost;
184    NM sup;
185    Node n;
186    Arc a;
187    Value k;
188
189    FlowMap fm;
190    PotMap pm;
191    bool b;
192    double x;
193    typename MCF::Value v;
194    typename MCF::Cost c;
195  };
196
197};
198
199
200// Check the feasibility of the given flow (primal soluiton)
201template < typename GR, typename LM, typename UM,
202           typename SM, typename FM >
203bool checkFlow( const GR& gr, const LM& lower, const UM& upper,
204                const SM& supply, const FM& flow,
205                SupplyType type = EQ )
206{
207  TEMPLATE_DIGRAPH_TYPEDEFS(GR);
208
209  for (ArcIt e(gr); e != INVALID; ++e) {
210    if (flow[e] < lower[e] || flow[e] > upper[e]) return false;
211  }
212
213  for (NodeIt n(gr); n != INVALID; ++n) {
214    typename SM::Value sum = 0;
215    for (OutArcIt e(gr, n); e != INVALID; ++e)
216      sum += flow[e];
217    for (InArcIt e(gr, n); e != INVALID; ++e)
218      sum -= flow[e];
219    bool b = (type ==  EQ && sum == supply[n]) ||
220             (type == GEQ && sum >= supply[n]) ||
221             (type == LEQ && sum <= supply[n]);
222    if (!b) return false;
223  }
224
225  return true;
226}
227
228// Check the feasibility of the given potentials (dual soluiton)
229// using the "Complementary Slackness" optimality condition
230template < typename GR, typename LM, typename UM,
231           typename CM, typename SM, typename FM, typename PM >
232bool checkPotential( const GR& gr, const LM& lower, const UM& upper,
233                     const CM& cost, const SM& supply, const FM& flow,
234                     const PM& pi, SupplyType type )
235{
236  TEMPLATE_DIGRAPH_TYPEDEFS(GR);
237
238  bool opt = true;
239  for (ArcIt e(gr); opt && e != INVALID; ++e) {
240    typename CM::Value red_cost =
241      cost[e] + pi[gr.source(e)] - pi[gr.target(e)];
242    opt = red_cost == 0 ||
243          (red_cost > 0 && flow[e] == lower[e]) ||
244          (red_cost < 0 && flow[e] == upper[e]);
245  }
246 
247  for (NodeIt n(gr); opt && n != INVALID; ++n) {
248    typename SM::Value sum = 0;
249    for (OutArcIt e(gr, n); e != INVALID; ++e)
250      sum += flow[e];
251    for (InArcIt e(gr, n); e != INVALID; ++e)
252      sum -= flow[e];
253    if (type != LEQ) {
254      opt = (pi[n] <= 0) && (sum == supply[n] || pi[n] == 0);
255    } else {
256      opt = (pi[n] >= 0) && (sum == supply[n] || pi[n] == 0);
257    }
258  }
259 
260  return opt;
261}
262
263// Check whether the dual cost is equal to the primal cost
264template < typename GR, typename LM, typename UM,
265           typename CM, typename SM, typename PM >
266bool checkDualCost( const GR& gr, const LM& lower, const UM& upper,
267                    const CM& cost, const SM& supply, const PM& pi,
268                    typename CM::Value total )
269{
270  TEMPLATE_DIGRAPH_TYPEDEFS(GR);
271
272  typename CM::Value dual_cost = 0;
273  SM red_supply(gr);
274  for (NodeIt n(gr); n != INVALID; ++n) {
275    red_supply[n] = supply[n];
276  }
277  for (ArcIt a(gr); a != INVALID; ++a) {
278    if (lower[a] != 0) {
279      dual_cost += lower[a] * cost[a];
280      red_supply[gr.source(a)] -= lower[a];
281      red_supply[gr.target(a)] += lower[a];
282    }
283  }
284 
285  for (NodeIt n(gr); n != INVALID; ++n) {
286    dual_cost -= red_supply[n] * pi[n];
287  }
288  for (ArcIt a(gr); a != INVALID; ++a) {
289    typename CM::Value red_cost =
290      cost[a] + pi[gr.source(a)] - pi[gr.target(a)];
291    dual_cost -= (upper[a] - lower[a]) * std::max(-red_cost, 0);
292  }
293 
294  return dual_cost == total;
295}
296
297// Run a minimum cost flow algorithm and check the results
298template < typename MCF, typename GR,
299           typename LM, typename UM,
300           typename CM, typename SM,
301           typename PT >
302void checkMcf( const MCF& mcf, PT mcf_result,
303               const GR& gr, const LM& lower, const UM& upper,
304               const CM& cost, const SM& supply,
305               PT result, bool optimal, typename CM::Value total,
306               const std::string &test_id = "",
307               SupplyType type = EQ )
308{
309  check(mcf_result == result, "Wrong result " + test_id);
310  if (optimal) {
311    typename GR::template ArcMap<typename SM::Value> flow(gr);
312    typename GR::template NodeMap<typename CM::Value> pi(gr);
313    mcf.flowMap(flow);
314    mcf.potentialMap(pi);
315    check(checkFlow(gr, lower, upper, supply, flow, type),
316          "The flow is not feasible " + test_id);
317    check(mcf.totalCost() == total, "The flow is not optimal " + test_id);
318    check(checkPotential(gr, lower, upper, cost, supply, flow, pi, type),
319          "Wrong potentials " + test_id);
320    check(checkDualCost(gr, lower, upper, cost, supply, pi, total),
321          "Wrong dual cost " + test_id);
322  }
323}
324
325template < typename MCF, typename Param >
326void runMcfGeqTests( Param param,
327                     const std::string &test_str = "",
328                     bool full_neg_cost_support = false )
329{
330  MCF mcf1(gr), mcf2(neg1_gr), mcf3(neg2_gr);
331 
332  // Basic tests
333  mcf1.upperMap(u).costMap(c).supplyMap(s1);
334  checkMcf(mcf1, mcf1.run(param), gr, l1, u, c, s1,
335           mcf1.OPTIMAL, true,     5240, test_str + "-1");
336  mcf1.stSupply(v, w, 27);
337  checkMcf(mcf1, mcf1.run(param), gr, l1, u, c, s2,
338           mcf1.OPTIMAL, true,     7620, test_str + "-2");
339  mcf1.lowerMap(l2).supplyMap(s1);
340  checkMcf(mcf1, mcf1.run(param), gr, l2, u, c, s1,
341           mcf1.OPTIMAL, true,     5970, test_str + "-3");
342  mcf1.stSupply(v, w, 27);
343  checkMcf(mcf1, mcf1.run(param), gr, l2, u, c, s2,
344           mcf1.OPTIMAL, true,     8010, test_str + "-4");
345  mcf1.reset().supplyMap(s1);
346  checkMcf(mcf1, mcf1.run(param), gr, l1, cu, cc, s1,
347           mcf1.OPTIMAL, true,       74, test_str + "-5");
348  mcf1.lowerMap(l2).stSupply(v, w, 27);
349  checkMcf(mcf1, mcf1.run(param), gr, l2, cu, cc, s2,
350           mcf1.OPTIMAL, true,       94, test_str + "-6");
351  mcf1.reset();
352  checkMcf(mcf1, mcf1.run(param), gr, l1, cu, cc, s3,
353           mcf1.OPTIMAL, true,        0, test_str + "-7");
354  mcf1.lowerMap(l2).upperMap(u);
355  checkMcf(mcf1, mcf1.run(param), gr, l2, u, cc, s3,
356           mcf1.INFEASIBLE, false,    0, test_str + "-8");
357  mcf1.lowerMap(l3).upperMap(u).costMap(c).supplyMap(s4);
358  checkMcf(mcf1, mcf1.run(param), gr, l3, u, c, s4,
359           mcf1.OPTIMAL, true,     6360, test_str + "-9");
360
361  // Tests for the GEQ form
362  mcf1.reset().upperMap(u).costMap(c).supplyMap(s5);
363  checkMcf(mcf1, mcf1.run(param), gr, l1, u, c, s5,
364           mcf1.OPTIMAL, true,     3530, test_str + "-10", GEQ);
365  mcf1.lowerMap(l2);
366  checkMcf(mcf1, mcf1.run(param), gr, l2, u, c, s5,
367           mcf1.OPTIMAL, true,     4540, test_str + "-11", GEQ);
368  mcf1.supplyMap(s6);
369  checkMcf(mcf1, mcf1.run(param), gr, l2, u, c, s6,
370           mcf1.INFEASIBLE, false,    0, test_str + "-12", GEQ);
371
372  // Tests with negative costs
373  mcf2.lowerMap(neg1_l1).costMap(neg1_c).supplyMap(neg1_s);
374  checkMcf(mcf2, mcf2.run(param), neg1_gr, neg1_l1, neg1_u1, neg1_c, neg1_s,
375           mcf2.UNBOUNDED, false,     0, test_str + "-13");
376  mcf2.upperMap(neg1_u2);
377  checkMcf(mcf2, mcf2.run(param), neg1_gr, neg1_l1, neg1_u2, neg1_c, neg1_s,
378           mcf2.OPTIMAL, true,   -40000, test_str + "-14");
379  mcf2.reset().lowerMap(neg1_l2).costMap(neg1_c).supplyMap(neg1_s);
380  checkMcf(mcf2, mcf2.run(param), neg1_gr, neg1_l2, neg1_u1, neg1_c, neg1_s,
381           mcf2.UNBOUNDED, false,     0, test_str + "-15");
382
383  mcf3.costMap(neg2_c).supplyMap(neg2_s);
384  if (full_neg_cost_support) {
385    checkMcf(mcf3, mcf3.run(param), neg2_gr, neg2_l, neg2_u, neg2_c, neg2_s,
386             mcf3.OPTIMAL, true,   -300, test_str + "-16", GEQ);
387  } else {
388    checkMcf(mcf3, mcf3.run(param), neg2_gr, neg2_l, neg2_u, neg2_c, neg2_s,
389             mcf3.UNBOUNDED, false,   0, test_str + "-17", GEQ);
390  }
391  mcf3.upperMap(neg2_u);
392  checkMcf(mcf3, mcf3.run(param), neg2_gr, neg2_l, neg2_u, neg2_c, neg2_s,
393           mcf3.OPTIMAL, true,     -300, test_str + "-18", GEQ);
394}
395
396template < typename MCF, typename Param >
397void runMcfLeqTests( Param param,
398                     const std::string &test_str = "" )
399{
400  // Tests for the LEQ form
401  MCF mcf1(gr);
402  mcf1.supplyType(mcf1.LEQ);
403  mcf1.upperMap(u).costMap(c).supplyMap(s6);
404  checkMcf(mcf1, mcf1.run(param), gr, l1, u, c, s6,
405           mcf1.OPTIMAL, true,   5080, test_str + "-19", LEQ);
406  mcf1.lowerMap(l2);
407  checkMcf(mcf1, mcf1.run(param), gr, l2, u, c, s6,
408           mcf1.OPTIMAL, true,   5930, test_str + "-20", LEQ);
409  mcf1.supplyMap(s5);
410  checkMcf(mcf1, mcf1.run(param), gr, l2, u, c, s5,
411           mcf1.INFEASIBLE, false,  0, test_str + "-21", LEQ);
412}
413
414
415int main()
416{
417  // Read the test networks
418  std::istringstream input(test_lgf);
419  DigraphReader<Digraph>(gr, input)
420    .arcMap("cost", c)
421    .arcMap("cap", u)
422    .arcMap("low1", l1)
423    .arcMap("low2", l2)
424    .arcMap("low3", l3)
425    .nodeMap("sup1", s1)
426    .nodeMap("sup2", s2)
427    .nodeMap("sup3", s3)
428    .nodeMap("sup4", s4)
429    .nodeMap("sup5", s5)
430    .nodeMap("sup6", s6)
431    .node("source", v)
432    .node("target", w)
433    .run();
434 
435  std::istringstream neg_inp1(test_neg1_lgf);
436  DigraphReader<Digraph>(neg1_gr, neg_inp1)
437    .arcMap("cost", neg1_c)
438    .arcMap("low1", neg1_l1)
439    .arcMap("low2", neg1_l2)
440    .nodeMap("sup", neg1_s)
441    .run();
442 
443  std::istringstream neg_inp2(test_neg2_lgf);
444  DigraphReader<Digraph>(neg2_gr, neg_inp2)
445    .arcMap("cost", neg2_c)
446    .nodeMap("sup", neg2_s)
447    .run();
448 
449  // Check the interface of NetworkSimplex
450  {
451    typedef concepts::Digraph GR;
452    checkConcept< McfClassConcept<GR, int, int>,
453                  NetworkSimplex<GR> >();
454    checkConcept< McfClassConcept<GR, double, double>,
455                  NetworkSimplex<GR, double> >();
456    checkConcept< McfClassConcept<GR, int, double>,
457                  NetworkSimplex<GR, int, double> >();
458  }
459
460  // Test NetworkSimplex
461  {
462    typedef NetworkSimplex<Digraph> MCF;
463    runMcfGeqTests<MCF>(MCF::FIRST_ELIGIBLE, "NS-FE", true);
464    runMcfLeqTests<MCF>(MCF::FIRST_ELIGIBLE, "NS-FE");
465    runMcfGeqTests<MCF>(MCF::BEST_ELIGIBLE,  "NS-BE", true);
466    runMcfLeqTests<MCF>(MCF::BEST_ELIGIBLE,  "NS-BE");
467    runMcfGeqTests<MCF>(MCF::BLOCK_SEARCH,   "NS-BS", true);
468    runMcfLeqTests<MCF>(MCF::BLOCK_SEARCH,   "NS-BS");
469    runMcfGeqTests<MCF>(MCF::CANDIDATE_LIST, "NS-CL", true);
470    runMcfLeqTests<MCF>(MCF::CANDIDATE_LIST, "NS-CL");
471    runMcfGeqTests<MCF>(MCF::ALTERING_LIST,  "NS-AL", true);
472    runMcfLeqTests<MCF>(MCF::ALTERING_LIST,  "NS-AL");
473  }
474
475  return 0;
476}
Note: See TracBrowser for help on using the repository browser.