gravatar
deba@inf.elte.hu
deba@inf.elte.hu
Fix GLPK tests (#213)
0 4 0
default
4 files changed with 14 insertions and 0 deletions:
↑ Collapse diff ↑
Ignore white space 6 line context
... ...
@@ -141,768 +141,772 @@
141 141
    }
142 142
  }
143 143

	
144 144
  void GlpkBase::_setColCoeffs(int ix, ExprIterator b,
145 145
                                     ExprIterator e) {
146 146

	
147 147
    std::vector<int> indexes;
148 148
    std::vector<Value> values;
149 149

	
150 150
    indexes.push_back(0);
151 151
    values.push_back(0);
152 152

	
153 153
    for(ExprIterator it = b; it != e; ++it) {
154 154
      indexes.push_back(it->first);
155 155
      values.push_back(it->second);
156 156
    }
157 157

	
158 158
    glp_set_mat_col(lp, ix, values.size() - 1,
159 159
                    &indexes.front(), &values.front());
160 160
  }
161 161

	
162 162
  void GlpkBase::_getColCoeffs(int ix, InsertIterator b) const {
163 163
    int length = glp_get_mat_col(lp, ix, 0, 0);
164 164

	
165 165
    std::vector<int> indexes(length + 1);
166 166
    std::vector<Value> values(length + 1);
167 167

	
168 168
    glp_get_mat_col(lp, ix, &indexes.front(), &values.front());
169 169

	
170 170
    for (int i = 1; i  <= length; ++i) {
171 171
      *b = std::make_pair(indexes[i], values[i]);
172 172
      ++b;
173 173
    }
174 174
  }
175 175

	
176 176
  void GlpkBase::_setCoeff(int ix, int jx, Value value) {
177 177

	
178 178
    if (glp_get_num_cols(lp) < glp_get_num_rows(lp)) {
179 179

	
180 180
      int length = glp_get_mat_row(lp, ix, 0, 0);
181 181

	
182 182
      std::vector<int> indexes(length + 2);
183 183
      std::vector<Value> values(length + 2);
184 184

	
185 185
      glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
186 186

	
187 187
      //The following code does not suppose that the elements of the
188 188
      //array indexes are sorted
189 189
      bool found = false;
190 190
      for (int i = 1; i  <= length; ++i) {
191 191
        if (indexes[i] == jx) {
192 192
          found = true;
193 193
          values[i] = value;
194 194
          break;
195 195
        }
196 196
      }
197 197
      if (!found) {
198 198
        ++length;
199 199
        indexes[length] = jx;
200 200
        values[length] = value;
201 201
      }
202 202

	
203 203
      glp_set_mat_row(lp, ix, length, &indexes.front(), &values.front());
204 204

	
205 205
    } else {
206 206

	
207 207
      int length = glp_get_mat_col(lp, jx, 0, 0);
208 208

	
209 209
      std::vector<int> indexes(length + 2);
210 210
      std::vector<Value> values(length + 2);
211 211

	
212 212
      glp_get_mat_col(lp, jx, &indexes.front(), &values.front());
213 213

	
214 214
      //The following code does not suppose that the elements of the
215 215
      //array indexes are sorted
216 216
      bool found = false;
217 217
      for (int i = 1; i <= length; ++i) {
218 218
        if (indexes[i] == ix) {
219 219
          found = true;
220 220
          values[i] = value;
221 221
          break;
222 222
        }
223 223
      }
224 224
      if (!found) {
225 225
        ++length;
226 226
        indexes[length] = ix;
227 227
        values[length] = value;
228 228
      }
229 229

	
230 230
      glp_set_mat_col(lp, jx, length, &indexes.front(), &values.front());
231 231
    }
232 232

	
233 233
  }
234 234

	
235 235
  GlpkBase::Value GlpkBase::_getCoeff(int ix, int jx) const {
236 236

	
237 237
    int length = glp_get_mat_row(lp, ix, 0, 0);
238 238

	
239 239
    std::vector<int> indexes(length + 1);
240 240
    std::vector<Value> values(length + 1);
241 241

	
242 242
    glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
243 243

	
244 244
    for (int i = 1; i  <= length; ++i) {
245 245
      if (indexes[i] == jx) {
246 246
        return values[i];
247 247
      }
248 248
    }
249 249

	
250 250
    return 0;
251 251
  }
252 252

	
253 253
  void GlpkBase::_setColLowerBound(int i, Value lo) {
254 254
    LEMON_ASSERT(lo != INF, "Invalid bound");
255 255

	
256 256
    int b = glp_get_col_type(lp, i);
257 257
    double up = glp_get_col_ub(lp, i);
258 258
    if (lo == -INF) {
259 259
      switch (b) {
260 260
      case GLP_FR:
261 261
      case GLP_LO:
262 262
        glp_set_col_bnds(lp, i, GLP_FR, lo, up);
263 263
        break;
264 264
      case GLP_UP:
265 265
        break;
266 266
      case GLP_DB:
267 267
      case GLP_FX:
268 268
        glp_set_col_bnds(lp, i, GLP_UP, lo, up);
269 269
        break;
270 270
      default:
271 271
        break;
272 272
      }
273 273
    } else {
274 274
      switch (b) {
275 275
      case GLP_FR:
276 276
      case GLP_LO:
277 277
        glp_set_col_bnds(lp, i, GLP_LO, lo, up);
278 278
        break;
279 279
      case GLP_UP:
280 280
      case GLP_DB:
281 281
      case GLP_FX:
282 282
        if (lo == up)
283 283
          glp_set_col_bnds(lp, i, GLP_FX, lo, up);
284 284
        else
285 285
          glp_set_col_bnds(lp, i, GLP_DB, lo, up);
286 286
        break;
287 287
      default:
288 288
        break;
289 289
      }
290 290
    }
291 291
  }
292 292

	
293 293
  GlpkBase::Value GlpkBase::_getColLowerBound(int i) const {
294 294
    int b = glp_get_col_type(lp, i);
295 295
    switch (b) {
296 296
    case GLP_LO:
297 297
    case GLP_DB:
298 298
    case GLP_FX:
299 299
      return glp_get_col_lb(lp, i);
300 300
    default:
301 301
      return -INF;
302 302
    }
303 303
  }
304 304

	
305 305
  void GlpkBase::_setColUpperBound(int i, Value up) {
306 306
    LEMON_ASSERT(up != -INF, "Invalid bound");
307 307

	
308 308
    int b = glp_get_col_type(lp, i);
309 309
    double lo = glp_get_col_lb(lp, i);
310 310
    if (up == INF) {
311 311
      switch (b) {
312 312
      case GLP_FR:
313 313
      case GLP_LO:
314 314
        break;
315 315
      case GLP_UP:
316 316
        glp_set_col_bnds(lp, i, GLP_FR, lo, up);
317 317
        break;
318 318
      case GLP_DB:
319 319
      case GLP_FX:
320 320
        glp_set_col_bnds(lp, i, GLP_LO, lo, up);
321 321
        break;
322 322
      default:
323 323
        break;
324 324
      }
325 325
    } else {
326 326
      switch (b) {
327 327
      case GLP_FR:
328 328
        glp_set_col_bnds(lp, i, GLP_UP, lo, up);
329 329
        break;
330 330
      case GLP_UP:
331 331
        glp_set_col_bnds(lp, i, GLP_UP, lo, up);
332 332
        break;
333 333
      case GLP_LO:
334 334
      case GLP_DB:
335 335
      case GLP_FX:
336 336
        if (lo == up)
337 337
          glp_set_col_bnds(lp, i, GLP_FX, lo, up);
338 338
        else
339 339
          glp_set_col_bnds(lp, i, GLP_DB, lo, up);
340 340
        break;
341 341
      default:
342 342
        break;
343 343
      }
344 344
    }
345 345

	
346 346
  }
