gravatar
kpeter (Peter Kovacs)
kpeter@inf.elte.hu
Extend min cost flow test file + check dual costs (#291)
0 1 0
default
1 file changed with 105 insertions and 48 deletions:
↑ Collapse diff ↑
Ignore white space 6 line context
... ...
@@ -174,7 +174,7 @@
174 174
           typename CM, typename SM, typename FM, typename PM >
175 175
bool checkPotential( const GR& gr, const LM& lower, const UM& upper,
176 176
                     const CM& cost, const SM& supply, const FM& flow, 
177
                     const PM& pi )
177
                     const PM& pi, SupplyType type )
178 178
{
179 179
  TEMPLATE_DIGRAPH_TYPEDEFS(GR);
180 180

	
... ...
@@ -193,12 +193,50 @@
193 193
      sum += flow[e];
194 194
    for (InArcIt e(gr, n); e != INVALID; ++e)
195 195
      sum -= flow[e];
196
    opt = (sum == supply[n]) || (pi[n] == 0);
196
    if (type != LEQ) {
197
      opt = (pi[n] <= 0) && (sum == supply[n] || pi[n] == 0);
198
    } else {
199
      opt = (pi[n] >= 0) && (sum == supply[n] || pi[n] == 0);
200
    }
197 201
  }
198 202
  
199 203
  return opt;
200 204
}
201 205

	
206
// Check whether the dual cost is equal to the primal cost
207
template < typename GR, typename LM, typename UM,
208
           typename CM, typename SM, typename PM >
209
bool checkDualCost( const GR& gr, const LM& lower, const UM& upper,
210
                    const CM& cost, const SM& supply, const PM& pi,
211
                    typename CM::Value total )
212
{
213
  TEMPLATE_DIGRAPH_TYPEDEFS(GR);
214

	
215
  typename CM::Value dual_cost = 0;
216
  SM red_supply(gr);
217
  for (NodeIt n(gr); n != INVALID; ++n) {
218
    red_supply[n] = supply[n];
219
  }
220
  for (ArcIt a(gr); a != INVALID; ++a) {
221
    if (lower[a] != 0) {
222
      dual_cost += lower[a] * cost[a];
223
      red_supply[gr.source(a)] -= lower[a];
224
      red_supply[gr.target(a)] += lower[a];
225
    }
226
  }
227
  
228
  for (NodeIt n(gr); n != INVALID; ++n) {
229
    dual_cost -= red_supply[n] * pi[n];
230
  }
231
  for (ArcIt a(gr); a != INVALID; ++a) {
232
    typename CM::Value red_cost =
233
      cost[a] + pi[gr.source(a)] - pi[gr.target(a)];
234
    dual_cost -= (upper[a] - lower[a]) * std::max(-red_cost, 0);
235
  }
236
  
237
  return dual_cost == total;
238
}
239

	
202 240
// Run a minimum cost flow algorithm and check the results
203 241
template < typename MCF, typename GR,
204 242
           typename LM, typename UM,
... ...
@@ -220,8 +258,10 @@
220 258
    check(checkFlow(gr, lower, upper, supply, flow, type),
221 259
          "The flow is not feasible " + test_id);
222 260
    check(mcf.totalCost() == total, "The flow is not optimal " + test_id);
223
    check(checkPotential(gr, lower, upper, cost, supply, flow, pi),
261
    check(checkPotential(gr, lower, upper, cost, supply, flow, pi, type),
224 262
          "Wrong potentials " + test_id);
263
    check(checkDualCost(gr, lower, upper, cost, supply, pi, total),
264
          "Wrong dual cost " + test_id);
225 265
  }
226 266
}
227 267

	
... ...
@@ -266,45 +306,56 @@
266 306
    .node("target", w)
267 307
    .run();
268 308
  
269
  // Build a test digraph for testing negative costs
270
  Digraph ngr;
271
  Node n1 = ngr.addNode();
272
  Node n2 = ngr.addNode();
