gravatar
deba@inf.elte.hu
deba@inf.elte.hu
Fix lp related errors and warnings (#241 and #242)
0 5 0
default
5 files changed with 25 insertions and 22 deletions:
↑ Collapse diff ↑
Ignore white space 384 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
 * Copyright (C) 2003-2008
5
 * Copyright (C) 2003-2009
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 <iostream>
20 20
#include <vector>
21 21
#include <cstring>
22 22

	
23 23
#include <lemon/cplex.h>
24 24

	
25 25
extern "C" {
26 26
#include <ilcplex/cplex.h>
27 27
}
28 28

	
29 29

	
30 30
///\file
31 31
///\brief Implementation of the LEMON-CPLEX lp solver interface.
32 32
namespace lemon {
33 33

	
34 34
  CplexEnv::LicenseError::LicenseError(int status) {
35 35
    if (!CPXgeterrorstring(0, status, _message)) {
36 36
      std::strcpy(_message, "Cplex unknown error");
37 37
    }
38 38
  }
39 39

	
40 40
  CplexEnv::CplexEnv() {
41 41
    int status;
42 42
    _cnt = new int;
43 43
    _env = CPXopenCPLEX(&status);
44 44
    if (_env == 0) {
45 45
      delete _cnt;
46 46
      _cnt = 0;
47 47
      throw LicenseError(status);
48 48
    }
49 49
  }
50 50

	
51 51
  CplexEnv::CplexEnv(const CplexEnv& other) {
52 52
    _env = other._env;
53 53
    _cnt = other._cnt;
54 54
    ++(*_cnt);
55 55
  }
56 56

	
57 57
  CplexEnv& CplexEnv::operator=(const CplexEnv& other) {
58 58
    _env = other._env;
59 59
    _cnt = other._cnt;
60 60
    ++(*_cnt);
61 61
    return *this;
62 62
  }
63 63

	
64 64
  CplexEnv::~CplexEnv() {
65 65
    --(*_cnt);
66 66
    if (*_cnt == 0) {
67 67
      delete _cnt;
68 68
      CPXcloseCPLEX(&_env);
69 69
    }
70 70
  }
71 71

	
72 72
  CplexBase::CplexBase() : LpBase() {
73 73
    int status;
74 74
    _prob = CPXcreateprob(cplexEnv(), &status, "Cplex problem");
75 75
  }
76 76

	
77 77
  CplexBase::CplexBase(const CplexEnv& env)
78 78
    : LpBase(), _env(env) {
79 79
    int status;
80 80
    _prob = CPXcreateprob(cplexEnv(), &status, "Cplex problem");
81 81
  }
82 82

	
83 83
  CplexBase::CplexBase(const CplexBase& cplex)
84 84
    : LpBase() {
85 85
    int status;
86 86
    _prob = CPXcloneprob(cplexEnv(), cplex._prob, &status);
87 87
    rows = cplex.rows;
88 88
    cols = cplex.cols;
89 89
  }
90 90

	
91 91
  CplexBase::~CplexBase() {
92 92
    CPXfreeprob(cplexEnv(),&_prob);
93 93
  }
94 94

	
95 95
  int CplexBase::_addCol() {
96 96
    int i = CPXgetnumcols(cplexEnv(), _prob);
97 97
    double lb = -INF, ub = INF;
98 98
    CPXnewcols(cplexEnv(), _prob, 1, 0, &lb, &ub, 0, 0);
99 99
    return i;
100 100
  }
101 101

	
102 102

	
103 103
  int CplexBase::_addRow() {
104 104
    int i = CPXgetnumrows(cplexEnv(), _prob);
105 105
    const double ub = INF;
106 106
    const char s = 'L';
107 107
    CPXnewrows(cplexEnv(), _prob, 1, &ub, &s, 0, 0);
108 108
    return i;
109 109
  }
110 110

	
111 111

	
112 112
  void CplexBase::_eraseCol(int i) {
113 113
    CPXdelcols(cplexEnv(), _prob, i, i);
114 114
  }
115 115

	
116 116
  void CplexBase::_eraseRow(int i) {
117 117
    CPXdelrows(cplexEnv(), _prob, i, i);
118 118
  }
119 119

	
120 120
  void CplexBase::_eraseColId(int i) {
121 121
    cols.eraseIndex(i);
122 122
    cols.shiftIndices(i);
123 123
  }
124 124
  void CplexBase::_eraseRowId(int i) {
125 125
    rows.eraseIndex(i);
126 126
    rows.shiftIndices(i);
127 127
  }
128 128

	
129 129
  void CplexBase::_getColName(int col, std::string &name) const {
130 130
    int size;
131 131
    CPXgetcolname(cplexEnv(), _prob, 0, 0, 0, &size, col, col);
132 132
    if (size == 0) {
133 133
      name.clear();
134 134
      return;
135 135
    }
136 136

	
137 137
    size *= -1;
138 138
    std::vector<char> buf(size);
139 139
    char *cname;
140 140
    int tmp;
141 141
    CPXgetcolname(cplexEnv(), _prob, &cname, &buf.front(), size,
142 142
                  &tmp, col, col);
143 143
    name = cname;
144 144
  }
145 145

	
146 146
  void CplexBase::_setColName(int col, const std::string &name) {
147 147
    char *cname;
148 148
    cname = const_cast<char*>(name.c_str());
149 149
    CPXchgcolname(cplexEnv(), _prob, 1, &col, &cname);
150 150
  }
151 151

	
152 152
  int CplexBase::_colByName(const std::string& name) const {
153 153
    int index;
154 154
    if (CPXgetcolindex(cplexEnv(), _prob,
155 155
                       const_cast<char*>(name.c_str()), &index) == 0) {
156 156
      return index;
157 157
    }
158 158
    return -1;
159 159
  }
160 160

	
161 161
  void CplexBase::_getRowName(int row, std::string &name) const {
162 162
    int size;
163 163
    CPXgetrowname(cplexEnv(), _prob, 0, 0, 0, &size, row, row);
164 164
    if (size == 0) {
165 165
      name.clear();
166 166
      return;
167 167
    }
168 168

	
169 169
    size *= -1;
170 170
    std::vector<char> buf(size);
171 171
    char *cname;
172 172
    int tmp;
173 173
    CPXgetrowname(cplexEnv(), _prob, &cname, &buf.front(), size,
174 174
                  &tmp, row, row);
175 175
    name = cname;
176 176
  }
177 177

	
178 178
  void CplexBase::_setRowName(int row, const std::string &name) {
179 179
    char *cname;
180 180
    cname = const_cast<char*>(name.c_str());
181 181
    CPXchgrowname(cplexEnv(), _prob, 1, &row, &cname);
182 182
  }
183 183

	
184 184
  int CplexBase::_rowByName(const std::string& name) const {
185 185
    int index;
186 186
    if (CPXgetrowindex(cplexEnv(), _prob,
187 187
                       const_cast<char*>(name.c_str()), &index) == 0) {
188 188
      return index;
189 189
    }
190 190
    return -1;
191 191
  }
192 192

	
193 193
  void CplexBase::_setRowCoeffs(int i, ExprIterator b,
194 194
                                      ExprIterator e)
195 195
  {
196 196
    std::vector<int> indices;
197 197
    std::vector<int> rowlist;
... ...
@@ -252,674 +252,674 @@
252 252
    CPXgetcols(cplexEnv(), _prob, &tmp1, &tmp2,
253 253
               &indices.front(), &values.front(),
254 254
               length, &tmp3, i, i);
255 255

	
256 256
    for (int i = 0; i < length; ++i) {
257 257
      *b = std::make_pair(indices[i], values[i]);
258 258
      ++b;
259 259
    }
260 260

	
261 261
  }
262 262

	
263 263
  void CplexBase::_setCoeff(int row, int col, Value value) {
264 264
    CPXchgcoef(cplexEnv(), _prob, row, col, value);
265 265
  }
266 266

	
267 267
  CplexBase::Value CplexBase::_getCoeff(int row, int col) const {
268 268
    CplexBase::Value value;
269 269
    CPXgetcoef(cplexEnv(), _prob, row, col, &value);
270 270
    return value;
271 271
  }
272 272

	
273 273
  void CplexBase::_setColLowerBound(int i, Value value) {
274 274
    const char s = 'L';
275 275
    CPXchgbds(cplexEnv(), _prob, 1, &i, &s, &value);
276 276
  }
277 277

	
278 278
  CplexBase::Value CplexBase::_getColLowerBound(int i) const {
279 279
    CplexBase::Value res;
280 280
    CPXgetlb(cplexEnv(), _prob, &res, i, i);
281 281
    return res <= -CPX_INFBOUND ? -INF : res;
282 282
  }
283 283

	
284 284
  void CplexBase::_setColUpperBound(int i, Value value)
285 285
  {
286 286
    const char s = 'U';
287 287
    CPXchgbds(cplexEnv(), _prob, 1, &i, &s, &value);
288 288
  }
289 289

	
290 290
  CplexBase::Value CplexBase::_getColUpperBound(int i) const {
291 291
    CplexBase::Value res;
292 292
    CPXgetub(cplexEnv(), _prob, &res, i, i);
293 293
    return res >= CPX_INFBOUND ? INF : res;
294 294
  }
295 295

	
296 296
  CplexBase::Value CplexBase::_getRowLowerBound(int i) const {
297 297
    char s;
298 298
    CPXgetsense(cplexEnv(), _prob, &s, i, i);
299 299
    CplexBase::Value res;
300 300

	
301 301
    switch (s) {
302 302
    case 'G':
303 303
    case 'R':
304 304
    case 'E':
305 305
      CPXgetrhs(cplexEnv(), _prob, &res, i, i);
306 306
      return res <= -CPX_INFBOUND ? -INF : res;
307 307
    default:
308 308
      return -INF;
309 309
    }
310 310
  }
311 311

	
312 312
  CplexBase::Value CplexBase::_getRowUpperBound(int i) const {
313 313
    char s;
314 314
    CPXgetsense(cplexEnv(), _prob, &s, i, i);
315 315
    CplexBase::Value res;
316 316

	
317 317
    switch (s) {
318 318
    case 'L':
319 319
    case 'E':
320 320
      CPXgetrhs(cplexEnv(), _prob, &res, i, i);
321 321
      return res >= CPX_INFBOUND ? INF : res;
322 322
    case 'R':
323 323
      CPXgetrhs(cplexEnv(), _prob, &res, i, i);
324 324
      {
325 325
        double rng;
326 326
        CPXgetrngval(cplexEnv(), _prob, &rng, i, i);
327 327
        res += rng;
328 328
      }
329 329
      return res >= CPX_INFBOUND ? INF : res;
330 330
    default:
331 331
      return INF;
332 332
    }
333 333
  }
334 334

	
335 335
  //This is easier to implement
336 336
  void CplexBase::_set_row_bounds(int i, Value lb, Value ub) {
337 337
    if (lb == -INF) {
338 338
      const char s = 'L';
339 339
      CPXchgsense(cplexEnv(), _prob, 1, &i, &s);
340 340
      CPXchgrhs(cplexEnv(), _prob, 1, &i, &ub);
341 341
    } else if (ub == INF) {
342 342
      const char s = 'G';
343 343
      CPXchgsense(cplexEnv(), _prob, 1, &i, &s);
344 344
      CPXchgrhs(cplexEnv(), _prob, 1, &i, &lb);
345 345
    } else if (lb == ub){
346 346
      const char s = 'E';
347 347
      CPXchgsense(cplexEnv(), _prob, 1, &i, &s);
348 348
      CPXchgrhs(cplexEnv(), _prob, 1, &i, &lb);
349 349
    } else {
350 350
      const char s = 'R';
351 351
      CPXchgsense(cplexEnv(), _prob, 1, &i, &s);
352 352
      CPXchgrhs(cplexEnv(), _prob, 1, &i, &lb);
353 353
      double len = ub - lb;
354 354
      CPXchgrngval(cplexEnv(), _prob, 1, &i, &len);
355 355
    }
356 356
  }
357 357

	
358 358
  void CplexBase::_setRowLowerBound(int i, Value lb)
359 359
  {
360 360
    LEMON_ASSERT(lb != INF, "Invalid bound");
361 361
    _set_row_bounds(i, lb, CplexBase::_getRowUpperBound(i));
362 362
  }
363 363

	
364 364
  void CplexBase::_setRowUpperBound(int i, Value ub)
365 365
  {
366 366

	
367 367
    LEMON_ASSERT(ub != -INF, "Invalid bound");
368 368
    _set_row_bounds(i, CplexBase::_getRowLowerBound(i), ub);
369 369
  }
370 370

	
371 371
  void CplexBase::_setObjCoeffs(ExprIterator b, ExprIterator e)
372 372
  {
373 373
    std::vector<int> indices;
374 374
    std::vector<Value> values;
375 375
    for(ExprIterator it=b; it!=e; ++it) {
376 376
      indices.push_back(it->first);
377 377
      values.push_back(it->second);
378 378
    }
379 379
    CPXchgobj(cplexEnv(), _prob, values.size(),
380 380
              &indices.front(), &values.front());
381 381

	
382 382
  }
383 383

	
384 384
  void CplexBase::_getObjCoeffs(InsertIterator b) const
385 385
  {
386 386
    int num = CPXgetnumcols(cplexEnv(), _prob);
387 387
    std::vector<Value> x(num);
388 388

	
389 389
    CPXgetobj(cplexEnv(), _prob, &x.front(), 0, num - 1);
390 390
    for (int i = 0; i < num; ++i) {
391 391
      if (x[i] != 0.0) {
392 392
        *b = std::make_pair(i, x[i]);
393 393
        ++b;
394 394
      }
395 395
    }
396 396
  }
397 397

	
398 398
  void CplexBase::_setObjCoeff(int i, Value obj_coef)
399 399
  {
400 400
    CPXchgobj(cplexEnv(), _prob, 1, &i, &obj_coef);
401 401
  }
402 402

	
403 403
  CplexBase::Value CplexBase::_getObjCoeff(int i) const
404 404
  {
405 405
    Value x;
406 406
    CPXgetobj(cplexEnv(), _prob, &x, i, i);
407 407
    return x;
408 408
  }
409 409

	
410 410
  void CplexBase::_setSense(CplexBase::Sense sense) {
411 411
    switch (sense) {
412 412
    case MIN:
413 413
      CPXchgobjsen(cplexEnv(), _prob, CPX_MIN);
414 414
      break;
415 415
    case MAX:
416 416
      CPXchgobjsen(cplexEnv(), _prob, CPX_MAX);
417 417
      break;
418 418
    }
419 419
  }
420 420

	
421 421
  CplexBase::Sense CplexBase::_getSense() const {
422 422
    switch (CPXgetobjsen(cplexEnv(), _prob)) {
423 423
    case CPX_MIN:
424 424
      return MIN;
425 425
    case CPX_MAX:
426 426
      return MAX;
427 427
    default:
428 428
      LEMON_ASSERT(false, "Invalid sense");
429 429
      return CplexBase::Sense();
430 430
    }
431 431
  }
432 432

	
433 433
  void CplexBase::_clear() {
434 434
    CPXfreeprob(cplexEnv(),&_prob);
435 435
    int status;
436 436
    _prob = CPXcreateprob(cplexEnv(), &status, "Cplex problem");
437 437
    rows.clear();
438 438
    cols.clear();
439 439
  }
440 440

	
441 441
  // CplexLp members
442 442

	
443 443
  CplexLp::CplexLp()
444
    : LpBase(), CplexBase(), LpSolver() {}
444
    : LpBase(), LpSolver(), CplexBase() {}
445 445

	
446 446
  CplexLp::CplexLp(const CplexEnv& env)
447
    : LpBase(), CplexBase(env), LpSolver() {}
447
    : LpBase(), LpSolver(), CplexBase(env) {}
448 448

	
449 449
  CplexLp::CplexLp(const CplexLp& other)
450
    : LpBase(), CplexBase(other), LpSolver() {}
450
    : LpBase(), LpSolver(), CplexBase(other) {}
451 451

	
452 452
  CplexLp::~CplexLp() {}
453 453

	
454 454
  CplexLp* CplexLp::newSolver() const { return new CplexLp; }
455 455
  CplexLp* CplexLp::cloneSolver() const {return new CplexLp(*this); }
456 456

	
457 457
  const char* CplexLp::_solverName() const { return "CplexLp"; }
458 458

	
459 459
  void CplexLp::_clear_temporals() {
460 460
    _col_status.clear();
461 461
    _row_status.clear();
462 462
    _primal_ray.clear();
463 463
    _dual_ray.clear();
464 464
  }
465 465

	
466 466
  // The routine returns zero unless an error occurred during the
467 467
  // optimization. Examples of errors include exhausting available
468 468
  // memory (CPXERR_NO_MEMORY) or encountering invalid data in the
469 469
  // CPLEX problem object (CPXERR_NO_PROBLEM). Exceeding a
470 470
  // user-specified CPLEX limit, or proving the model infeasible or
471 471
  // unbounded, are not considered errors. Note that a zero return
472 472
  // value does not necessarily mean that a solution exists. Use query
473 473
  // routines CPXsolninfo, CPXgetstat, and CPXsolution to obtain
474 474
  // further information about the status of the optimization.
475 475
  CplexLp::SolveExitStatus CplexLp::convertStatus(int status) {
476 476
#if CPX_VERSION >= 800
477 477
    if (status == 0) {
478 478
      switch (CPXgetstat(cplexEnv(), _prob)) {
479 479
      case CPX_STAT_OPTIMAL:
480 480
      case CPX_STAT_INFEASIBLE:
481 481
      case CPX_STAT_UNBOUNDED:
482 482
        return SOLVED;
483 483
      default:
484 484
        return UNSOLVED;
485 485
      }
486 486
    } else {
487 487
      return UNSOLVED;
488 488
    }
489 489
#else
490 490
    if (status == 0) {
491 491
      //We want to exclude some cases
492 492
      switch (CPXgetstat(cplexEnv(), _prob)) {
493 493
      case CPX_OBJ_LIM:
494 494
      case CPX_IT_LIM_FEAS:
495 495
      case CPX_IT_LIM_INFEAS:
496 496
      case CPX_TIME_LIM_FEAS:
497 497
      case CPX_TIME_LIM_INFEAS:
498 498
        return UNSOLVED;
499 499
      default:
500 500
        return SOLVED;
501 501
      }
502 502
    } else {
503 503
      return UNSOLVED;
504 504
    }
505 505
#endif
506 506
  }
507 507

	
508 508
  CplexLp::SolveExitStatus CplexLp::_solve() {
509 509
    _clear_temporals();
510 510
    return convertStatus(CPXlpopt(cplexEnv(), _prob));
511 511
  }
512 512

	
513 513
  CplexLp::SolveExitStatus CplexLp::solvePrimal() {
514 514
    _clear_temporals();
515 515
    return convertStatus(CPXprimopt(cplexEnv(), _prob));
516 516
  }
517 517

	
518 518
  CplexLp::SolveExitStatus CplexLp::solveDual() {
519 519
    _clear_temporals();
520 520
    return convertStatus(CPXdualopt(cplexEnv(), _prob));
521 521
  }
522 522

	
523 523
  CplexLp::SolveExitStatus CplexLp::solveBarrier() {
524 524
    _clear_temporals();
525 525
    return convertStatus(CPXbaropt(cplexEnv(), _prob));
526 526
  }
527 527

	
528 528
  CplexLp::Value CplexLp::_getPrimal(int i) const {
529 529
    Value x;
530 530
    CPXgetx(cplexEnv(), _prob, &x, i, i);
531 531
    return x;
532 532
  }
533 533

	
534 534
  CplexLp::Value CplexLp::_getDual(int i) const {
535 535
    Value y;
536 536
    CPXgetpi(cplexEnv(), _prob, &y, i, i);
537 537
    return y;
538 538
  }
539 539

	
540 540
  CplexLp::Value CplexLp::_getPrimalValue() const {
541 541
    Value objval;
542 542
    CPXgetobjval(cplexEnv(), _prob, &objval);
543 543
    return objval;
544 544
  }
545 545

	
546 546
  CplexLp::VarStatus CplexLp::_getColStatus(int i) const {
547 547
    if (_col_status.empty()) {
548 548
      _col_status.resize(CPXgetnumcols(cplexEnv(), _prob));
549 549
      CPXgetbase(cplexEnv(), _prob, &_col_status.front(), 0);
550 550
    }
551 551
    switch (_col_status[i]) {
552 552
    case CPX_BASIC:
553 553
      return BASIC;
554 554
    case CPX_FREE_SUPER:
555 555
      return FREE;
556 556
    case CPX_AT_LOWER:
557 557
      return LOWER;
558 558
    case CPX_AT_UPPER:
559 559
      return UPPER;
560 560
    default:
561 561
      LEMON_ASSERT(false, "Wrong column status");
562 562
      return CplexLp::VarStatus();
563 563
    }
564 564
  }
565 565

	
566 566
  CplexLp::VarStatus CplexLp::_getRowStatus(int i) const {
567 567
    if (_row_status.empty()) {
568 568
      _row_status.resize(CPXgetnumrows(cplexEnv(), _prob));
569 569
      CPXgetbase(cplexEnv(), _prob, 0, &_row_status.front());
570 570
    }
571 571
    switch (_row_status[i]) {
572 572
    case CPX_BASIC:
573 573
      return BASIC;
574 574
    case CPX_AT_LOWER:
575 575
      {
576 576
        char s;
577 577
        CPXgetsense(cplexEnv(), _prob, &s, i, i);
578 578
        return s != 'L' ? LOWER : UPPER;
579 579
      }
580 580
    case CPX_AT_UPPER:
581 581
      return UPPER;
582 582
    default:
583 583
      LEMON_ASSERT(false, "Wrong row status");
584 584
      return CplexLp::VarStatus();
585 585
    }
586 586
  }
587 587

	
588 588
  CplexLp::Value CplexLp::_getPrimalRay(int i) const {
589 589
    if (_primal_ray.empty()) {
590 590
      _primal_ray.resize(CPXgetnumcols(cplexEnv(), _prob));
591 591
      CPXgetray(cplexEnv(), _prob, &_primal_ray.front());
592 592
    }
593 593
    return _primal_ray[i];
594 594
  }
595 595

	
596 596
  CplexLp::Value CplexLp::_getDualRay(int i) const {
597 597
    if (_dual_ray.empty()) {
598 598

	
599 599
    }
600 600
    return _dual_ray[i];
601 601
  }
602 602

	
603 603
  //7.5-os cplex statusai (Vigyazat: a 9.0-asei masok!)
604 604
  // This table lists the statuses, returned by the CPXgetstat()
605 605
  // routine, for solutions to LP problems or mixed integer problems. If
606 606
  // no solution exists, the return value is zero.
607 607

	
608 608
  // For Simplex, Barrier
609 609
  // 1          CPX_OPTIMAL
610 610
  //          Optimal solution found
611 611
  // 2          CPX_INFEASIBLE
612 612
  //          Problem infeasible
613 613
  // 3    CPX_UNBOUNDED
614 614
  //          Problem unbounded
615 615
  // 4          CPX_OBJ_LIM
616 616
  //          Objective limit exceeded in Phase II
617 617
  // 5          CPX_IT_LIM_FEAS
618 618
  //          Iteration limit exceeded in Phase II
619 619
  // 6          CPX_IT_LIM_INFEAS
620 620
  //          Iteration limit exceeded in Phase I
621 621
  // 7          CPX_TIME_LIM_FEAS
622 622
  //          Time limit exceeded in Phase II
623 623
  // 8          CPX_TIME_LIM_INFEAS
624 624
  //          Time limit exceeded in Phase I
625 625
  // 9          CPX_NUM_BEST_FEAS
626 626
  //          Problem non-optimal, singularities in Phase II
627 627
  // 10         CPX_NUM_BEST_INFEAS
628 628
  //          Problem non-optimal, singularities in Phase I
629 629
  // 11         CPX_OPTIMAL_INFEAS
630 630
  //          Optimal solution found, unscaled infeasibilities
631 631
  // 12         CPX_ABORT_FEAS
632 632
  //          Aborted in Phase II
633 633
  // 13         CPX_ABORT_INFEAS
634 634
  //          Aborted in Phase I
635 635
  // 14          CPX_ABORT_DUAL_INFEAS
636 636
  //          Aborted in barrier, dual infeasible
637 637
  // 15          CPX_ABORT_PRIM_INFEAS
638 638
  //          Aborted in barrier, primal infeasible
639 639
  // 16          CPX_ABORT_PRIM_DUAL_INFEAS
640 640
  //          Aborted in barrier, primal and dual infeasible
641 641
  // 17          CPX_ABORT_PRIM_DUAL_FEAS
642 642
  //          Aborted in barrier, primal and dual feasible
643 643
  // 18          CPX_ABORT_CROSSOVER
644 644
  //          Aborted in crossover
645 645
  // 19          CPX_INForUNBD
646 646
  //          Infeasible or unbounded
647 647
  // 20   CPX_PIVOT
648 648
  //       User pivot used
649 649
  //
650 650
  //     Ezeket hova tegyem:
651 651
  // ??case CPX_ABORT_DUAL_INFEAS
652 652
  // ??case CPX_ABORT_CROSSOVER
653 653
  // ??case CPX_INForUNBD
654 654
  // ??case CPX_PIVOT
655 655

	
656 656
  //Some more interesting stuff:
657 657

	
658 658
  // CPX_PARAM_PROBMETHOD  1062  int  LPMETHOD
659 659
  // 0 Automatic
660 660
  // 1 Primal Simplex
661 661
  // 2 Dual Simplex
662 662
  // 3 Network Simplex
663 663
  // 4 Standard Barrier
664 664
  // Default: 0
665 665
  // Description: Method for linear optimization.
666 666
  // Determines which algorithm is used when CPXlpopt() (or "optimize"
667 667
  // in the Interactive Optimizer) is called. Currently the behavior of
668 668
  // the "Automatic" setting is that CPLEX simply invokes the dual
669 669
  // simplex method, but this capability may be expanded in the future
670 670
  // so that CPLEX chooses the method based on problem characteristics
671 671
#if CPX_VERSION < 900
672 672
  void statusSwitch(CPXENVptr cplexEnv(),int& stat){
673 673
    int lpmethod;
674 674
    CPXgetintparam (cplexEnv(),CPX_PARAM_PROBMETHOD,&lpmethod);
675 675
    if (lpmethod==2){
676 676
      if (stat==CPX_UNBOUNDED){
677 677
        stat=CPX_INFEASIBLE;
678 678
      }
679 679
      else{
680 680
        if (stat==CPX_INFEASIBLE)
681 681
          stat=CPX_UNBOUNDED;
682 682
      }
683 683
    }
684 684
  }
685 685
#else
686 686
  void statusSwitch(CPXENVptr,int&){}
687 687
#endif
688 688

	
689 689
  CplexLp::ProblemType CplexLp::_getPrimalType() const {
690 690
    // Unboundedness not treated well: the following is from cplex 9.0 doc
691 691
    // About Unboundedness
692 692

	
693 693
    // The treatment of models that are unbounded involves a few
694 694
    // subtleties. Specifically, a declaration of unboundedness means that
695 695
    // ILOG CPLEX has determined that the model has an unbounded
696 696
    // ray. Given any feasible solution x with objective z, a multiple of
697 697
    // the unbounded ray can be added to x to give a feasible solution
698 698
    // with objective z-1 (or z+1 for maximization models). Thus, if a
699 699
    // feasible solution exists, then the optimal objective is
700 700
    // unbounded. Note that ILOG CPLEX has not necessarily concluded that
701 701
    // a feasible solution exists. Users can call the routine CPXsolninfo
702 702
    // to determine whether ILOG CPLEX has also concluded that the model
703 703
    // has a feasible solution.
704 704

	
705 705
    int stat = CPXgetstat(cplexEnv(), _prob);
706 706
#if CPX_VERSION >= 800
707 707
    switch (stat)
708 708
      {
709 709
      case CPX_STAT_OPTIMAL:
710 710
        return OPTIMAL;
711 711
      case CPX_STAT_UNBOUNDED:
712 712
        return UNBOUNDED;
713 713
      case CPX_STAT_INFEASIBLE:
714 714
        return INFEASIBLE;
715 715
      default:
716 716
        return UNDEFINED;
717 717
      }
718 718
#else
719 719
    statusSwitch(cplexEnv(),stat);
720 720
    //CPXgetstat(cplexEnv(), _prob);
721 721
    //printf("A primal status: %d, CPX_OPTIMAL=%d \n",stat,CPX_OPTIMAL);
722 722
    switch (stat) {
723 723
    case 0:
724 724
      return UNDEFINED; //Undefined
725 725
    case CPX_OPTIMAL://Optimal
726 726
      return OPTIMAL;
727 727
    case CPX_UNBOUNDED://Unbounded
728 728
      return INFEASIBLE;//In case of dual simplex
729 729
      //return UNBOUNDED;
730 730
    case CPX_INFEASIBLE://Infeasible
731 731
      //    case CPX_IT_LIM_INFEAS:
732 732
      //     case CPX_TIME_LIM_INFEAS:
733 733
      //     case CPX_NUM_BEST_INFEAS:
734 734
      //     case CPX_OPTIMAL_INFEAS:
735 735
      //     case CPX_ABORT_INFEAS:
736 736
      //     case CPX_ABORT_PRIM_INFEAS:
737 737
      //     case CPX_ABORT_PRIM_DUAL_INFEAS:
738 738
      return UNBOUNDED;//In case of dual simplex
739 739
      //return INFEASIBLE;
740 740
      //     case CPX_OBJ_LIM:
741 741
      //     case CPX_IT_LIM_FEAS:
742 742
      //     case CPX_TIME_LIM_FEAS:
743 743
      //     case CPX_NUM_BEST_FEAS:
744 744
      //     case CPX_ABORT_FEAS:
745 745
      //     case CPX_ABORT_PRIM_DUAL_FEAS:
746 746
      //       return FEASIBLE;
747 747
    default:
748 748
      return UNDEFINED; //Everything else comes here
749 749
      //FIXME error
750 750
    }
751 751
#endif
752 752
  }
753 753

	
754 754
  //9.0-as cplex verzio statusai
755 755
  // CPX_STAT_ABORT_DUAL_OBJ_LIM
756 756
  // CPX_STAT_ABORT_IT_LIM
757 757
  // CPX_STAT_ABORT_OBJ_LIM
758 758
  // CPX_STAT_ABORT_PRIM_OBJ_LIM
759 759
  // CPX_STAT_ABORT_TIME_LIM
760 760
  // CPX_STAT_ABORT_USER
761 761
  // CPX_STAT_FEASIBLE_RELAXED
762 762
  // CPX_STAT_INFEASIBLE
763 763
  // CPX_STAT_INForUNBD
764 764
  // CPX_STAT_NUM_BEST
765 765
  // CPX_STAT_OPTIMAL
766 766
  // CPX_STAT_OPTIMAL_FACE_UNBOUNDED
767 767
  // CPX_STAT_OPTIMAL_INFEAS
768 768
  // CPX_STAT_OPTIMAL_RELAXED
769 769
  // CPX_STAT_UNBOUNDED
770 770

	
771 771
  CplexLp::ProblemType CplexLp::_getDualType() const {
772 772
    int stat = CPXgetstat(cplexEnv(), _prob);
773 773
#if CPX_VERSION >= 800
774 774
    switch (stat) {
775 775
    case CPX_STAT_OPTIMAL:
776 776
      return OPTIMAL;
777 777
    case CPX_STAT_UNBOUNDED:
778 778
      return INFEASIBLE;
779 779
    default:
780 780
      return UNDEFINED;
781 781
    }
782 782
#else
783 783
    statusSwitch(cplexEnv(),stat);
784 784
    switch (stat) {
785 785
    case 0:
786 786
      return UNDEFINED; //Undefined
787 787
    case CPX_OPTIMAL://Optimal
788 788
      return OPTIMAL;
789 789
    case CPX_UNBOUNDED:
790 790
      return INFEASIBLE;
791 791
    default:
792 792
      return UNDEFINED; //Everything else comes here
793 793
      //FIXME error
794 794
    }
795 795
#endif
796 796
  }
797 797

	
798 798
  // CplexMip members
799 799

	
800 800
  CplexMip::CplexMip()
801
    : LpBase(), CplexBase(), MipSolver() {
801
    : LpBase(), MipSolver(), CplexBase() {
802 802

	
803 803
#if CPX_VERSION < 800
804 804
    CPXchgprobtype(cplexEnv(),  _prob, CPXPROB_MIP);
805 805
#else
806 806
    CPXchgprobtype(cplexEnv(),  _prob, CPXPROB_MILP);
807 807
#endif
808 808
  }
809 809

	
810 810
  CplexMip::CplexMip(const CplexEnv& env)
811
    : LpBase(), CplexBase(env), MipSolver() {
811
    : LpBase(), MipSolver(), CplexBase(env) {
812 812

	
813 813
#if CPX_VERSION < 800
814 814
    CPXchgprobtype(cplexEnv(),  _prob, CPXPROB_MIP);
815 815
#else
816 816
    CPXchgprobtype(cplexEnv(),  _prob, CPXPROB_MILP);
817 817
#endif
818 818

	
819 819
  }
820 820

	
821 821
  CplexMip::CplexMip(const CplexMip& other)
822
    : LpBase(), CplexBase(other), MipSolver() {}
822
    : LpBase(), MipSolver(), CplexBase(other) {}
823 823

	
824 824
  CplexMip::~CplexMip() {}
825 825

	
826 826
  CplexMip* CplexMip::newSolver() const { return new CplexMip; }
827 827
  CplexMip* CplexMip::cloneSolver() const {return new CplexMip(*this); }
828 828

	
829 829
  const char* CplexMip::_solverName() const { return "CplexMip"; }
830 830

	
831 831
  void CplexMip::_setColType(int i, CplexMip::ColTypes col_type) {
832 832

	
833 833
    // Note If a variable is to be changed to binary, a call to CPXchgbds
834 834
    // should also be made to change the bounds to 0 and 1.
835 835

	
836 836
    switch (col_type){
837 837
    case INTEGER: {
838 838
      const char t = 'I';
839 839
      CPXchgctype (cplexEnv(), _prob, 1, &i, &t);
840 840
    } break;
841 841
    case REAL: {
842 842
      const char t = 'C';
843 843
      CPXchgctype (cplexEnv(), _prob, 1, &i, &t);
844 844
    } break;
845 845
    default:
846 846
      break;
847 847
    }
848 848
  }
849 849

	
850 850
  CplexMip::ColTypes CplexMip::_getColType(int i) const {
851 851
    char t;
852 852
    CPXgetctype (cplexEnv(), _prob, &t, i, i);
853 853
    switch (t) {
854 854
    case 'I':
855 855
      return INTEGER;
856 856
    case 'C':
857 857
      return REAL;
858 858
    default:
859 859
      LEMON_ASSERT(false, "Invalid column type");
860 860
      return ColTypes();
861 861
    }
862 862

	
863 863
  }
864 864

	
865 865
  CplexMip::SolveExitStatus CplexMip::_solve() {
866 866
    int status;
867 867
    status = CPXmipopt (cplexEnv(), _prob);
868 868
    if (status==0)
869 869
      return SOLVED;
870 870
    else
871 871
      return UNSOLVED;
872 872

	
873 873
  }
874 874

	
875 875

	
876 876
  CplexMip::ProblemType CplexMip::_getType() const {
877 877

	
878 878
    int stat = CPXgetstat(cplexEnv(), _prob);
879 879

	
880 880
    //Fortunately, MIP statuses did not change for cplex 8.0
881 881
    switch (stat) {
882 882
    case CPXMIP_OPTIMAL:
883 883
      // Optimal integer solution has been found.
884 884
    case CPXMIP_OPTIMAL_TOL:
885 885
      // Optimal soluton with the tolerance defined by epgap or epagap has
886 886
      // been found.
887 887
      return OPTIMAL;
888 888
      //This also exists in later issues
889 889
      //    case CPXMIP_UNBOUNDED:
890 890
      //return UNBOUNDED;
891 891
      case CPXMIP_INFEASIBLE:
892 892
        return INFEASIBLE;
893 893
    default:
894 894
      return UNDEFINED;
895 895
    }
896 896
    //Unboundedness not treated well: the following is from cplex 9.0 doc
897 897
    // About Unboundedness
898 898

	
899 899
    // The treatment of models that are unbounded involves a few
900 900
    // subtleties. Specifically, a declaration of unboundedness means that
901 901
    // ILOG CPLEX has determined that the model has an unbounded
902 902
    // ray. Given any feasible solution x with objective z, a multiple of
903 903
    // the unbounded ray can be added to x to give a feasible solution
904 904
    // with objective z-1 (or z+1 for maximization models). Thus, if a
905 905
    // feasible solution exists, then the optimal objective is
906 906
    // unbounded. Note that ILOG CPLEX has not necessarily concluded that
907 907
    // a feasible solution exists. Users can call the routine CPXsolninfo
908 908
    // to determine whether ILOG CPLEX has also concluded that the model
909 909
    // has a feasible solution.
910 910
  }
911 911

	
912 912
  CplexMip::Value CplexMip::_getSol(int i) const {
913 913
    Value x;
914 914
    CPXgetmipx(cplexEnv(), _prob, &x, i, i);
915 915
    return x;
916 916
  }
917 917

	
918 918
  CplexMip::Value CplexMip::_getSolValue() const {
919 919
    Value objval;
920 920
    CPXgetmipobjval(cplexEnv(), _prob, &objval);
921 921
    return objval;
922 922
  }
923 923

	
924 924
} //namespace lemon
925 925

	
Ignore white space 384 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
 * Copyright (C) 2003-2008
5
 * Copyright (C) 2003-2009
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_CPLEX_H
20 20
#define LEMON_CPLEX_H
21 21

	
22 22
///\file
23 23
///\brief Header of the LEMON-CPLEX lp solver interface.
24 24

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

	
27 27
struct cpxenv;
28 28
struct cpxlp;
29 29

	
30 30
namespace lemon {
31 31

	
32 32
  /// \brief Reference counted wrapper around cpxenv pointer
33 33
  ///
34 34
  /// The cplex uses environment object which is responsible for
35 35
  /// checking the proper license usage. This class provides a simple
36 36
  /// interface for share the environment object between different
37 37
  /// problems.
38 38
  class CplexEnv {
39 39
    friend class CplexBase;
40 40
  private:
41 41
    cpxenv* _env;
42 42
    mutable int* _cnt;
43 43

	
44 44
  public:
45 45

	
46 46
    /// \brief This exception is thrown when the license check is not
47 47
    /// sufficient
48 48
    class LicenseError : public Exception {
49 49
      friend class CplexEnv;
50 50
    private:
51 51

	
52 52
      LicenseError(int status);
53 53
      char _message[510];
54 54

	
55 55
    public:
56 56

	
57 57
      /// The short error message
58 58
      virtual const char* what() const throw() {
59 59
        return _message;
60 60
      }
61 61
    };
62 62

	
63 63
    /// Constructor
64 64
    CplexEnv();
65 65
    /// Shallow copy constructor
66 66
    CplexEnv(const CplexEnv&);
67 67
    /// Shallow assignement
68 68
    CplexEnv& operator=(const CplexEnv&);
69 69
    /// Destructor
70 70
    virtual ~CplexEnv();
71 71

	
72 72
  protected:
73 73

	
74 74
    cpxenv* cplexEnv() { return _env; }
75 75
    const cpxenv* cplexEnv() const { return _env; }
76 76
  };
77 77

	
78 78
  /// \brief Base interface for the CPLEX LP and MIP solver
79 79
  ///
80 80
  /// This class implements the common interface of the CPLEX LP and
81
  /// MIP solvers.  
81
  /// MIP solvers.
82 82
  /// \ingroup lp_group
83 83
  class CplexBase : virtual public LpBase {
84 84
  protected:
85 85

	
86 86
    CplexEnv _env;
87 87
    cpxlp* _prob;
88 88

	
89 89
    CplexBase();
90 90
    CplexBase(const CplexEnv&);
91 91
    CplexBase(const CplexBase &);
92 92
    virtual ~CplexBase();
93 93

	
94 94
    virtual int _addCol();
95 95
    virtual int _addRow();
96 96

	
97 97
    virtual void _eraseCol(int i);
98 98
    virtual void _eraseRow(int i);
99 99

	
100 100
    virtual void _eraseColId(int i);
101 101
    virtual void _eraseRowId(int i);
102 102

	
103 103
    virtual void _getColName(int col, std::string& name) const;
104 104
    virtual void _setColName(int col, const std::string& name);
105 105
    virtual int _colByName(const std::string& name) const;
106 106

	
107 107
    virtual void _getRowName(int row, std::string& name) const;
108 108
    virtual void _setRowName(int row, const std::string& name);
109 109
    virtual int _rowByName(const std::string& name) const;
110 110

	
111 111
    virtual void _setRowCoeffs(int i, ExprIterator b, ExprIterator e);
112 112
    virtual void _getRowCoeffs(int i, InsertIterator b) const;
113 113

	
114 114
    virtual void _setColCoeffs(int i, ExprIterator b, ExprIterator e);
115 115
    virtual void _getColCoeffs(int i, InsertIterator b) const;
116 116

	
117 117
    virtual void _setCoeff(int row, int col, Value value);
118 118
    virtual Value _getCoeff(int row, int col) const;
119 119

	
120 120
    virtual void _setColLowerBound(int i, Value value);
121 121
    virtual Value _getColLowerBound(int i) const;
122 122

	
123 123
    virtual void _setColUpperBound(int i, Value value);
124 124
    virtual Value _getColUpperBound(int i) const;
125 125

	
126 126
  private:
127 127
    void _set_row_bounds(int i, Value lb, Value ub);
128 128
  protected:
129 129

	
130 130
    virtual void _setRowLowerBound(int i, Value value);
131 131
    virtual Value _getRowLowerBound(int i) const;
132 132

	
133 133
    virtual void _setRowUpperBound(int i, Value value);
134 134
    virtual Value _getRowUpperBound(int i) const;
135 135

	
136 136
    virtual void _setObjCoeffs(ExprIterator b, ExprIterator e);
137 137
    virtual void _getObjCoeffs(InsertIterator b) const;
138 138

	
139 139
    virtual void _setObjCoeff(int i, Value obj_coef);
140 140
    virtual Value _getObjCoeff(int i) const;
141 141

	
142 142
    virtual void _setSense(Sense sense);
143 143
    virtual Sense _getSense() const;
144 144

	
145 145
    virtual void _clear();
146 146

	
147 147
  public:
148 148

	
149 149
    /// Returns the used \c CplexEnv instance
150 150
    const CplexEnv& env() const { return _env; }
151 151
    ///
152 152
    const cpxenv* cplexEnv() const { return _env.cplexEnv(); }
153 153

	
154 154
    cpxlp* cplexLp() { return _prob; }
155 155
    const cpxlp* cplexLp() const { return _prob; }
156 156

	
157 157
  };
158 158

	
159 159
  /// \brief Interface for the CPLEX LP solver
160 160
  ///
161 161
  /// This class implements an interface for the CPLEX LP solver.
162 162
  ///\ingroup lp_group
163 163
  class CplexLp : public LpSolver, public CplexBase {
164 164
  public:
165 165
    /// \e
166 166
    CplexLp();
167 167
    /// \e
168 168
    CplexLp(const CplexEnv&);
169 169
    /// \e
170 170
    CplexLp(const CplexLp&);
171 171
    /// \e
172 172
    virtual ~CplexLp();
173 173

	
174 174
    /// \e
175 175
    virtual CplexLp* cloneSolver() const;
176 176
    /// \e
177 177
    virtual CplexLp* newSolver() const;
178 178

	
179 179
  private:
180 180

	
181 181
    // these values cannot retrieved element by element
182 182
    mutable std::vector<int> _col_status;
183 183
    mutable std::vector<int> _row_status;
184 184

	
185 185
    mutable std::vector<Value> _primal_ray;
186 186
    mutable std::vector<Value> _dual_ray;
187 187

	
188 188
    void _clear_temporals();
189 189

	
190 190
    SolveExitStatus convertStatus(int status);
191 191

	
192 192
  protected:
193 193

	
194 194
    virtual const char* _solverName() const;
195 195

	
196 196
    virtual SolveExitStatus _solve();
197 197
    virtual Value _getPrimal(int i) const;
198 198
    virtual Value _getDual(int i) const;
199 199
    virtual Value _getPrimalValue() const;
200 200

	
201 201
    virtual VarStatus _getColStatus(int i) const;
202 202
    virtual VarStatus _getRowStatus(int i) const;
203 203

	
204 204
    virtual Value _getPrimalRay(int i) const;
205 205
    virtual Value _getDualRay(int i) const;
206 206

	
207 207
    virtual ProblemType _getPrimalType() const;
208 208
    virtual ProblemType _getDualType() const;
209 209

	
210 210
  public:
211 211

	
212 212
    /// Solve with primal simplex method
213 213
    SolveExitStatus solvePrimal();
214 214

	
215 215
    /// Solve with dual simplex method
216 216
    SolveExitStatus solveDual();
217 217

	
218 218
    /// Solve with barrier method
219 219
    SolveExitStatus solveBarrier();
220 220

	
221 221
  };
222 222

	
223 223
  /// \brief Interface for the CPLEX MIP solver
224 224
  ///
225 225
  /// This class implements an interface for the CPLEX MIP solver.
226 226
  ///\ingroup lp_group
227 227
  class CplexMip : public MipSolver, public CplexBase {
228 228
  public:
229 229
    /// \e
230 230
    CplexMip();
231 231
    /// \e
232 232
    CplexMip(const CplexEnv&);
233 233
    /// \e
234 234
    CplexMip(const CplexMip&);
235 235
    /// \e
236 236
    virtual ~CplexMip();
237 237

	
238
    /// \e
239
    virtual CplexMip* cloneSolver() const;
240
    /// \e
241
    virtual CplexMip* newSolver() const;
242

	
238 243
  protected:
239 244

	
240
    virtual CplexMip* _cloneSolver() const;
241
    virtual CplexMip* _newSolver() const;
242 245

	
243 246
    virtual const char* _solverName() const;
244 247

	
245 248
    virtual ColTypes _getColType(int col) const;
246 249
    virtual void _setColType(int col, ColTypes col_type);
247 250

	
248 251
    virtual SolveExitStatus _solve();
249 252
    virtual ProblemType _getType() const;
250 253
    virtual Value _getSol(int i) const;
251 254
    virtual Value _getSolValue() const;
252 255

	
253 256
  };
254 257

	
255 258
} //END OF NAMESPACE LEMON
256 259

	
257 260
#endif //LEMON_CPLEX_H
258 261

	
Ignore white space 384 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
 * Copyright (C) 2003-2008
5
 * Copyright (C) 2003-2009
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
///\file
20 20
///\brief Implementation of the LEMON GLPK LP and MIP solver interface.
21 21

	
22 22
#include <lemon/glpk.h>
23 23
#include <glpk.h>
24 24

	
25 25
#include <lemon/assert.h>
26 26

	
27 27
namespace lemon {
28 28

	
29 29
  // GlpkBase members
30 30

	
31 31
  GlpkBase::GlpkBase() : LpBase() {
32 32
    lp = glp_create_prob();
33 33
    glp_create_index(lp);
34 34
  }
35 35

	
36 36
  GlpkBase::GlpkBase(const GlpkBase &other) : LpBase() {
37 37
    lp = glp_create_prob();
38 38
    glp_copy_prob(lp, other.lp, GLP_ON);
39 39
    glp_create_index(lp);
40 40
    rows = other.rows;
41 41
    cols = other.cols;
42 42
  }
43 43

	
44 44
  GlpkBase::~GlpkBase() {
45 45
    glp_delete_prob(lp);
46 46
  }
47 47

	
48 48
  int GlpkBase::_addCol() {
49 49
    int i = glp_add_cols(lp, 1);
50 50
    glp_set_col_bnds(lp, i, GLP_FR, 0.0, 0.0);
51 51
    return i;
52 52
  }
53 53

	
54 54
  int GlpkBase::_addRow() {
55 55
    int i = glp_add_rows(lp, 1);
56 56
    glp_set_row_bnds(lp, i, GLP_FR, 0.0, 0.0);
57 57
    return i;
58 58
  }
59 59

	
60 60
  void GlpkBase::_eraseCol(int i) {
61 61
    int ca[2];
62 62
    ca[1] = i;
63 63
    glp_del_cols(lp, 1, ca);
64 64
  }
65 65

	
66 66
  void GlpkBase::_eraseRow(int i) {
67 67
    int ra[2];
68 68
    ra[1] = i;
69 69
    glp_del_rows(lp, 1, ra);
70 70
  }
71 71

	
72 72
  void GlpkBase::_eraseColId(int i) {
73 73
    cols.eraseIndex(i);
74 74
    cols.shiftIndices(i);
75 75
  }
76 76

	
77 77
  void GlpkBase::_eraseRowId(int i) {
78 78
    rows.eraseIndex(i);
79 79
    rows.shiftIndices(i);
80 80
  }
81 81

	
82 82
  void GlpkBase::_getColName(int c, std::string& name) const {
83 83
    const char *str = glp_get_col_name(lp, c);
84 84
    if (str) name = str;
85 85
    else name.clear();
86 86
  }
87 87

	
88 88
  void GlpkBase::_setColName(int c, const std::string & name) {
89 89
    glp_set_col_name(lp, c, const_cast<char*>(name.c_str()));
90 90

	
91 91
  }
92 92

	
93 93
  int GlpkBase::_colByName(const std::string& name) const {
94 94
    int k = glp_find_col(lp, const_cast<char*>(name.c_str()));
95 95
    return k > 0 ? k : -1;
96 96
  }
97 97

	
98 98
  void GlpkBase::_getRowName(int r, std::string& name) const {
99 99
    const char *str = glp_get_row_name(lp, r);
100 100
    if (str) name = str;
101 101
    else name.clear();
102 102
  }
103 103

	
104 104
  void GlpkBase::_setRowName(int r, const std::string & name) {
105 105
    glp_set_row_name(lp, r, const_cast<char*>(name.c_str()));
106 106

	
107 107
  }
108 108

	
109 109
  int GlpkBase::_rowByName(const std::string& name) const {
110 110
    int k = glp_find_row(lp, const_cast<char*>(name.c_str()));
111 111
    return k > 0 ? k : -1;
112 112
  }
113 113

	
114 114
  void GlpkBase::_setRowCoeffs(int i, ExprIterator b, ExprIterator e) {
115 115
    std::vector<int> indexes;
116 116
    std::vector<Value> values;
117 117

	
118 118
    indexes.push_back(0);
119 119
    values.push_back(0);
120 120

	
121 121
    for(ExprIterator it = b; it != e; ++it) {
122 122
      indexes.push_back(it->first);
123 123
      values.push_back(it->second);
124 124
    }
125 125

	
126 126
    glp_set_mat_row(lp, i, values.size() - 1,
127 127
                    &indexes.front(), &values.front());
128 128
  }
129 129

	
130 130
  void GlpkBase::_getRowCoeffs(int ix, InsertIterator b) const {
131 131
    int length = glp_get_mat_row(lp, ix, 0, 0);
132 132

	
133 133
    std::vector<int> indexes(length + 1);
134 134
    std::vector<Value> values(length + 1);
135 135

	
136 136
    glp_get_mat_row(lp, ix, &indexes.front(), &values.front());
137 137

	
138 138
    for (int i = 1; i <= length; ++i) {
139 139
      *b = std::make_pair(indexes[i], values[i]);
140 140
      ++b;
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) {
... ...
@@ -342,617 +342,617 @@
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 525
  void GlpkBase::freeEnv() {
526 526
    glp_free_env();
527 527
  }
528 528

	
529 529
  GlpkBase::FreeEnvHelper GlpkBase::freeEnvHelper;
530 530

	
531 531
  // GlpkLp members
532 532

	
533 533
  GlpkLp::GlpkLp()
534
    : LpBase(), GlpkBase(), LpSolver() {
534
    : LpBase(), LpSolver(), GlpkBase() {
535 535
    messageLevel(MESSAGE_NO_OUTPUT);
536 536
  }
537 537

	
538 538
  GlpkLp::GlpkLp(const GlpkLp& other)
539
    : LpBase(other), GlpkBase(other), LpSolver(other) {
539
    : LpBase(other), LpSolver(other), GlpkBase(other) {
540 540
    messageLevel(MESSAGE_NO_OUTPUT);
541 541
  }
542 542

	
543 543
  GlpkLp* GlpkLp::newSolver() const { return new GlpkLp; }
544 544
  GlpkLp* GlpkLp::cloneSolver() const { return new GlpkLp(*this); }
545 545

	
546 546
  const char* GlpkLp::_solverName() const { return "GlpkLp"; }
547 547

	
548 548
  void GlpkLp::_clear_temporals() {
549 549
    _primal_ray.clear();
550 550
    _dual_ray.clear();
551 551
  }
552 552

	
553 553
  GlpkLp::SolveExitStatus GlpkLp::_solve() {
554 554
    return solvePrimal();
555 555
  }
556 556

	
557 557
  GlpkLp::SolveExitStatus GlpkLp::solvePrimal() {
558 558
    _clear_temporals();
559 559

	
560 560
    glp_smcp smcp;
561 561
    glp_init_smcp(&smcp);
562 562

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

	
578 578
    if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
579 579
    return SOLVED;
580 580
  }
581 581

	
582 582
  GlpkLp::SolveExitStatus GlpkLp::solveDual() {
583 583
    _clear_temporals();
584 584

	
585 585
    glp_smcp smcp;
586 586
    glp_init_smcp(&smcp);
587 587

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

	
604 604
    if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
605 605
    return SOLVED;
606 606
  }
607 607

	
608 608
  GlpkLp::Value GlpkLp::_getPrimal(int i) const {
609 609
    return glp_get_col_prim(lp, i);
610 610
  }
611 611

	
612 612
  GlpkLp::Value GlpkLp::_getDual(int i) const {
613 613
    return glp_get_row_dual(lp, i);
614 614
  }
615 615

	
616 616
  GlpkLp::Value GlpkLp::_getPrimalValue() const {
617 617
    return glp_get_obj_val(lp);
618 618
  }
619 619

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

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

	
656 656
  GlpkLp::Value GlpkLp::_getPrimalRay(int i) const {
657 657
    if (_primal_ray.empty()) {
658 658
      int row_num = glp_get_num_rows(lp);
659 659
      int col_num = glp_get_num_cols(lp);
660 660

	
661 661
      _primal_ray.resize(col_num + 1, 0.0);
662 662

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

	
670 670
        bool negate = glp_get_obj_dir(lp) == GLP_MAX;
671 671

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

	
683 683
        std::vector<int> ray_indexes(row_num + 1);
684 684
        std::vector<Value> ray_values(row_num + 1);
685 685
        int ray_length = glp_eval_tab_col(lp, index, &ray_indexes.front(),
686 686
                                          &ray_values.front());
687 687

	
688 688
        for (int i = 1; i <= ray_length; ++i) {
689 689
          if (ray_indexes[i] > row_num) {
690 690
            _primal_ray[ray_indexes[i] - row_num] = ray_values[i];
691 691
          }
692 692
        }
693 693

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

	
708 708
  GlpkLp::Value GlpkLp::_getDualRay(int i) const {
709 709
    if (_dual_ray.empty()) {
710 710
      int row_num = glp_get_num_rows(lp);
711 711

	
712 712
      _dual_ray.resize(row_num + 1, 0.0);
713 713

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

	
720 720
                     "Wrong dual ray");
721 721

	
722 722
        int idx;
723 723
        bool negate = false;
724 724

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

	
738 738
        _dual_ray[idx] = negate ?  - 1.0 : 1.0;
739 739

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

	
770 770
        glp_btran(lp, &_dual_ray.front());
771 771

	
772 772
        for (int i = 1; i <= row_num; ++i) {
773 773
          _dual_ray[i] /= glp_get_rii(lp, i);
774 774
        }
775 775
      }
776 776
    }
777 777
    return _dual_ray[i];
778 778
  }
779 779

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

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

	
822 822
  void GlpkLp::presolver(bool b) {
823 823
    lpx_set_int_parm(lp, LPX_K_PRESOL, b ? 1 : 0);
824 824
  }
825 825

	
826 826
  void GlpkLp::messageLevel(MessageLevel m) {
827 827
    _message_level = m;
828 828
  }
829 829

	
830 830
  // GlpkMip members
831 831

	
832 832
  GlpkMip::GlpkMip()
833
    : LpBase(), GlpkBase(), MipSolver() {
833
    : LpBase(), MipSolver(), GlpkBase() {
834 834
    messageLevel(MESSAGE_NO_OUTPUT);
835 835
  }
836 836

	
837 837
  GlpkMip::GlpkMip(const GlpkMip& other)
838
    : LpBase(), GlpkBase(other), MipSolver() {
838
    : LpBase(), MipSolver(), GlpkBase(other) {
839 839
    messageLevel(MESSAGE_NO_OUTPUT);
840 840
  }
841 841

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

	
853 853
  GlpkMip::ColTypes GlpkMip::_getColType(int i) const {
854 854
    switch (glp_get_col_kind(lp, i)) {
855 855
    case GLP_IV:
856 856
    case GLP_BV:
857 857
      return INTEGER;
858 858
    default:
859 859
      return REAL;
860 860
    }
861 861

	
862 862
  }
863 863

	
864 864
  GlpkMip::SolveExitStatus GlpkMip::_solve() {
865 865
    glp_smcp smcp;
866 866
    glp_init_smcp(&smcp);
867 867

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

	
884 884
    if (glp_simplex(lp, &smcp) != 0) return UNSOLVED;
885 885
    if (glp_get_status(lp) != GLP_OPT) return SOLVED;
886 886

	
887 887
    glp_iocp iocp;
888 888
    glp_init_iocp(&iocp);
889 889

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

	
905 905
    if (glp_intopt(lp, &iocp) != 0) return UNSOLVED;
906 906
    return SOLVED;
907 907
  }
908 908

	
909 909

	
910 910
  GlpkMip::ProblemType GlpkMip::_getType() const {
911 911
    switch (glp_get_status(lp)) {
912 912
    case GLP_OPT:
913 913
      switch (glp_mip_status(lp)) {
914 914
      case GLP_UNDEF:
915 915
        return UNDEFINED;
916 916
      case GLP_NOFEAS:
917 917
        return INFEASIBLE;
918 918
      case GLP_FEAS:
919 919
        return FEASIBLE;
920 920
      case GLP_OPT:
921 921
        return OPTIMAL;
922 922
      default:
923 923
        LEMON_ASSERT(false, "Wrong problem type.");
924 924
        return GlpkMip::ProblemType();
925 925
      }
926 926
    case GLP_NOFEAS:
927 927
      return INFEASIBLE;
928 928
    case GLP_INFEAS:
929 929
    case GLP_FEAS:
930 930
      if (glp_get_dual_stat(lp) == GLP_NOFEAS) {
931 931
        return UNBOUNDED;
932 932
      } else {
933 933
        return UNDEFINED;
934 934
      }
935 935
    default:
936 936
      LEMON_ASSERT(false, "Wrong problem type.");
937 937
      return GlpkMip::ProblemType();
938 938
    }
939 939
  }
940 940

	
941 941
  GlpkMip::Value GlpkMip::_getSol(int i) const {
942 942
    return glp_mip_col_val(lp, i);
943 943
  }
944 944

	
945 945
  GlpkMip::Value GlpkMip::_getSolValue() const {
946 946
    return glp_mip_obj_val(lp);
947 947
  }
948 948

	
949 949
  GlpkMip* GlpkMip::newSolver() const { return new GlpkMip; }
950 950
  GlpkMip* GlpkMip::cloneSolver() const {return new GlpkMip(*this); }
951 951

	
952 952
  const char* GlpkMip::_solverName() const { return "GlpkMip"; }
953 953

	
954 954
  void GlpkMip::messageLevel(MessageLevel m) {
955 955
    _message_level = m;
956 956
  }
957 957

	
958 958
} //END OF NAMESPACE LEMON
Ignore white space 384 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
 * Copyright (C) 2003-2008
5
 * Copyright (C) 2003-2009
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
    //Test for clone/new
201 201
    LP* lpnew = lp.newSolver();
202 202
    LP* lpclone = lp.cloneSolver();
203 203
    delete lpnew;
204 204
    delete lpclone;
205 205

	
206 206
  }
207 207

	
208 208
  {
209 209
    LP::DualExpr e,f,g;
210 210
    LP::Row p1 = INVALID, p2 = INVALID, p3 = INVALID,
211 211
      p4 = INVALID, p5 = INVALID;
212 212

	
213 213
    e[p1]=2;
214 214
    e[p1]+=2;
215 215
    e[p1]-=2;
216 216

	
217 217
    e=p1;
218 218
    e=f;
219 219

	
220 220
    e+=p1;
221 221
    e+=f;
222 222

	
223 223
    e-=p1;
224 224
    e-=f;
225 225

	
226 226
    e*=2;
227 227
    e*=2.2;
228 228
    e/=2;
229 229
    e/=2.2;
230 230

	
231 231
    e=((p1+p2)+(p1-p2)+
232 232
       (p1+f)+(f+p1)+(f+g)+
233 233
       (p1-f)+(f-p1)+(f-g)+
234 234
       2.2*f+f*2.2+f/2.2+
235 235
       2*f+f*2+f/2+
236 236
       2.2*p1+p1*2.2+p1/2.2+
237 237
       2*p1+p1*2+p1/2
238 238
       );
239 239
  }
240 240

	
241 241
}
242 242

	
243 243
void solveAndCheck(LpSolver& lp, LpSolver::ProblemType stat,
244 244
                   double exp_opt) {
245 245
  using std::string;
246 246
  lp.solve();
247 247

	
248 248
  std::ostringstream buf;
249 249
  buf << "PrimalType should be: " << int(stat) << int(lp.primalType());
250 250

	
251 251
  check(lp.primalType()==stat, buf.str());
252 252

	
253 253
  if (stat ==  LpSolver::OPTIMAL) {
254 254
    std::ostringstream sbuf;
255 255
    sbuf << "Wrong optimal value (" << lp.primal() <<") with "
256 256
         << lp.solverName() <<"\n     the right optimum is " << exp_opt;
257 257
    check(std::abs(lp.primal()-exp_opt) < 1e-3, sbuf.str());
258 258
  }
259 259
}
260 260

	
261 261
void aTest(LpSolver & lp)
262 262
{
263 263
  typedef LpSolver LP;
264 264

	
265 265
 //The following example is very simple
266 266

	
267 267
  typedef LpSolver::Row Row;
268 268
  typedef LpSolver::Col Col;
269 269

	
270 270

	
271 271
  Col x1 = lp.addCol();
272 272
  Col x2 = lp.addCol();
273 273

	
274 274

	
275 275
  //Constraints
276 276
  Row upright=lp.addRow(x1+2*x2 <=1);
277 277
  lp.addRow(x1+x2 >=-1);
278 278
  lp.addRow(x1-x2 <=1);
279 279
  lp.addRow(x1-x2 >=-1);
280 280
  //Nonnegativity of the variables
281 281
  lp.colLowerBound(x1, 0);
282 282
  lp.colLowerBound(x2, 0);
283 283
  //Objective function
284 284
  lp.obj(x1+x2);
285 285

	
286 286
  lp.sense(lp.MAX);
287 287

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

	
304 304
  lp.row(upright, x1+x2 <=1);
305 305
  e = lp.row(upright);
306 306
  check(e[x1] == 1, "The first coefficient should 1.");
307 307
  check(e[x2] == 1, "The second coefficient should 1.");
308 308

	
309 309
  LpSolver::DualExpr de = lp.col(x1);
310 310
  check(  de[upright] == 1, "The first coefficient should 1.");
311 311

	
312 312
  LpSolver* clp = lp.cloneSolver();
313 313

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

	
324 324
  check(lp.rowLowerBound(upright)==-LpSolver::INF,
325 325
        "The lower bound for the first row should be -infty.");
326 326
  check(lp.rowUpperBound(upright)==1,
327 327
        "The upper bound for the first row should be 1.");
328 328
  e = clp->row(upright);
329 329
  check(e[x1] == 1, "The first coefficient should 1.");
330 330
  check(e[x2] == 1, "The second coefficient should 1.");
331 331

	
332 332
  de = clp->col(x1);
333 333
  check(de[upright] == 1, "The first coefficient should 1.");
334 334

	
335 335
  delete clp;
336 336

	
337 337
  //Maximization of x1+x2
338 338
  //over the triangle with vertices (0,0) (0,1) (1,0)
339 339
  double expected_opt=1;
340 340
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
341 341

	
342 342
  //Minimization
343 343
  lp.sense(lp.MIN);
344 344
  expected_opt=0;
345 345
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
346 346

	
347 347
  //Vertex (-1,0) instead of (0,0)
348 348
  lp.colLowerBound(x1, -LpSolver::INF);
349 349
  expected_opt=-1;
350 350
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
351 351

	
352 352
  //Erase one constraint and return to maximization
353 353
  lp.erase(upright);
354 354
  lp.sense(lp.MAX);
355 355
  expected_opt=LpSolver::INF;
356 356
  solveAndCheck(lp, LpSolver::UNBOUNDED, expected_opt);
357 357

	
358 358
  //Infeasibilty
359 359
  lp.addRow(x1+x2 <=-2);
360 360
  solveAndCheck(lp, LpSolver::INFEASIBLE, expected_opt);
361 361

	
362 362
}
363 363

	
364 364
template<class LP>
365 365
void cloneTest()
366 366
{
367 367
  //Test for clone/new
368
  
368

	
369 369
  LP* lp = new LP();
370 370
  LP* lpnew = lp->newSolver();
371 371
  LP* lpclone = lp->cloneSolver();
372 372
  delete lp;
373 373
  delete lpnew;
374 374
  delete lpclone;
375 375
}
376 376

	
377 377
int main()
378 378
{
379 379
  LpSkeleton lp_skel;
380 380
  lpTest(lp_skel);
381 381

	
382 382
#ifdef HAVE_GLPK
383 383
  {
384 384
    GlpkLp lp_glpk1,lp_glpk2;
385 385
    lpTest(lp_glpk1);
386 386
    aTest(lp_glpk2);
387 387
    cloneTest<GlpkLp>();
388 388
  }
389 389
#endif
390 390

	
391 391
#ifdef HAVE_CPLEX
392 392
  try {
393 393
    CplexLp lp_cplex1,lp_cplex2;
394 394
    lpTest(lp_cplex1);
395 395
    aTest(lp_cplex2);
396
    cloneTest<CplexLp>();
396 397
  } catch (CplexEnv::LicenseError& error) {
397 398
#ifdef LEMON_FORCE_CPLEX_CHECK
398 399
    check(false, error.what());
399 400
#else
400 401
    std::cerr << error.what() << std::endl;
401 402
    std::cerr << "Cplex license check failed, lp check skipped" << std::endl;
402 403
#endif
403 404
  }
404
    cloneTest<CplexLp>();
405 405
#endif
406 406

	
407 407
#ifdef HAVE_SOPLEX
408 408
  {
409 409
    SoplexLp lp_soplex1,lp_soplex2;
410 410
    lpTest(lp_soplex1);
411 411
    aTest(lp_soplex2);
412 412
    cloneTest<SoplexLp>();
413 413
  }
414 414
#endif
415 415

	
416 416
#ifdef HAVE_CLP
417 417
  {
418 418
    ClpLp lp_clp1,lp_clp2;
419 419
    lpTest(lp_clp1);
420 420
    aTest(lp_clp2);
421 421
    cloneTest<ClpLp>();
422 422
  }
423 423
#endif
424 424

	
425 425
  return 0;
426 426
}
Ignore white space 384 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
 * Copyright (C) 2003-2008
5
 * Copyright (C) 2003-2009
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
template<class MIP>
110 110
void cloneTest()
111 111
{
112
  
112

	
113 113
  MIP* mip = new MIP();
114 114
  MIP* mipnew = mip->newSolver();
115 115
  MIP* mipclone = mip->cloneSolver();
116 116
  delete mip;
117 117
  delete mipnew;
118 118
  delete mipclone;
119 119
}
120 120

	
121 121
int main()
122 122
{
123 123

	
124 124
#ifdef HAVE_GLPK
125 125
  {
126 126
    GlpkMip mip1;
127 127
    aTest(mip1);
128 128
    cloneTest<GlpkMip>();
129 129
  }
130 130
#endif
131 131

	
132 132
#ifdef HAVE_CPLEX
133 133
  try {
134 134
    CplexMip mip2;
135 135
    aTest(mip2);
136
    cloneTest<CplexMip>();
136 137
  } catch (CplexEnv::LicenseError& error) {
137 138
#ifdef LEMON_FORCE_CPLEX_CHECK
138 139
    check(false, error.what());
139 140
#else
140 141
    std::cerr << error.what() << std::endl;
141 142
    std::cerr << "Cplex license check failed, lp check skipped" << std::endl;
142 143
#endif
143 144
  }
144
  cloneTest<CplexMip>();
145 145
#endif
146 146

	
147 147
  return 0;
148 148

	
149 149
}
0 comments (0 inline)