347 347

	
348 348
  GlpkBase::Value GlpkBase::_getColUpperBound(int i) const {
349 349
    int b = glp_get_col_type(lp, i);
350 350
      switch (b) {
351 351
      case GLP_UP:
352 352
      case GLP_DB:
353 353
      case GLP_FX:
354 354
        return glp_get_col_ub(lp, i);
355 355
      default:
356 356
        return INF;
357 357
      }
358 358
  }
359 359

	
360 360
  void GlpkBase::_setRowLowerBound(int i, Value lo) {
361 361
    LEMON_ASSERT(lo != INF, "Invalid bound");
362 362

	
363 363
    int b = glp_get_row_type(lp, i);
364 364
    double up = glp_get_row_ub(lp, i);
365 365
    if (lo == -INF) {
366 366
      switch (b) {
367 367
      case GLP_FR:
368 368
      case GLP_LO:
369 369
        glp_set_row_bnds(lp, i, GLP_FR, lo, up);
370 370
        break;
371 371
      case GLP_UP:
372 372
        break;
373 373
      case GLP_DB:
374 374
      case GLP_FX:
375 375
        glp_set_row_bnds(lp, i, GLP_UP, lo, up);
376 376
        break;
377 377
      default:
378 378
        break;
379 379
      }
380 380
    } else {
381 381
      switch (b) {
382 382
      case GLP_FR:
383 383
      case GLP_LO:
384 384
        glp_set_row_bnds(lp, i, GLP_LO, lo, up);
385 385
        break;
386 386
      case GLP_UP:
387 387
      case GLP_DB:
388 388
      case GLP_FX:
389 389
        if (lo == up)
390 390
          glp_set_row_bnds(lp, i, GLP_FX, lo, up);
391 391
        else
392 392
          glp_set_row_bnds(lp, i, GLP_DB, lo, up);
393 393
        break;
394 394
      default:
395 395
        break;
396 396
      }
397 397
    }
398 398

	
399 399
  }
400 400

	
401 401
  GlpkBase::Value GlpkBase::_getRowLowerBound(int i) const {
402 402
    int b = glp_get_row_type(lp, i);
403 403
    switch (b) {
404 404
    case GLP_LO:
405 405
    case GLP_DB:
406 406
    case GLP_FX:
407 407
      return glp_get_row_lb(lp, i);
408 408
    default:
409 409
      return -INF;
410 410
    }
411 411
  }
412 412

	
413 413
  void GlpkBase::_setRowUpperBound(int i, Value up) {
414 414
    LEMON_ASSERT(up != -INF, "Invalid bound");
415 415

	
416 416
    int b = glp_get_row_type(lp, i);
417 417
    double lo = glp_get_row_lb(lp, i);
418 418
    if (up == INF) {
419 419
      switch (b) {
420 420
      case GLP_FR:
421 421
      case GLP_LO:
422 422
        break;
423 423
      case GLP_UP:
424 424
        glp_set_row_bnds(lp, i, GLP_FR, lo, up);
425 425
        break;
426 426
      case GLP_DB:
427 427
      case GLP_FX:
428 428
        glp_set_row_bnds(lp, i, GLP_LO, lo, up);
429 429
        break;
430 430
      default:
431 431
        break;
432 432
      }
433 433
    } else {
434 434
      switch (b) {
435 435
      case GLP_FR:
436 436
        glp_set_row_bnds(lp, i, GLP_UP, lo, up);
437 437
        break;
438 438
      case GLP_UP:
439 439
        glp_set_row_bnds(lp, i, GLP_UP, lo, up);
440 440
        break;
441 441
      case GLP_LO:
442 442
      case GLP_DB:
443 443
      case GLP_FX:
444 444
        if (lo == up)
445 445
          glp_set_row_bnds(lp, i, GLP_FX, lo, up);
446 446
        else
447 447
          glp_set_row_bnds(lp, i, GLP_DB, lo, up);
448 448
        break;
449 449
      default:
450 450
        break;
451 451
      }
452 452
    }
453 453
  }
454 454

	
455 455
  GlpkBase::Value GlpkBase::_getRowUpperBound(int i) const {
456 456
    int b = glp_get_row_type(lp, i);
457 457
    switch (b) {
458 458
    case GLP_UP:
459 459
    case GLP_DB:
460 460
    case GLP_FX:
461 461
      return glp_get_row_ub(lp, i);
462 462
    default:
463 463
      return INF;
464 464
    }
465 465
  }
466 466

	
467 467
  void GlpkBase::_setObjCoeffs(ExprIterator b, ExprIterator e) {
468 468
    for (int i = 1; i <= glp_get_num_cols(lp); ++i) {
469 469
      glp_set_obj_coef(lp, i, 0.0);
470 470
    }
471 471
    for (ExprIterator it = b; it != e; ++it) {
472 472
      glp_set_obj_coef(lp, it->first, it->second);
473 473
    }
474 474
  }
475 475

	
476 476
  void GlpkBase::_getObjCoeffs(InsertIterator b) const {
477 477
    for (int i = 1; i <= glp_get_num_cols(lp); ++i) {
478 478
      Value val = glp_get_obj_coef(lp, i);
479 479
      if (val != 0.0) {
480 480
        *b = std::make_pair(i, val);
481 481
        ++b;
482 482
      }
483 483
    }
484 484
  }
485 485

	
486 486
  void GlpkBase::_setObjCoeff(int i, Value obj_coef) {
487 487
    //i = 0 means the constant term (shift)
488 488
    glp_set_obj_coef(lp, i, obj_coef);
489 489
  }
490 490

	
491 491
  GlpkBase::Value GlpkBase::_getObjCoeff(int i) const {
492 492
    //i = 0 means the constant term (shift)
493 493
    return glp_get_obj_coef(lp, i);
494 494
  }
495 495

	
496 496
  void GlpkBase::_setSense(GlpkBase::Sense sense) {
497 497
    switch (sense) {
498 498
    case MIN:
499 499
      glp_set_obj_dir(lp, GLP_MIN);
500 500
      break;
501 501
    case MAX:
502 502
      glp_set_obj_dir(lp, GLP_MAX);
503 503
      break;
504 504
    }
505 505
  }
506 506

	
507 507
  GlpkBase::Sense GlpkBase::_getSense() const {
508 508
    switch(glp_get_obj_dir(lp)) {
509 509
    case GLP_MIN:
510 510
      return MIN;
511 511
    case GLP_MAX:
512 512
      return MAX;
513 513
    default:
514 514
      LEMON_ASSERT(false, "Wrong sense");
515 515
      return GlpkBase::Sense();
516 516
    }
517 517
  }
518 518

	
519 519
  void GlpkBase::_clear() {
520 520
    glp_erase_prob(lp);
521 521
    rows.clear();
522 522
    cols.clear();
523 523
  }
524 524

	
525
  void GlpkBase::freeEnv() {
526
    glp_free_env();
527
  }
528

	
525 529
  // GlpkLp members
526 530

	
527 531
  GlpkLp::GlpkLp()
528 532
    : LpBase(), GlpkBase(), LpSolver() {
529 533
    messageLevel(MESSAGE_NO_OUTPUT);
530 534
  }
531 535

	
532 536
  GlpkLp::GlpkLp(const GlpkLp& other)
533 537
    : LpBase(other), GlpkBase(other), LpSolver(other) {
534 538
    messageLevel(MESSAGE_NO_OUTPUT);
535 539
  }
536 540

	
537 541
  GlpkLp* GlpkLp::_newSolver() const { return new GlpkLp; }
538 542
  GlpkLp* GlpkLp::_cloneSolver() const { return new GlpkLp(*this); }
539 543

	
540 544
  const char* GlpkLp::_solverName() const { return "GlpkLp"; }
