test/min_cost_flow_test.cc
changeset 687 6c408d864fa1
parent 662 e3d9bff447ed
child 689 111698359429
equal deleted inserted replaced
5:107c8e3f7229 6:a0f49214daf0
    16  *
    16  *
    17  */
    17  */
    18 
    18 
    19 #include <iostream>
    19 #include <iostream>
    20 #include <fstream>
    20 #include <fstream>
       
    21 #include <limits>
    21 
    22 
    22 #include <lemon/list_graph.h>
    23 #include <lemon/list_graph.h>
    23 #include <lemon/lgf_reader.h>
    24 #include <lemon/lgf_reader.h>
    24 
    25 
    25 #include <lemon/network_simplex.h>
    26 #include <lemon/network_simplex.h>
    31 
    32 
    32 using namespace lemon;
    33 using namespace lemon;
    33 
    34 
    34 char test_lgf[] =
    35 char test_lgf[] =
    35   "@nodes\n"
    36   "@nodes\n"
    36   "label  sup1 sup2 sup3 sup4 sup5\n"
    37   "label  sup1 sup2 sup3 sup4 sup5 sup6\n"
    37   "    1    20   27    0   20   30\n"
    38   "    1    20   27    0   30   20   30\n"
    38   "    2    -4    0    0   -8   -3\n"
    39   "    2    -4    0    0    0   -8   -3\n"
    39   "    3     0    0    0    0    0\n"
    40   "    3     0    0    0    0    0    0\n"
    40   "    4     0    0    0    0    0\n"
    41   "    4     0    0    0    0    0    0\n"
    41   "    5     9    0    0    6   11\n"
    42   "    5     9    0    0    0    6   11\n"
    42   "    6    -6    0    0   -5   -6\n"
    43   "    6    -6    0    0    0   -5   -6\n"
    43   "    7     0    0    0    0    0\n"
    44   "    7     0    0    0    0    0    0\n"
    44   "    8     0    0    0    0    3\n"
    45   "    8     0    0    0    0    0    3\n"
    45   "    9     3    0    0    0    0\n"
    46   "    9     3    0    0    0    0    0\n"
    46   "   10    -2    0    0   -7   -2\n"
    47   "   10    -2    0    0    0   -7   -2\n"
    47   "   11     0    0    0  -10    0\n"
    48   "   11     0    0    0    0  -10    0\n"
    48   "   12   -20  -27    0  -30  -20\n"
    49   "   12   -20  -27    0  -30  -30  -20\n"
    49   "\n"
    50   "\n"                
    50   "@arcs\n"
    51   "@arcs\n"
    51   "       cost  cap low1 low2\n"
    52   "       cost  cap low1 low2 low3\n"
    52   " 1  2    70   11    0    8\n"
    53   " 1  2    70   11    0    8    8\n"
    53   " 1  3   150    3    0    1\n"
    54   " 1  3   150    3    0    1    0\n"
    54   " 1  4    80   15    0    2\n"
    55   " 1  4    80   15    0    2    2\n"
    55   " 2  8    80   12    0    0\n"
    56   " 2  8    80   12    0    0    0\n"
    56   " 3  5   140    5    0    3\n"
    57   " 3  5   140    5    0    3    1\n"
    57   " 4  6    60   10    0    1\n"
    58   " 4  6    60   10    0    1    0\n"
    58   " 4  7    80    2    0    0\n"
    59   " 4  7    80    2    0    0    0\n"
    59   " 4  8   110    3    0    0\n"
    60   " 4  8   110    3    0    0    0\n"
    60   " 5  7    60   14    0    0\n"
    61   " 5  7    60   14    0    0    0\n"
    61   " 5 11   120   12    0    0\n"
    62   " 5 11   120   12    0    0    0\n"
    62   " 6  3     0    3    0    0\n"
    63   " 6  3     0    3    0    0    0\n"
    63   " 6  9   140    4    0    0\n"
    64   " 6  9   140    4    0    0    0\n"
    64   " 6 10    90    8    0    0\n"
    65   " 6 10    90    8    0    0    0\n"
    65   " 7  1    30    5    0    0\n"
    66   " 7  1    30    5    0    0   -5\n"
    66   " 8 12    60   16    0    4\n"
    67   " 8 12    60   16    0    4    3\n"
    67   " 9 12    50    6    0    0\n"
    68   " 9 12    50    6    0    0    0\n"
    68   "10 12    70   13    0    5\n"
    69   "10 12    70   13    0    5    2\n"
    69   "10  2   100    7    0    0\n"
    70   "10  2   100    7    0    0    0\n"
    70   "10  7    60   10    0    0\n"
    71   "10  7    60   10    0    0   -3\n"
    71   "11 10    20   14    0    6\n"
    72   "11 10    20   14    0    6  -20\n"
    72   "12 11    30   10    0    0\n"
    73   "12 11    30   10    0    0  -10\n"
    73   "\n"
    74   "\n"
    74   "@attributes\n"
    75   "@attributes\n"
    75   "source 1\n"
    76   "source 1\n"
    76   "target 12\n";
    77   "target 12\n";
    77 
    78 
    78 
    79 
    79 enum ProblemType {
    80 enum SupplyType {
    80   EQ,
    81   EQ,
    81   GEQ,
    82   GEQ,
    82   LEQ
    83   LEQ
    83 };
    84 };
    84 
    85 
    96       MCF mcf(g);
    97       MCF mcf(g);
    97 
    98 
    98       b = mcf.reset()
    99       b = mcf.reset()
    99              .lowerMap(lower)
   100              .lowerMap(lower)
   100              .upperMap(upper)
   101              .upperMap(upper)
   101              .capacityMap(upper)
       
   102              .boundMaps(lower, upper)
       
   103              .costMap(cost)
   102              .costMap(cost)
   104              .supplyMap(sup)
   103              .supplyMap(sup)
   105              .stSupply(n, n, k)
   104              .stSupply(n, n, k)
   106              .flowMap(flow)
   105              .flowMap(flow)
   107              .potentialMap(pot)
   106              .potentialMap(pot)
   110       const MCF& const_mcf = mcf;
   109       const MCF& const_mcf = mcf;
   111 
   110 
   112       const typename MCF::FlowMap &fm = const_mcf.flowMap();
   111       const typename MCF::FlowMap &fm = const_mcf.flowMap();
   113       const typename MCF::PotentialMap &pm = const_mcf.potentialMap();
   112       const typename MCF::PotentialMap &pm = const_mcf.potentialMap();
   114 
   113 
   115       v = const_mcf.totalCost();
   114       c = const_mcf.totalCost();
   116       double x = const_mcf.template totalCost<double>();
   115       double x = const_mcf.template totalCost<double>();
   117       v = const_mcf.flow(a);
   116       v = const_mcf.flow(a);
   118       v = const_mcf.potential(n);
   117       c = const_mcf.potential(n);
       
   118       
       
   119       v = const_mcf.INF;
   119 
   120 
   120       ignore_unused_variable_warning(fm);
   121       ignore_unused_variable_warning(fm);
   121       ignore_unused_variable_warning(pm);
   122       ignore_unused_variable_warning(pm);
   122       ignore_unused_variable_warning(x);
   123       ignore_unused_variable_warning(x);
   123     }
   124     }
   135     const NM &sup;
   136     const NM &sup;
   136     const Node &n;
   137     const Node &n;
   137     const Arc &a;
   138     const Arc &a;
   138     const Flow &k;
   139     const Flow &k;
   139     Flow v;
   140     Flow v;
       
   141     Cost c;
   140     bool b;
   142     bool b;
   141 
   143 
   142     typename MCF::FlowMap &flow;
   144     typename MCF::FlowMap &flow;
   143     typename MCF::PotentialMap &pot;
   145     typename MCF::PotentialMap &pot;
   144   };
   146   };
   149 // Check the feasibility of the given flow (primal soluiton)
   151 // Check the feasibility of the given flow (primal soluiton)
   150 template < typename GR, typename LM, typename UM,
   152 template < typename GR, typename LM, typename UM,
   151            typename SM, typename FM >
   153            typename SM, typename FM >
   152 bool checkFlow( const GR& gr, const LM& lower, const UM& upper,
   154 bool checkFlow( const GR& gr, const LM& lower, const UM& upper,
   153                 const SM& supply, const FM& flow,
   155                 const SM& supply, const FM& flow,
   154                 ProblemType type = EQ )
   156                 SupplyType type = EQ )
   155 {
   157 {
   156   TEMPLATE_DIGRAPH_TYPEDEFS(GR);
   158   TEMPLATE_DIGRAPH_TYPEDEFS(GR);
   157 
   159 
   158   for (ArcIt e(gr); e != INVALID; ++e) {
   160   for (ArcIt e(gr); e != INVALID; ++e) {
   159     if (flow[e] < lower[e] || flow[e] > upper[e]) return false;
   161     if (flow[e] < lower[e] || flow[e] > upper[e]) return false;
   206 }
   208 }
   207 
   209 
   208 // Run a minimum cost flow algorithm and check the results
   210 // Run a minimum cost flow algorithm and check the results
   209 template < typename MCF, typename GR,
   211 template < typename MCF, typename GR,
   210            typename LM, typename UM,
   212            typename LM, typename UM,
   211            typename CM, typename SM >
   213            typename CM, typename SM,
   212 void checkMcf( const MCF& mcf, bool mcf_result,
   214            typename PT >
       
   215 void checkMcf( const MCF& mcf, PT mcf_result,
   213                const GR& gr, const LM& lower, const UM& upper,
   216                const GR& gr, const LM& lower, const UM& upper,
   214                const CM& cost, const SM& supply,
   217                const CM& cost, const SM& supply,
   215                bool result, typename CM::Value total,
   218                PT result, bool optimal, typename CM::Value total,
   216                const std::string &test_id = "",
   219                const std::string &test_id = "",
   217                ProblemType type = EQ )
   220                SupplyType type = EQ )
   218 {
   221 {
   219   check(mcf_result == result, "Wrong result " + test_id);
   222   check(mcf_result == result, "Wrong result " + test_id);
   220   if (result) {
   223   if (optimal) {
   221     check(checkFlow(gr, lower, upper, supply, mcf.flowMap(), type),
   224     check(checkFlow(gr, lower, upper, supply, mcf.flowMap(), type),
   222           "The flow is not feasible " + test_id);
   225           "The flow is not feasible " + test_id);
   223     check(mcf.totalCost() == total, "The flow is not optimal " + test_id);
   226     check(mcf.totalCost() == total, "The flow is not optimal " + test_id);
   224     check(checkPotential(gr, lower, upper, cost, supply, mcf.flowMap(),
   227     check(checkPotential(gr, lower, upper, cost, supply, mcf.flowMap(),
   225                          mcf.potentialMap()),
   228                          mcf.potentialMap()),
   242   typedef ListDigraph Digraph;
   245   typedef ListDigraph Digraph;
   243   DIGRAPH_TYPEDEFS(ListDigraph);
   246   DIGRAPH_TYPEDEFS(ListDigraph);
   244 
   247 
   245   // Read the test digraph
   248   // Read the test digraph
   246   Digraph gr;
   249   Digraph gr;
   247   Digraph::ArcMap<int> c(gr), l1(gr), l2(gr), u(gr);
   250   Digraph::ArcMap<int> c(gr), l1(gr), l2(gr), l3(gr), u(gr);
   248   Digraph::NodeMap<int> s1(gr), s2(gr), s3(gr), s4(gr), s5(gr);
   251   Digraph::NodeMap<int> s1(gr), s2(gr), s3(gr), s4(gr), s5(gr), s6(gr);
   249   ConstMap<Arc, int> cc(1), cu(std::numeric_limits<int>::max());
   252   ConstMap<Arc, int> cc(1), cu(std::numeric_limits<int>::max());
   250   Node v, w;
   253   Node v, w;
   251 
   254 
   252   std::istringstream input(test_lgf);
   255   std::istringstream input(test_lgf);
   253   DigraphReader<Digraph>(gr, input)
   256   DigraphReader<Digraph>(gr, input)
   254     .arcMap("cost", c)
   257     .arcMap("cost", c)
   255     .arcMap("cap", u)
   258     .arcMap("cap", u)
   256     .arcMap("low1", l1)
   259     .arcMap("low1", l1)
   257     .arcMap("low2", l2)
   260     .arcMap("low2", l2)
       
   261     .arcMap("low3", l3)
   258     .nodeMap("sup1", s1)
   262     .nodeMap("sup1", s1)
   259     .nodeMap("sup2", s2)
   263     .nodeMap("sup2", s2)
   260     .nodeMap("sup3", s3)
   264     .nodeMap("sup3", s3)
   261     .nodeMap("sup4", s4)
   265     .nodeMap("sup4", s4)
   262     .nodeMap("sup5", s5)
   266     .nodeMap("sup5", s5)
       
   267     .nodeMap("sup6", s6)
   263     .node("source", v)
   268     .node("source", v)
   264     .node("target", w)
   269     .node("target", w)
   265     .run();
   270     .run();
       
   271   
       
   272   // Build a test digraph for testing negative costs
       
   273   Digraph ngr;
       
   274   Node n1 = ngr.addNode();
       
   275   Node n2 = ngr.addNode();
       
   276   Node n3 = ngr.addNode();
       
   277   Node n4 = ngr.addNode();
       
   278   Node n5 = ngr.addNode();
       
   279   Node n6 = ngr.addNode();
       
   280   Node n7 = ngr.addNode();
       
   281   
       
   282   Arc a1 = ngr.addArc(n1, n2);
       
   283   Arc a2 = ngr.addArc(n1, n3);
       
   284   Arc a3 = ngr.addArc(n2, n4);
       
   285   Arc a4 = ngr.addArc(n3, n4);
       
   286   Arc a5 = ngr.addArc(n3, n2);
       
   287   Arc a6 = ngr.addArc(n5, n3);
       
   288   Arc a7 = ngr.addArc(n5, n6);
       
   289   Arc a8 = ngr.addArc(n6, n7);
       
   290   Arc a9 = ngr.addArc(n7, n5);
       
   291   
       
   292   Digraph::ArcMap<int> nc(ngr), nl1(ngr, 0), nl2(ngr, 0);
       
   293   ConstMap<Arc, int> nu1(std::numeric_limits<int>::max()), nu2(5000);
       
   294   Digraph::NodeMap<int> ns(ngr, 0);
       
   295   
       
   296   nl2[a7] =  1000;
       
   297   nl2[a8] = -1000;
       
   298   
       
   299   ns[n1] =  100;
       
   300   ns[n4] = -100;
       
   301   
       
   302   nc[a1] =  100;
       
   303   nc[a2] =   30;
       
   304   nc[a3] =   20;
       
   305   nc[a4] =   80;
       
   306   nc[a5] =   50;
       
   307   nc[a6] =   10;
       
   308   nc[a7] =   80;
       
   309   nc[a8] =   30;
       
   310   nc[a9] = -120;
   266 
   311 
   267   // A. Test NetworkSimplex with the default pivot rule
   312   // A. Test NetworkSimplex with the default pivot rule
   268   {
   313   {
   269     NetworkSimplex<Digraph> mcf(gr);
   314     NetworkSimplex<Digraph> mcf(gr);
   270 
   315 
   271     // Check the equality form
   316     // Check the equality form
   272     mcf.upperMap(u).costMap(c);
   317     mcf.upperMap(u).costMap(c);
   273     checkMcf(mcf, mcf.supplyMap(s1).run(),
   318     checkMcf(mcf, mcf.supplyMap(s1).run(),
   274              gr, l1, u, c, s1, true,  5240, "#A1");
   319              gr, l1, u, c, s1, mcf.OPTIMAL, true,   5240, "#A1");
   275     checkMcf(mcf, mcf.stSupply(v, w, 27).run(),
   320     checkMcf(mcf, mcf.stSupply(v, w, 27).run(),
   276              gr, l1, u, c, s2, true,  7620, "#A2");
   321              gr, l1, u, c, s2, mcf.OPTIMAL, true,   7620, "#A2");
   277     mcf.lowerMap(l2);
   322     mcf.lowerMap(l2);
   278     checkMcf(mcf, mcf.supplyMap(s1).run(),
   323     checkMcf(mcf, mcf.supplyMap(s1).run(),
   279              gr, l2, u, c, s1, true,  5970, "#A3");
   324              gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#A3");
   280     checkMcf(mcf, mcf.stSupply(v, w, 27).run(),
   325     checkMcf(mcf, mcf.stSupply(v, w, 27).run(),
   281              gr, l2, u, c, s2, true,  8010, "#A4");
   326              gr, l2, u, c, s2, mcf.OPTIMAL, true,   8010, "#A4");
   282     mcf.reset();
   327     mcf.reset();
   283     checkMcf(mcf, mcf.supplyMap(s1).run(),
   328     checkMcf(mcf, mcf.supplyMap(s1).run(),
   284              gr, l1, cu, cc, s1, true,  74, "#A5");
   329              gr, l1, cu, cc, s1, mcf.OPTIMAL, true,   74, "#A5");
   285     checkMcf(mcf, mcf.lowerMap(l2).stSupply(v, w, 27).run(),
   330     checkMcf(mcf, mcf.lowerMap(l2).stSupply(v, w, 27).run(),
   286              gr, l2, cu, cc, s2, true,  94, "#A6");
   331              gr, l2, cu, cc, s2, mcf.OPTIMAL, true,   94, "#A6");
   287     mcf.reset();
   332     mcf.reset();
   288     checkMcf(mcf, mcf.run(),
   333     checkMcf(mcf, mcf.run(),
   289              gr, l1, cu, cc, s3, true,   0, "#A7");
   334              gr, l1, cu, cc, s3, mcf.OPTIMAL, true,    0, "#A7");
   290     checkMcf(mcf, mcf.boundMaps(l2, u).run(),
   335     checkMcf(mcf, mcf.lowerMap(l2).upperMap(u).run(),
   291              gr, l2, u, cc, s3, false,   0, "#A8");
   336              gr, l2, u, cc, s3, mcf.INFEASIBLE, false, 0, "#A8");
       
   337     mcf.reset().lowerMap(l3).upperMap(u).costMap(c).supplyMap(s4);
       
   338     checkMcf(mcf, mcf.run(),
       
   339              gr, l3, u, c, s4, mcf.OPTIMAL, true,   6360, "#A9");
   292 
   340 
   293     // Check the GEQ form
   341     // Check the GEQ form
   294     mcf.reset().upperMap(u).costMap(c).supplyMap(s4);
   342     mcf.reset().upperMap(u).costMap(c).supplyMap(s5);
   295     checkMcf(mcf, mcf.run(),
   343     checkMcf(mcf, mcf.run(),
   296              gr, l1, u, c, s4, true,  3530, "#A9", GEQ);
   344              gr, l1, u, c, s5, mcf.OPTIMAL, true,   3530, "#A10", GEQ);
   297     mcf.problemType(mcf.GEQ);
   345     mcf.supplyType(mcf.GEQ);
   298     checkMcf(mcf, mcf.lowerMap(l2).run(),
   346     checkMcf(mcf, mcf.lowerMap(l2).run(),
   299              gr, l2, u, c, s4, true,  4540, "#A10", GEQ);
   347              gr, l2, u, c, s5, mcf.OPTIMAL, true,   4540, "#A11", GEQ);
   300     mcf.problemType(mcf.CARRY_SUPPLIES).supplyMap(s5);
   348     mcf.supplyType(mcf.CARRY_SUPPLIES).supplyMap(s6);
   301     checkMcf(mcf, mcf.run(),
   349     checkMcf(mcf, mcf.run(),
   302              gr, l2, u, c, s5, false,    0, "#A11", GEQ);
   350              gr, l2, u, c, s6, mcf.INFEASIBLE, false,  0, "#A12", GEQ);
   303 
   351 
   304     // Check the LEQ form
   352     // Check the LEQ form
   305     mcf.reset().problemType(mcf.LEQ);
   353     mcf.reset().supplyType(mcf.LEQ);
   306     mcf.upperMap(u).costMap(c).supplyMap(s5);
   354     mcf.upperMap(u).costMap(c).supplyMap(s6);
   307     checkMcf(mcf, mcf.run(),
   355     checkMcf(mcf, mcf.run(),
   308              gr, l1, u, c, s5, true,  5080, "#A12", LEQ);
   356              gr, l1, u, c, s6, mcf.OPTIMAL, true,   5080, "#A13", LEQ);
   309     checkMcf(mcf, mcf.lowerMap(l2).run(),
   357     checkMcf(mcf, mcf.lowerMap(l2).run(),
   310              gr, l2, u, c, s5, true,  5930, "#A13", LEQ);
   358              gr, l2, u, c, s6, mcf.OPTIMAL, true,   5930, "#A14", LEQ);
   311     mcf.problemType(mcf.SATISFY_DEMANDS).supplyMap(s4);
   359     mcf.supplyType(mcf.SATISFY_DEMANDS).supplyMap(s5);
   312     checkMcf(mcf, mcf.run(),
   360     checkMcf(mcf, mcf.run(),
   313              gr, l2, u, c, s4, false,    0, "#A14", LEQ);
   361              gr, l2, u, c, s5, mcf.INFEASIBLE, false,  0, "#A15", LEQ);
       
   362 
       
   363     // Check negative costs
       
   364     NetworkSimplex<Digraph> nmcf(ngr);
       
   365     nmcf.lowerMap(nl1).costMap(nc).supplyMap(ns);
       
   366     checkMcf(nmcf, nmcf.run(),
       
   367       ngr, nl1, nu1, nc, ns, nmcf.UNBOUNDED, false,    0, "#A16");
       
   368     checkMcf(nmcf, nmcf.upperMap(nu2).run(),
       
   369       ngr, nl1, nu2, nc, ns, nmcf.OPTIMAL, true,  -40000, "#A17");
       
   370     nmcf.reset().lowerMap(nl2).costMap(nc).supplyMap(ns);
       
   371     checkMcf(nmcf, nmcf.run(),
       
   372       ngr, nl2, nu1, nc, ns, nmcf.UNBOUNDED, false,    0, "#A18");
   314   }
   373   }
   315 
   374 
   316   // B. Test NetworkSimplex with each pivot rule
   375   // B. Test NetworkSimplex with each pivot rule
   317   {
   376   {
   318     NetworkSimplex<Digraph> mcf(gr);
   377     NetworkSimplex<Digraph> mcf(gr);
   319     mcf.supplyMap(s1).costMap(c).capacityMap(u).lowerMap(l2);
   378     mcf.supplyMap(s1).costMap(c).upperMap(u).lowerMap(l2);
   320 
   379 
   321     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::FIRST_ELIGIBLE),
   380     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::FIRST_ELIGIBLE),
   322              gr, l2, u, c, s1, true,  5970, "#B1");
   381              gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B1");
   323     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::BEST_ELIGIBLE),
   382     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::BEST_ELIGIBLE),
   324              gr, l2, u, c, s1, true,  5970, "#B2");
   383              gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B2");
   325     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::BLOCK_SEARCH),
   384     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::BLOCK_SEARCH),
   326              gr, l2, u, c, s1, true,  5970, "#B3");
   385              gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B3");
   327     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::CANDIDATE_LIST),
   386     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::CANDIDATE_LIST),
   328              gr, l2, u, c, s1, true,  5970, "#B4");
   387              gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B4");
   329     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::ALTERING_LIST),
   388     checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::ALTERING_LIST),
   330              gr, l2, u, c, s1, true,  5970, "#B5");
   389              gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B5");
   331   }
   390   }
   332 
   391 
   333   return 0;
   392   return 0;
   334 }
   393 }