273
  Node n3 = ngr.addNode();
274
  Node n4 = ngr.addNode();
275
  Node n5 = ngr.addNode();
276
  Node n6 = ngr.addNode();
277
  Node n7 = ngr.addNode();
309
  // Build test digraphs with negative costs
310
  Digraph neg_gr;
311
  Node n1 = neg_gr.addNode();
312
  Node n2 = neg_gr.addNode();
313
  Node n3 = neg_gr.addNode();
314
  Node n4 = neg_gr.addNode();
315
  Node n5 = neg_gr.addNode();
316
  Node n6 = neg_gr.addNode();
317
  Node n7 = neg_gr.addNode();
278 318
  
279
  Arc a1 = ngr.addArc(n1, n2);
280
  Arc a2 = ngr.addArc(n1, n3);
281
  Arc a3 = ngr.addArc(n2, n4);
282
  Arc a4 = ngr.addArc(n3, n4);
283
  Arc a5 = ngr.addArc(n3, n2);
284
  Arc a6 = ngr.addArc(n5, n3);
285
  Arc a7 = ngr.addArc(n5, n6);
286
  Arc a8 = ngr.addArc(n6, n7);
287
  Arc a9 = ngr.addArc(n7, n5);
319
  Arc a1 = neg_gr.addArc(n1, n2);
320
  Arc a2 = neg_gr.addArc(n1, n3);
321
  Arc a3 = neg_gr.addArc(n2, n4);
322
  Arc a4 = neg_gr.addArc(n3, n4);
323
  Arc a5 = neg_gr.addArc(n3, n2);
324
  Arc a6 = neg_gr.addArc(n5, n3);
325
  Arc a7 = neg_gr.addArc(n5, n6);
326
  Arc a8 = neg_gr.addArc(n6, n7);
327
  Arc a9 = neg_gr.addArc(n7, n5);
288 328
  
289
  Digraph::ArcMap<int> nc(ngr), nl1(ngr, 0), nl2(ngr, 0);
290
  ConstMap<Arc, int> nu1(std::numeric_limits<int>::max()), nu2(5000);
291
  Digraph::NodeMap<int> ns(ngr, 0);
329
  Digraph::ArcMap<int> neg_c(neg_gr), neg_l1(neg_gr, 0), neg_l2(neg_gr, 0);
330
  ConstMap<Arc, int> neg_u1(std::numeric_limits<int>::max()), neg_u2(5000);
331
  Digraph::NodeMap<int> neg_s(neg_gr, 0);
292 332
  
293
  nl2[a7] =  1000;
294
  nl2[a8] = -1000;
333
  neg_l2[a7] =  1000;
334
  neg_l2[a8] = -1000;
295 335
  
296
  ns[n1] =  100;
297
  ns[n4] = -100;
336
  neg_s[n1] =  100;
337
  neg_s[n4] = -100;
298 338
  
299
  nc[a1] =  100;
300
  nc[a2] =   30;
301
  nc[a3] =   20;
302
  nc[a4] =   80;
303
  nc[a5] =   50;
304
  nc[a6] =   10;
305
  nc[a7] =   80;
306
  nc[a8] =   30;
307
  nc[a9] = -120;
339
  neg_c[a1] =  100;
340
  neg_c[a2] =   30;
341
  neg_c[a3] =   20;
342
  neg_c[a4] =   80;
343
  neg_c[a5] =   50;
344
  neg_c[a6] =   10;
345
  neg_c[a7] =   80;
346
  neg_c[a8] =   30;
347
  neg_c[a9] = -120;
348

	
349
  Digraph negs_gr;
350
  Digraph::NodeMap<int> negs_s(negs_gr);
351
  Digraph::ArcMap<int> negs_c(negs_gr);
352
  ConstMap<Arc, int> negs_l(0), negs_u(1000);
353
  n1 = negs_gr.addNode();
354
  n2 = negs_gr.addNode();
355
  negs_s[n1] = 100;