541 545

	
542 546
  void GlpkLp::_clear_temporals() {
543 547
    _primal_ray.clear();
544 548
    _dual_ray.clear();
545 549
  }
546 550

	
547 551
  GlpkLp::SolveExitStatus GlpkLp::_solve() {
548 552
    return solvePrimal();
549 553
  }
550 554

	
551 555
  GlpkLp::SolveExitStatus GlpkLp::solvePrimal() {
552 556
    _clear_temporals();
553 557

	
554 558
    glp_smcp smcp;
555 559
    glp_init_smcp(&smcp);
556 560

	
557 561
    switch (_message_level) {
558 562
    case MESSAGE_NO_OUTPUT:
559 563
      smcp.msg_lev = GLP_MSG_OFF;
560 564
      break;
561 565
    case MESSAGE_ERROR_MESSAGE:
562 566
      smcp.msg_lev = GLP_MSG_ERR;
563 567
      break;
564 568
    case MESSAGE_NORMAL_OUTPUT:
565 569
      smcp.msg_lev = GLP_MSG_ON;
566 570
      break;
567 571
    case MESSAGE_FULL_OUTPUT:
568 572
      smcp.msg_lev = GLP_MSG_ALL;
569 573
      break;
570 574
    }
571 575

	
572 576
    if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
573 577
    return SOLVED;
574 578
  }
575 579

	
576 580
  GlpkLp::SolveExitStatus GlpkLp::solveDual() {
577 581
    _clear_temporals();
578 582

	
579 583
    glp_smcp smcp;
580 584
    glp_init_smcp(&smcp);
581 585

	
582 586
    switch (_message_level) {
583 587
    case MESSAGE_NO_OUTPUT:
584 588
      smcp.msg_lev = GLP_MSG_OFF;
585 589
      break;
586 590
    case MESSAGE_ERROR_MESSAGE:
587 591
      smcp.msg_lev = GLP_MSG_ERR;
588 592
      break;
589 593
    case MESSAGE_NORMAL_OUTPUT:
590 594
      smcp.msg_lev = GLP_MSG_ON;
591 595
      break;
592 596
    case MESSAGE_FULL_OUTPUT:
593 597
      smcp.msg_lev = GLP_MSG_ALL;
594 598
      break;
595 599
    }
596 600
    smcp.meth = GLP_DUAL;
597 601

	
598 602
    if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
599 603
    return SOLVED;
600 604
  }
601 605

	
602 606
  GlpkLp::Value GlpkLp::_getPrimal(int i) const {
603 607
    return glp_get_col_prim(lp, i);
604 608
  }
605 609

	
606 610
  GlpkLp::Value GlpkLp::_getDual(int i) const {
607 611
    return glp_get_row_dual(lp, i);
608 612
  }
609 613

	
610 614
  GlpkLp::Value GlpkLp::_getPrimalValue() const {
611 615
    return glp_get_obj_val(lp);
612 616
  }
613 617

	
614 618
  GlpkLp::VarStatus GlpkLp::_getColStatus(int i) const {
615 619
    switch (glp_get_col_stat(lp, i)) {
616 620
    case GLP_BS:
617 621
      return BASIC;
618 622
    case GLP_UP:
619 623
      return UPPER;
620 624
    case GLP_LO:
621 625
      return LOWER;
622 626
    case GLP_NF:
623 627
      return FREE;
624 628
    case GLP_NS:
625 629
      return FIXED;
626 630
    default:
627 631
      LEMON_ASSERT(false, "Wrong column status");
628 632
      return GlpkLp::VarStatus();
629 633
    }
630 634
  }
631 635

	
632 636
  GlpkLp::VarStatus GlpkLp::_getRowStatus(int i) const {
633 637
    switch (glp_get_row_stat(lp, i)) {
634 638
    case GLP_BS:
635 639
      return BASIC;
636 640
    case GLP_UP:
637 641
      return UPPER;
638 642
    case GLP_LO:
639 643
      return LOWER;
640 644
    case GLP_NF:
641 645
      return FREE;
642 646
    case GLP_NS:
643 647
      return FIXED;
644 648
    default:
645 649
      LEMON_ASSERT(false, "Wrong row status");
646 650
      return GlpkLp::VarStatus();
647 651
    }
648 652
  }
649 653

	
650 654
  GlpkLp::Value GlpkLp::_getPrimalRay(int i) const {
651 655
    if (_primal_ray.empty()) {
652 656
      int row_num = glp_get_num_rows(lp);
653 657
      int col_num = glp_get_num_cols(lp);
654 658

	
655 659
      _primal_ray.resize(col_num + 1, 0.0);
656 660

	
657 661
      int index = glp_get_unbnd_ray(lp);
658 662
      if (index != 0) {
659 663
        // The primal ray is found in primal simplex second phase
660 664
        LEMON_ASSERT((index <= row_num ? glp_get_row_stat(lp, index) :
661 665
                      glp_get_col_stat(lp, index - row_num)) != GLP_BS,
662 666
                     "Wrong primal ray");
663 667

	
664 668
        bool negate = glp_get_obj_dir(lp) == GLP_MAX;
665 669

	
666 670
        if (index > row_num) {
667 671
          _primal_ray[index - row_num] = 1.0;
668 672
          if (glp_get_col_dual(lp, index - row_num) > 0) {
669 673
            negate = !negate;
670 674
          }
671 675
        } else {
672 676
          if (glp_get_row_dual(lp, index) > 0) {
673 677
            negate = !negate;
674 678
          }
675 679
        }
676 680

	
677 681
        std::vector<int> ray_indexes(row_num + 1);
678 682
        std::vector<Value> ray_values(row_num + 1);
679 683
        int ray_length = glp_eval_tab_col(lp, index, &ray_indexes.front(),
680 684
                                          &ray_values.front());
681 685

	
682 686
        for (int i = 1; i <= ray_length; ++i) {
683 687
          if (ray_indexes[i] > row_num) {
684 688
            _primal_ray[ray_indexes[i] - row_num] = ray_values[i];
685 689
          }
686 690
        }
687 691

	
688 692
        if (negate) {
689 693
          for (int i = 1; i <= col_num; ++i) {
690 694
            _primal_ray[i] = - _primal_ray[i];
691 695
          }
692 696
        }
693 697
      } else {
694 698
        for (int i = 1; i <= col_num; ++i) {
695 699
          _primal_ray[i] = glp_get_col_prim(lp, i);
696 700
        }
697 701
      }
698 702
    }
699 703
    return _primal_ray[i];
700 704
  }
701 705

	
702 706
  GlpkLp::Value GlpkLp::_getDualRay(int i) const {
703 707
    if (_dual_ray.empty()) {
704 708
      int row_num = glp_get_num_rows(lp);
705 709

	
706 710
      _dual_ray.resize(row_num + 1, 0.0);
707 711

	
708 712
      int index = glp_get_unbnd_ray(lp);
709 713
      if (index != 0) {
710 714
        // The dual ray is found in dual simplex second phase
711 715
        LEMON_ASSERT((index <= row_num ? glp_get_row_stat(lp, index) :
712 716
                      glp_get_col_stat(lp, index - row_num)) == GLP_BS,
713 717

	
714 718
                     "Wrong dual ray");
715 719

	
716 720
        int idx;
717 721
        bool negate = false;
718 722

	
719 723
        if (index > row_num) {
720 724
          idx = glp_get_col_bind(lp, index - row_num);
721 725
          if (glp_get_col_prim(lp, index - row_num) >
722 726
              glp_get_col_ub(lp, index - row_num)) {
723 727
            negate = true;
724 728
          }
725 729
        } else {
726 730
          idx = glp_get_row_bind(lp, index);
727 731
          if (glp_get_row_prim(lp, index) > glp_get_row_ub(lp, index)) {
728 732
            negate = true;
729 733
          }
730 734
        }
731 735

	
732 736
        _dual_ray[idx] = negate ?  - 1.0 : 1.0;
733 737

	
734 738
        glp_btran(lp, &_dual_ray.front());
735 739
      } else {
736 740
        double eps = 1e-7;
737 741
        // The dual ray is found in primal simplex first phase
738 742
        // We assume that the glpk minimizes the slack to get feasible solution
739 743
        for (int i = 1; i <= row_num; ++i) {
740 744
          int index = glp_get_bhead(lp, i);
741 745
          if (index <= row_num) {
742 746
            double res = glp_get_row_prim(lp, index);
743 747
            if (res > glp_get_row_ub(lp, index) + eps) {
744 748
              _dual_ray[i] = -1;
745 749
            } else if (res < glp_get_row_lb(lp, index) - eps) {
746 750
              _dual_ray[i] = 1;
747 751
            } else {
748 752
              _dual_ray[i] = 0;
749 753
            }
750 754
            _dual_ray[i] *= glp_get_rii(lp, index);
751 755
          } else {
752 756
            double res = glp_get_col_prim(lp, index - row_num);
753 757
            if (res > glp_get_col_ub(lp, index - row_num) + eps) {
754 758
              _dual_ray[i] = -1;
755 759
            } else if (res < glp_get_col_lb(lp, index - row_num) - eps) {
756 760
              _dual_ray[i] = 1;
757 761
            } else {
758 762
              _dual_ray[i] = 0;
759 763
            }
760 764
            _dual_ray[i] /= glp_get_sjj(lp, index - row_num);
761 765
          }
762 766
        }
763 767

	
764 768
        glp_btran(lp, &_dual_ray.front());
765 769

	
766 770
        for (int i = 1; i <= row_num; ++i) {
767 771
          _dual_ray[i] /= glp_get_rii(lp, i);
768 772
        }
769 773
      }
770 774
    }
771 775
    return _dual_ray[i];
772 776
  }
773 777

	
774 778
  GlpkLp::ProblemType GlpkLp::_getPrimalType() const {
775 779
    if (glp_get_status(lp) == GLP_OPT)
776 780
      return OPTIMAL;
777 781
    switch (glp_get_prim_stat(lp)) {
778 782
    case GLP_UNDEF:
779 783
      return UNDEFINED;
780 784
    case GLP_FEAS:
781 785
    case GLP_INFEAS:
782 786
      if (glp_get_dual_stat(lp) == GLP_NOFEAS) {
783 787
        return UNBOUNDED;
784 788
      } else {
785 789
        return UNDEFINED;
786 790
      }
787 791
    case GLP_NOFEAS:
788 792
      return INFEASIBLE;
789 793
    default:
790 794
      LEMON_ASSERT(false, "Wrong primal type");
791 795
      return  GlpkLp::ProblemType();
792 796
    }
793 797
  }
794 798

	
795 799
  GlpkLp::ProblemType GlpkLp::_getDualType() const {
796 800
    if (glp_get_status(lp) == GLP_OPT)
797 801
      return OPTIMAL;
798 802
    switch (glp_get_dual_stat(lp)) {
799 803
    case GLP_UNDEF:
800 804
      return UNDEFINED;
801 805
    case GLP_FEAS:
802 806
    case GLP_INFEAS:
803 807
      if (glp_get_prim_stat(lp) == GLP_NOFEAS) {
804 808
        return UNBOUNDED;
805 809
      } else {
806 810
        return UNDEFINED;
807 811
      }
808 812
    case GLP_NOFEAS:
809 813
      return INFEASIBLE;
810 814
    default:
811 815
      LEMON_ASSERT(false, "Wrong primal type");
812 816
      return  GlpkLp::ProblemType();
813 817
    }
814 818
  }
815 819

	
816 820
  void GlpkLp::presolver(bool b) {
817 821
    lpx_set_int_parm(lp, LPX_K_PRESOL, b ? 1 : 0);
818 822
  }
819 823

	
820 824
  void GlpkLp::messageLevel(MessageLevel m) {
821 825
    _message_level = m;
822 826
  }
823 827

	
824 828
  // GlpkMip members
825 829

	
826 830
  GlpkMip::GlpkMip()
827 831
    : LpBase(), GlpkBase(), MipSolver() {
828 832
    messageLevel(MESSAGE_NO_OUTPUT);
829 833
  }
830 834

	
831 835
  GlpkMip::GlpkMip(const GlpkMip& other)
832 836
    : LpBase(), GlpkBase(other), MipSolver() {
833 837
    messageLevel(MESSAGE_NO_OUTPUT);
834 838
  }
835 839

	
836 840
  void GlpkMip::_setColType(int i, GlpkMip::ColTypes col_type) {
837 841
    switch (col_type) {
838 842
    case INTEGER:
839 843
      glp_set_col_kind(lp, i, GLP_IV);
840 844
      break;
841 845
    case REAL:
842 846
      glp_set_col_kind(lp, i, GLP_CV);
843 847
      break;
844 848
    }
845 849
  }
846 850

	
847 851
  GlpkMip::ColTypes GlpkMip::_getColType(int i) const {
848 852
    switch (glp_get_col_kind(lp, i)) {
849 853
    case GLP_IV:
850 854
    case GLP_BV:
851 855
      return INTEGER;
852 856
    default:
853 857
      return REAL;
854 858
    }
855 859

	
856 860
  }
857 861

	
858 862
  GlpkMip::SolveExitStatus GlpkMip::_solve() {
859 863
    glp_smcp smcp;
860 864
    glp_init_smcp(&smcp);
861 865

	
862 866
    switch (_message_level) {
863 867
    case MESSAGE_NO_OUTPUT:
864 868
      smcp.msg_lev = GLP_MSG_OFF;
865 869
      break;
866 870
    case MESSAGE_ERROR_MESSAGE:
867 871
      smcp.msg_lev = GLP_MSG_ERR;
868 872
      break;
869 873
    case MESSAGE_NORMAL_OUTPUT:
870 874
      smcp.msg_lev = GLP_MSG_ON;
871 875
      break;
872 876
    case MESSAGE_FULL_OUTPUT:
873 877
      smcp.msg_lev = GLP_MSG_ALL;
874 878
      break;
875 879
    }
876 880
    smcp.meth = GLP_DUAL;
877 881

	
878 882
    if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
879 883
    if (glp_get_status(lp) != GLP_OPT) return SOLVED;
880 884

	
881 885
    glp_iocp iocp;
882 886
    glp_init_iocp(&iocp);
883 887

	
884 888
    switch (_message_level) {
885 889
    case MESSAGE_NO_OUTPUT:
886 890
      iocp.msg_lev = GLP_MSG_OFF;
887 891
      break;
888 892
    case MESSAGE_ERROR_MESSAGE:
889 893
      iocp.msg_lev = GLP_MSG_ERR;
890 894
      break;
891 895
    case MESSAGE_NORMAL_OUTPUT:
892 896
      iocp.msg_lev = GLP_MSG_ON;
893 897
      break;
894 898
    case MESSAGE_FULL_OUTPUT:
895 899
      iocp.msg_lev = GLP_MSG_ALL;
896 900
      break;
897 901
    }
898 902

	
899 903
    if (glp_intopt(lp, &iocp) != 0) return UNSOLVED;
900 904
    return SOLVED;
901 905
  }