356
  negs_s[n2] = -300;
357
  negs_c[negs_gr.addArc(n1, n2)] = -1;
358

	
308 359

	
309 360
  // A. Test NetworkSimplex with the default pivot rule
310 361
  {
... ...
@@ -342,7 +393,7 @@
342 393
    mcf.supplyType(mcf.GEQ);
343 394
    checkMcf(mcf, mcf.lowerMap(l2).run(),
344 395
             gr, l2, u, c, s5, mcf.OPTIMAL, true,   4540, "#A11", GEQ);
345
    mcf.supplyType(mcf.CARRY_SUPPLIES).supplyMap(s6);
396
    mcf.supplyMap(s6);
346 397
    checkMcf(mcf, mcf.run(),
347 398
             gr, l2, u, c, s6, mcf.INFEASIBLE, false,  0, "#A12", GEQ);
348 399

	
... ...
@@ -353,20 +404,26 @@
353 404
             gr, l1, u, c, s6, mcf.OPTIMAL, true,   5080, "#A13", LEQ);
354 405
    checkMcf(mcf, mcf.lowerMap(l2).run(),
355 406
             gr, l2, u, c, s6, mcf.OPTIMAL, true,   5930, "#A14", LEQ);
356
    mcf.supplyType(mcf.SATISFY_DEMANDS).supplyMap(s5);
407
    mcf.supplyMap(s5);
357 408
    checkMcf(mcf, mcf.run(),
358 409
             gr, l2, u, c, s5, mcf.INFEASIBLE, false,  0, "#A15", LEQ);
359 410

	
360 411
    // Check negative costs
361
    NetworkSimplex<Digraph> nmcf(ngr);
362
    nmcf.lowerMap(nl1).costMap(nc).supplyMap(ns);
363
    checkMcf(nmcf, nmcf.run(),
364
      ngr, nl1, nu1, nc, ns, nmcf.UNBOUNDED, false,    0, "#A16");
365
    checkMcf(nmcf, nmcf.upperMap(nu2).run(),
366
      ngr, nl1, nu2, nc, ns, nmcf.OPTIMAL, true,  -40000, "#A17");
367
    nmcf.reset().lowerMap(nl2).costMap(nc).supplyMap(ns);
368
    checkMcf(nmcf, nmcf.run(),
369
      ngr, nl2, nu1, nc, ns, nmcf.UNBOUNDED, false,    0, "#A18");
412
    NetworkSimplex<Digraph> neg_mcf(neg_gr);
413
    neg_mcf.lowerMap(neg_l1).costMap(neg_c).supplyMap(neg_s);
414
    checkMcf(neg_mcf, neg_mcf.run(), neg_gr, neg_l1, neg_u1,
415
      neg_c, neg_s, neg_mcf.UNBOUNDED, false,    0, "#A16");
416
    neg_mcf.upperMap(neg_u2);
417
    checkMcf(neg_mcf, neg_mcf.run(), neg_gr, neg_l1, neg_u2,
418
      neg_c, neg_s, neg_mcf.OPTIMAL, true,  -40000, "#A17");
419
    neg_mcf.reset().lowerMap(neg_l2).costMap(neg_c).supplyMap(neg_s);
420
    checkMcf(neg_mcf, neg_mcf.run(), neg_gr, neg_l2, neg_u1,
421
      neg_c, neg_s, neg_mcf.UNBOUNDED, false,    0, "#A18");
422
      
423
    NetworkSimplex<Digraph> negs_mcf(negs_gr);
424
    negs_mcf.costMap(negs_c).supplyMap(negs_s);
425
    checkMcf(negs_mcf, negs_mcf.run(), negs_gr, negs_l, negs_u,
426
      negs_c, negs_s, negs_mcf.OPTIMAL, true, -300, "#A19", GEQ);
370 427
  }
371 428

	
372 429
  // B. Test NetworkSimplex with each pivot rule
0 comments (0 inline)