902 906

	
903 907

	
904 908
  GlpkMip::ProblemType GlpkMip::_getType() const {
905 909
    switch (glp_get_status(lp)) {
906 910
    case GLP_OPT:
907 911
      switch (glp_mip_status(lp)) {
908 912
      case GLP_UNDEF:
Ignore white space 6 line context
1 1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
2 2
 *
3 3
 * This file is a part of LEMON, a generic C++ optimization library.
4 4
 *
5 5
 * Copyright (C) 2003-2008
6 6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 8
 *
9 9
 * Permission to use, modify and distribute this software is granted
10 10
 * provided that this copyright notice appears in all copies. For
11 11
 * precise terms see the accompanying LICENSE file.
12 12
 *
13 13
 * This software is provided "AS IS" with no warranty of any kind,
14 14
 * express or implied, and with no claim as to its suitability for any
15 15
 * purpose.
16 16
 *
17 17
 */
18 18

	
19 19
#ifndef LEMON_GLPK_H
20 20
#define LEMON_GLPK_H
21 21

	
22 22
///\file
23 23
///\brief Header of the LEMON-GLPK lp solver interface.
24 24
///\ingroup lp_group
25 25

	
26 26
#include <lemon/lp_base.h>
27 27

	
28 28
// forward declaration
29 29
#ifndef _GLP_PROB
30 30
#define _GLP_PROB
31 31
typedef struct { double _prob; } glp_prob;
32 32
/* LP/MIP problem object */
33 33
#endif
34 34

	
35 35
namespace lemon {
36 36

	
37 37

	
38 38
  /// \brief Base interface for the GLPK LP and MIP solver
39 39
  ///
40 40
  /// This class implements the common interface of the GLPK LP and MIP solver.
41 41
  /// \ingroup lp_group
42 42
  class GlpkBase : virtual public LpBase {
43 43
  protected:
44 44

	
45 45
    typedef glp_prob LPX;
46 46
    glp_prob* lp;
47 47

	
48 48
    GlpkBase();
49 49
    GlpkBase(const GlpkBase&);
50 50
    virtual ~GlpkBase();
51 51

	
52 52
  protected:
53 53

	
54 54
    virtual int _addCol();
55 55
    virtual int _addRow();
56 56

	
57 57
    virtual void _eraseCol(int i);
58 58
    virtual void _eraseRow(int i);
59 59

	
60 60
    virtual void _eraseColId(int i);
61 61
    virtual void _eraseRowId(int i);
62 62

	
63 63
    virtual void _getColName(int col, std::string& name) const;
64 64
    virtual void _setColName(int col, const std::string& name);
65 65
    virtual int _colByName(const std::string& name) const;
66 66

	
67 67
    virtual void _getRowName(int row, std::string& name) const;
68 68
    virtual void _setRowName(int row, const std::string& name);
69 69
    virtual int _rowByName(const std::string& name) const;
70 70

	
71 71
    virtual void _setRowCoeffs(int i, ExprIterator b, ExprIterator e);
72 72
    virtual void _getRowCoeffs(int i, InsertIterator b) const;
73 73

	
74 74
    virtual void _setColCoeffs(int i, ExprIterator b, ExprIterator e);
75 75
    virtual void _getColCoeffs(int i, InsertIterator b) const;
76 76

	
77 77
    virtual void _setCoeff(int row, int col, Value value);
78 78
    virtual Value _getCoeff(int row, int col) const;
79 79

	
80 80
    virtual void _setColLowerBound(int i, Value value);
81 81
    virtual Value _getColLowerBound(int i) const;
82 82

	
83 83
    virtual void _setColUpperBound(int i, Value value);
84 84
    virtual Value _getColUpperBound(int i) const;
85 85

	
86 86
    virtual void _setRowLowerBound(int i, Value value);
87 87
    virtual Value _getRowLowerBound(int i) const;
88 88

	
89 89
    virtual void _setRowUpperBound(int i, Value value);
90 90
    virtual Value _getRowUpperBound(int i) const;
91 91

	
92 92
    virtual void _setObjCoeffs(ExprIterator b, ExprIterator e);
93 93
    virtual void _getObjCoeffs(InsertIterator b) const;
94 94

	
95 95
    virtual void _setObjCoeff(int i, Value obj_coef);
96 96
    virtual Value _getObjCoeff(int i) const;
97 97

	
98 98
    virtual void _setSense(Sense);
99 99
    virtual Sense _getSense() const;
100 100

	
101 101
    virtual void _clear();
102 102

	
103 103
  public:
104 104

	
105
    /// \brief Deallocates the globally allocated memory of GLPK.
106

	
107
    /// Deallocates the globally allocated memory of GLPK.  \note
108
    /// Usually, it do not have to be called, because the GLPK use
109
    /// only a small amount of global memory, and it is deallocated
110
    /// automatically at the end of program.
111
    static void freeEnv();
112

	
105 113
    ///Pointer to the underlying GLPK data structure.
106 114
    LPX *lpx() {return lp;}
107 115
    ///Const pointer to the underlying GLPK data structure.
108 116
    const LPX *lpx() const {return lp;}
109 117

	
110 118
    ///Returns the constraint identifier understood by GLPK.
111 119
    int lpxRow(Row r) const { return rows(id(r)); }
112 120

	
113 121
    ///Returns the variable identifier understood by GLPK.
114 122
    int lpxCol(Col c) const { return cols(id(c)); }
115 123

	
116 124
  };
117 125

	
118 126
  /// \brief Interface for the GLPK LP solver
119 127
  ///
120 128
  /// This class implements an interface for the GLPK LP solver.
121 129
  ///\ingroup lp_group
122 130
  class GlpkLp : public GlpkBase, public LpSolver {
123 131
  public:
124 132

	
125 133
    ///\e
126 134
    GlpkLp();
127 135
    ///\e
128 136
    GlpkLp(const GlpkLp&);
129 137

	
130 138
  private:
131 139

	
132 140
    mutable std::vector<double> _primal_ray;
133 141
    mutable std::vector<double> _dual_ray;
134 142

	
135 143
    void _clear_temporals();
136 144

	
137 145
  protected:
138 146

	
139 147
    virtual GlpkLp* _cloneSolver() const;
140 148
    virtual GlpkLp* _newSolver() const;
141 149

	
142 150
    virtual const char* _solverName() const;
143 151

	
144 152
    virtual SolveExitStatus _solve();
145 153
    virtual Value _getPrimal(int i) const;
146 154
    virtual Value _getDual(int i) const;
147 155

	
148 156
    virtual Value _getPrimalValue() const;
149 157

	
150 158
    virtual VarStatus _getColStatus(int i) const;
151 159
    virtual VarStatus _getRowStatus(int i) const;
152 160

	
153 161
    virtual Value _getPrimalRay(int i) const;
154 162
    virtual Value _getDualRay(int i) const;
155 163

	
156 164
    ///\todo It should be clarified
157 165
    ///
158 166
    virtual ProblemType _getPrimalType() const;
159 167
    virtual ProblemType _getDualType() const;
160 168

	
161 169
  public:
162 170

	
163 171
    ///Solve with primal simplex
164 172
    SolveExitStatus solvePrimal();
165 173

	
166 174
    ///Solve with dual simplex
167 175
    SolveExitStatus solveDual();
168 176

	
169 177
    ///Turns on or off the presolver
170 178

	
171 179
    ///Turns on (\c b is \c true) or off (\c b is \c false) the presolver
172 180
    ///
173 181
    ///The presolver is off by default.
174 182
    void presolver(bool b);
175 183

	
176 184
    ///Enum for \c messageLevel() parameter
177 185
    enum MessageLevel {
178 186
      /// no output (default value)
179 187
      MESSAGE_NO_OUTPUT = 0,
180 188
      /// error messages only
181 189
      MESSAGE_ERROR_MESSAGE = 1,
182 190
      /// normal output
183 191
      MESSAGE_NORMAL_OUTPUT = 2,
184 192
      /// full output (includes informational messages)
185 193
      MESSAGE_FULL_OUTPUT = 3
186 194
    };
187 195

	
188 196
  private:
189 197

	
190 198
    MessageLevel _message_level;
191 199

	
192 200
  public:
193 201

	
194 202
    ///Set the verbosity of the messages
195 203

	
196 204
    ///Set the verbosity of the messages
197 205
    ///
198 206
    ///\param m is the level of the messages output by the solver routines.
199 207
    void messageLevel(MessageLevel m);
200 208
  };
201 209

	
202 210
  /// \brief Interface for the GLPK MIP solver
203 211
  ///
204 212
  /// This class implements an interface for the GLPK MIP solver.
205 213
  ///\ingroup lp_group
206 214
  class GlpkMip : public GlpkBase, public MipSolver {
207 215
  public:
208 216

	
209 217
    ///\e
210 218
    GlpkMip();
211 219
    ///\e
212 220
    GlpkMip(const GlpkMip&);
213 221

	
214 222
  protected:
215 223

	
216 224
    virtual GlpkMip* _cloneSolver() const;
217 225
    virtual GlpkMip* _newSolver() const;
218 226

	
219 227
    virtual const char* _solverName() const;
220 228

	
221 229
    virtual ColTypes _getColType(int col) const;
222 230
    virtual void _setColType(int col, ColTypes col_type);
223 231

	
224 232
    virtual SolveExitStatus _solve();
225 233
    virtual ProblemType _getType() const;
226 234
    virtual Value _getSol(int i) const;
227 235
    virtual Value _getSolValue() const;
228 236

	
229 237
    ///Enum for \c messageLevel() parameter
230 238
    enum MessageLevel {
231 239
      /// no output (default value)
232 240
      MESSAGE_NO_OUTPUT = 0,
233 241
      /// error messages only
234 242
      MESSAGE_ERROR_MESSAGE = 1,
235 243
      /// normal output
236 244
      MESSAGE_NORMAL_OUTPUT = 2,
237 245
      /// full output (includes informational messages)
238 246
      MESSAGE_FULL_OUTPUT = 3
239 247
    };
240 248

	
241 249
  private:
242 250

	
243 251
    MessageLevel _message_level;
244 252

	
245 253
  public:
246 254

	
247 255
    ///Set the verbosity of the messages
248 256

	
249 257
    ///Set the verbosity of the messages
250 258
    ///
251 259
    ///\param m is the level of the messages output by the solver routines.
252 260
    void messageLevel(MessageLevel m);
253 261
  };
254 262

	
255 263

	
256 264
} //END OF NAMESPACE LEMON
257 265

	
258 266
#endif //LEMON_GLPK_H
259 267

	
Ignore white space 6 line context
1 1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
2 2
 *
3 3
 * This file is a part of LEMON, a generic C++ optimization library.
4 4
 *
5 5
 * Copyright (C) 2003-2008
6 6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 8
 *
9 9
 * Permission to use, modify and distribute this software is granted
10 10
 * provided that this copyright notice appears in all copies. For
11 11
 * precise terms see the accompanying LICENSE file.
12 12
 *
13 13
 * This software is provided "AS IS" with no warranty of any kind,
14 14
 * express or implied, and with no claim as to its suitability for any
15 15
 * purpose.
16 16
 *
17 17
 */
18 18

	
19 19
#include <sstream>
20 20
#include <lemon/lp_skeleton.h>
21 21
#include "test_tools.h"
22 22
#include <lemon/tolerance.h>
23 23

	
24 24
#ifdef HAVE_CONFIG_H
25 25
#include <lemon/config.h>
26 26
#endif
27 27

	
28 28
#ifdef HAVE_GLPK
29 29
#include <lemon/glpk.h>
30 30
#endif
31 31

	
32 32
#ifdef HAVE_CPLEX
33 33
#include <lemon/cplex.h>
34 34
#endif
35 35

	
36 36
#ifdef HAVE_SOPLEX
37 37
#include <lemon/soplex.h>
38 38
#endif
39 39

	
40 40
#ifdef HAVE_CLP
41 41
#include <lemon/clp.h>
42 42
#endif
43 43

	
44 44
using namespace lemon;
45 45

	
46 46
void lpTest(LpSolver& lp)
47 47
{
48 48

	
49 49
  typedef LpSolver LP;
50 50

	
51 51
  std::vector<LP::Col> x(10);
52 52
  //  for(int i=0;i<10;i++) x.push_back(lp.addCol());
53 53
  lp.addColSet(x);
54 54
  lp.colLowerBound(x,1);
55 55
  lp.colUpperBound(x,1);
56 56
  lp.colBounds(x,1,2);
57 57

	
58 58
  std::vector<LP::Col> y(10);
59 59
  lp.addColSet(y);
60 60

	
61 61
  lp.colLowerBound(y,1);
62 62
  lp.colUpperBound(y,1);
63 63
  lp.colBounds(y,1,2);
64 64

	
65 65
  std::map<int,LP::Col> z;
66 66

	
67 67
  z.insert(std::make_pair(12,INVALID));
68 68
  z.insert(std::make_pair(2,INVALID));
69 69
  z.insert(std::make_pair(7,INVALID));
70 70
  z.insert(std::make_pair(5,INVALID));
71 71

	
72 72
  lp.addColSet(z);
73 73

	
74 74
  lp.colLowerBound(z,1);
75 75
  lp.colUpperBound(z,1);
76 76
  lp.colBounds(z,1,2);
77 77

	
78 78
  {
79 79
    LP::Expr e,f,g;
80 80
    LP::Col p1,p2,p3,p4,p5;
81 81
    LP::Constr c;
82 82

	
83 83
    p1=lp.addCol();
84 84
    p2=lp.addCol();
85 85
    p3=lp.addCol();
86 86
    p4=lp.addCol();
87 87
    p5=lp.addCol();
88 88

	
89 89
    e[p1]=2;
90 90
    *e=12;
91 91
    e[p1]+=2;
92 92
    *e+=12;
93 93
    e[p1]-=2;
94 94
    *e-=12;
95 95

	
96 96
    e=2;
97 97
    e=2.2;
98 98
    e=p1;
99 99
    e=f;
100 100

	
101 101
    e+=2;
102 102
    e+=2.2;
103 103
    e+=p1;
104 104
    e+=f;
105 105

	
106 106
    e-=2;
107 107
    e-=2.2;
108 108
    e-=p1;
109 109
    e-=f;
110 110

	
111 111
    e*=2;
112 112
    e*=2.2;
113 113
    e/=2;
114 114
    e/=2.2;
115 115

	
116 116
    e=((p1+p2)+(p1-p2)+(p1+12)+(12+p1)+(p1-12)+(12-p1)+
117 117
       (f+12)+(12+f)+(p1+f)+(f+p1)+(f+g)+
118 118
       (f-12)+(12-f)+(p1-f)+(f-p1)+(f-g)+
119 119
       2.2*f+f*2.2+f/2.2+
120 120
       2*f+f*2+f/2+
121 121
       2.2*p1+p1*2.2+p1/2.2+
122 122
       2*p1+p1*2+p1/2
123 123
       );
124 124

	
125 125

	
126 126
    c = (e  <= f  );
127 127
    c = (e  <= 2.2);
128 128
    c = (e  <= 2  );
129 129
    c = (e  <= p1 );
130 130
    c = (2.2<= f  );
131 131
    c = (2  <= f  );
132 132
    c = (p1 <= f  );
133 133
    c = (p1 <= p2 );
134 134
    c = (p1 <= 2.2);
135 135
    c = (p1 <= 2  );
136 136
    c = (2.2<= p2 );
137 137
    c = (2  <= p2 );
138 138

	
139 139
    c = (e  >= f  );
140 140
    c = (e  >= 2.2);
141 141
    c = (e  >= 2  );
142 142
    c = (e  >= p1 );
143 143
    c = (2.2>= f  );
144 144
    c = (2  >= f  );
145 145
    c = (p1 >= f  );
146 146
    c = (p1 >= p2 );
147 147
    c = (p1 >= 2.2);
148 148
    c = (p1 >= 2  );
149 149
    c = (2.2>= p2 );
150 150
    c = (2  >= p2 );
151 151

	
152 152
    c = (e  == f  );
153 153
    c = (e  == 2.2);
154 154
    c = (e  == 2  );
155 155
    c = (e  == p1 );
156 156
    c = (2.2== f  );
157 157
    c = (2  == f  );
158 158
    c = (p1 == f  );
159 159
    //c = (p1 == p2 );
160 160
    c = (p1 == 2.2);
161 161
    c = (p1 == 2  );
162 162
    c = (2.2== p2 );
163 163
    c = (2  == p2 );
164 164

	
165 165
    c = ((2 <= e) <= 3);
166 166
    c = ((2 <= p1) <= 3);
167 167

	
168 168
    c = ((2 >= e) >= 3);
169 169
    c = ((2 >= p1) >= 3);
170 170

	
171 171
    e[x[3]]=2;
172 172
    e[x[3]]=4;
173 173
    e[x[3]]=1;
174 174
    *e=12;
175 175

	
176 176
    lp.addRow(-LP::INF,e,23);
177 177
    lp.addRow(-LP::INF,3.0*(x[1]+x[2]/2)-x[3],23);
178 178
    lp.addRow(-LP::INF,3.0*(x[1]+x[2]*2-5*x[3]+12-x[4]/3)+2*x[4]-4,23);
179 179

	
180 180
    lp.addRow(x[1]+x[3]<=x[5]-3);
181 181
    lp.addRow((-7<=x[1]+x[3]-12)<=3);
182 182
    lp.addRow(x[1]<=x[5]);
183 183

	
184 184
    std::ostringstream buf;
185 185

	
186 186

	
187 187
    e=((p1+p2)+(p1-0.99*p2));
188 188
    //e.prettyPrint(std::cout);
189 189
    //(e<=2).prettyPrint(std::cout);
190 190
    double tolerance=0.001;
191 191
    e.simplify(tolerance);
192 192
    buf << "Coeff. of p2 should be 0.01";
193 193
    check(e[p2]>0, buf.str());
194 194

	
195 195
    tolerance=0.02;
196 196
    e.simplify(tolerance);
197 197
    buf << "Coeff. of p2 should be 0";
198 198
    check(const_cast<const LpSolver::Expr&>(e)[p2]==0, buf.str());
199 199

	
200 200

	
201 201
  }
202 202

	
203 203
  {
204 204
    LP::DualExpr e,f,g;
205 205
    LP::Row p1 = INVALID, p2 = INVALID, p3 = INVALID,
206 206
      p4 = INVALID, p5 = INVALID;
207 207

	
208 208
    e[p1]=2;
209 209
    e[p1]+=2;
210 210
    e[p1]-=2;
211 211

	
212 212
    e=p1;
213 213
    e=f;
214 214

	
215 215
    e+=p1;
216 216
    e+=f;
217 217

	
218 218
    e-=p1;
219 219
    e-=f;
220 220

	
221 221
    e*=2;
222 222
    e*=2.2;
223 223
    e/=2;
224 224
    e/=2.2;
225 225

	
226 226
    e=((p1+p2)+(p1-p2)+
227 227
       (p1+f)+(f+p1)+(f+g)+
228 228
       (p1-f)+(f-p1)+(f-g)+
229 229
       2.2*f+f*2.2+f/2.2+
230 230
       2*f+f*2+f/2+
231 231
       2.2*p1+p1*2.2+p1/2.2+
232 232
       2*p1+p1*2+p1/2
233 233
       );
234 234
  }
235 235

	
236 236
}
237 237

	
238 238
void solveAndCheck(LpSolver& lp, LpSolver::ProblemType stat,
239 239
                   double exp_opt) {
240 240
  using std::string;
241 241
  lp.solve();
242 242

	
243 243
  std::ostringstream buf;
244 244
  buf << "PrimalType should be: " << int(stat) << int(lp.primalType());
245 245

	
246 246
  check(lp.primalType()==stat, buf.str());
247 247

	
248 248
  if (stat ==  LpSolver::OPTIMAL) {
249 249
    std::ostringstream sbuf;
250 250
    sbuf << "Wrong optimal value: the right optimum is " << exp_opt;
251 251
    check(std::abs(lp.primal()-exp_opt) < 1e-3, sbuf.str());
252 252
  }
253 253
}
254 254

	
255 255
void aTest(LpSolver & lp)
256 256
{
257 257
  typedef LpSolver LP;
258 258

	
259 259
 //The following example is very simple
260 260

	
261 261
  typedef LpSolver::Row Row;
262 262
  typedef LpSolver::Col Col;
263 263

	
264 264

	
265 265
  Col x1 = lp.addCol();
266 266
  Col x2 = lp.addCol();
267 267

	
268 268

	
269 269
  //Constraints
270 270
  Row upright=lp.addRow(x1+2*x2 <=1);
271 271
  lp.addRow(x1+x2 >=-1);
272 272
  lp.addRow(x1-x2 <=1);
273 273
  lp.addRow(x1-x2 >=-1);
274 274
  //Nonnegativity of the variables
275 275
  lp.colLowerBound(x1, 0);
276 276
  lp.colLowerBound(x2, 0);
277 277
  //Objective function
278 278
  lp.obj(x1+x2);
279 279

	
280 280
  lp.sense(lp.MAX);
281 281

	
282 282
  //Testing the problem retrieving routines
283 283
  check(lp.objCoeff(x1)==1,"First term should be 1 in the obj function!");
284 284
  check(lp.sense() == lp.MAX,"This is a maximization!");
285 285
  check(lp.coeff(upright,x1)==1,"The coefficient in question is 1!");
286 286
  check(lp.colLowerBound(x1)==0,
287 287
        "The lower bound for variable x1 should be 0.");
288 288
  check(lp.colUpperBound(x1)==LpSolver::INF,
289 289
        "The upper bound for variable x1 should be infty.");
290 290
  check(lp.rowLowerBound(upright) == -LpSolver::INF,
291 291
        "The lower bound for the first row should be -infty.");
292 292
  check(lp.rowUpperBound(upright)==1,
293 293
        "The upper bound for the first row should be 1.");
294 294
  LpSolver::Expr e = lp.row(upright);
295 295
  check(e[x1] == 1, "The first coefficient should 1.");
296 296
  check(e[x2] == 2, "The second coefficient should 1.");
297 297

	
298 298
  lp.row(upright, x1+x2 <=1);
299 299
  e = lp.row(upright);
300 300
  check(e[x1] == 1, "The first coefficient should 1.");
301 301
  check(e[x2] == 1, "The second coefficient should 1.");
302 302

	
303 303
  LpSolver::DualExpr de = lp.col(x1);
304 304
  check(  de[upright] == 1, "The first coefficient should 1.");
305 305

	
306 306
  LpSolver* clp = lp.cloneSolver();
307 307

	
308 308
  //Testing the problem retrieving routines
309 309
  check(clp->objCoeff(x1)==1,"First term should be 1 in the obj function!");
310 310
  check(clp->sense() == clp->MAX,"This is a maximization!");
311 311
  check(clp->coeff(upright,x1)==1,"The coefficient in question is 1!");
312 312
  //  std::cout<<lp.colLowerBound(x1)<<std::endl;
313 313
  check(clp->colLowerBound(x1)==0,
314 314
        "The lower bound for variable x1 should be 0.");
315 315
  check(clp->colUpperBound(x1)==LpSolver::INF,
316 316
        "The upper bound for variable x1 should be infty.");
317 317

	
318 318
  check(lp.rowLowerBound(upright)==-LpSolver::INF,
319 319
        "The lower bound for the first row should be -infty.");
320 320
  check(lp.rowUpperBound(upright)==1,
321 321
        "The upper bound for the first row should be 1.");
322 322
  e = clp->row(upright);
323 323
  check(e[x1] == 1, "The first coefficient should 1.");
324 324
  check(e[x2] == 1, "The second coefficient should 1.");
325 325

	
326 326
  de = clp->col(x1);
327 327
  check(de[upright] == 1, "The first coefficient should 1.");
328 328

	
329 329
  delete clp;
330 330

	
331 331
  //Maximization of x1+x2
332 332
  //over the triangle with vertices (0,0) (0,1) (1,0)
333 333
  double expected_opt=1;
334 334
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
335 335

	
336 336
  //Minimization
337 337
  lp.sense(lp.MIN);
338 338
  expected_opt=0;
339 339
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
340 340

	
341 341
  //Vertex (-1,0) instead of (0,0)
342 342
  lp.colLowerBound(x1, -LpSolver::INF);
343 343
  expected_opt=-1;
344 344
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
345 345

	
346 346
  //Erase one constraint and return to maximization
347 347
  lp.erase(upright);
348 348
  lp.sense(lp.MAX);
349 349
  expected_opt=LpSolver::INF;
350 350
  solveAndCheck(lp, LpSolver::UNBOUNDED, expected_opt);
351 351

	
352 352
  //Infeasibilty
353 353
  lp.addRow(x1+x2 <=-2);
354 354
  solveAndCheck(lp, LpSolver::INFEASIBLE, expected_opt);
355 355

	
356 356
}
357 357

	
358 358
int main()
359 359
{
360 360
  LpSkeleton lp_skel;
361 361
  lpTest(lp_skel);
362 362

	
363 363
#ifdef HAVE_GLPK
364 364
  {
365 365
    GlpkLp lp_glpk1,lp_glpk2;
366 366
    lpTest(lp_glpk1);
367 367
    aTest(lp_glpk2);
368 368
  }
369
  GlpkLp::freeEnv();
369 370
#endif
370 371

	
371 372
#ifdef HAVE_CPLEX
372 373
  try {
373 374
    CplexLp lp_cplex1,lp_cplex2;
374 375
    lpTest(lp_cplex1);
375 376
    aTest(lp_cplex2);
376 377
  } catch (CplexEnv::LicenseError& error) {
377 378
#ifdef LEMON_FORCE_CPLEX_CHECK
378 379
    check(false, error.what());
379 380
#else
380 381
    std::cerr << error.what() << std::endl;
381 382
    std::cerr << "Cplex license check failed, lp check skipped" << std::endl;
382 383
#endif
383 384
  }
384 385
#endif
385 386

	
386 387
#ifdef HAVE_SOPLEX
387 388
  {
388 389
    SoplexLp lp_soplex1,lp_soplex2;
389 390
    lpTest(lp_soplex1);
390 391
    aTest(lp_soplex2);
391 392
  }
392 393
#endif
393 394

	
394 395
#ifdef HAVE_CLP
395 396
  {
396 397
    ClpLp lp_clp1,lp_clp2;
397 398
    lpTest(lp_clp1);
398 399
    aTest(lp_clp2);
399 400
  }
400 401
#endif
401 402

	
402 403
  return 0;
403 404
}
Ignore white space 768 line context
1 1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
2 2
 *
3 3
 * This file is a part of LEMON, a generic C++ optimization library.
4 4
 *
5 5
 * Copyright (C) 2003-2008
6 6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 8
 *
9 9
 * Permission to use, modify and distribute this software is granted
10 10
 * provided that this copyright notice appears in all copies. For
11 11
 * precise terms see the accompanying LICENSE file.
12 12
 *
13 13
 * This software is provided "AS IS" with no warranty of any kind,
14 14
 * express or implied, and with no claim as to its suitability for any
15 15
 * purpose.
16 16
 *
17 17
 */
18 18

	
19 19
#include "test_tools.h"
20 20

	
21 21

	
22 22
#ifdef HAVE_CONFIG_H
23 23
#include <lemon/config.h>
24 24
#endif
25 25

	
26 26
#ifdef HAVE_CPLEX
27 27
#include <lemon/cplex.h>
28 28
#endif
29 29

	
30 30
#ifdef HAVE_GLPK
31 31
#include <lemon/glpk.h>
32 32
#endif
33 33

	
34 34

	
35 35
using namespace lemon;
36 36

	
37 37
void solveAndCheck(MipSolver& mip, MipSolver::ProblemType stat,
38 38
                   double exp_opt) {
39 39
  using std::string;
40 40

	
41 41
  mip.solve();
42 42
  //int decimal,sign;
43 43
  std::ostringstream buf;
44 44
  buf << "Type should be: " << int(stat)<<" and it is "<<int(mip.type());
45 45

	
46 46

	
47 47
  //  itoa(stat,buf1, 10);
48 48
  check(mip.type()==stat, buf.str());
49 49

	
50 50
  if (stat ==  MipSolver::OPTIMAL) {
51 51
    std::ostringstream sbuf;
52 52
    buf << "Wrong optimal value: the right optimum is " << exp_opt;
53 53
    check(std::abs(mip.solValue()-exp_opt) < 1e-3, sbuf.str());
54 54
    //+ecvt(exp_opt,2)
55 55
  }
56 56
}
57 57

	
58 58
void aTest(MipSolver& mip)
59 59
{
60 60
 //The following example is very simple
61 61

	
62 62

	
63 63
  typedef MipSolver::Row Row;
64 64
  typedef MipSolver::Col Col;
65 65

	
66 66

	
67 67

	
68 68
  Col x1 = mip.addCol();
69 69
  Col x2 = mip.addCol();
70 70

	
71 71

	
72 72
  //Objective function
73 73
  mip.obj(x1);
74 74

	
75 75
  mip.max();
76 76

	
77 77

	
78 78
  //Unconstrained optimization
79 79
  mip.solve();
80 80
  //Check it out!
81 81

	
82 82
  //Constraints
83 83
  mip.addRow(2*x1+x2 <=2);
84 84
  mip.addRow(x1-2*x2 <=0);
85 85

	
86 86
  //Nonnegativity of the variable x1
87 87
  mip.colLowerBound(x1, 0);
88 88

	
89 89
  //Maximization of x1
90 90
  //over the triangle with vertices (0,0),(4/5,2/5),(0,2)
91 91
  double expected_opt=4.0/5.0;
92 92
  solveAndCheck(mip, MipSolver::OPTIMAL, expected_opt);
93 93

	
94 94
  //Restrict x2 to integer
95 95
  mip.colType(x2,MipSolver::INTEGER);
96 96
  expected_opt=1.0/2.0;
97 97
  solveAndCheck(mip, MipSolver::OPTIMAL, expected_opt);
98 98

	
99 99

	
100 100
  //Restrict both to integer
101 101
  mip.colType(x1,MipSolver::INTEGER);
102 102
  expected_opt=0;
103 103
  solveAndCheck(mip, MipSolver::OPTIMAL, expected_opt);
104 104

	
105 105

	
106 106

	
107 107
}
108 108

	
109 109

	
110 110
int main()
111 111
{
112 112

	
113 113
#ifdef HAVE_GLPK
114 114
  {
115 115
    GlpkMip mip1;
116 116
    aTest(mip1);
117 117
  }
118
  GlpkLp::freeEnv();
118 119
#endif
119 120

	
120 121
#ifdef HAVE_CPLEX
121 122
  try {
122 123
    CplexMip mip2;
123 124
    aTest(mip2);
124 125
  } catch (CplexEnv::LicenseError& error) {
125 126
#ifdef LEMON_FORCE_CPLEX_CHECK
126 127
    check(false, error.what());
127 128
#else
128 129
    std::cerr << error.what() << std::endl;
129 130
    std::cerr << "Cplex license check failed, lp check skipped" << std::endl;
130 131
#endif
131 132
  }
132 133
#endif
133 134

	
134 135
  return 0;
135 136

	
136 137
}
0 comments (0 inline)