test/lp_test.cc
author Alpar Juttner <alpar@cs.elte.hu>
Fri, 23 Mar 2018 16:09:06 +0100
branch1.3
changeset 1390 f0f15d07bf51
parent 1270 dceba191c00d
child 1337 4add05447ca0
permissions -rw-r--r--
Merge bugfix #609 to branch 1.3
deba@481
     1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
deba@481
     2
 *
deba@481
     3
 * This file is a part of LEMON, a generic C++ optimization library.
deba@481
     4
 *
alpar@1270
     5
 * Copyright (C) 2003-2013
deba@481
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
deba@481
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
deba@481
     8
 *
deba@481
     9
 * Permission to use, modify and distribute this software is granted
deba@481
    10
 * provided that this copyright notice appears in all copies. For
deba@481
    11
 * precise terms see the accompanying LICENSE file.
deba@481
    12
 *
deba@481
    13
 * This software is provided "AS IS" with no warranty of any kind,
deba@481
    14
 * express or implied, and with no claim as to its suitability for any
deba@481
    15
 * purpose.
deba@481
    16
 *
deba@481
    17
 */
deba@481
    18
deba@481
    19
#include <sstream>
deba@481
    20
#include <lemon/lp_skeleton.h>
deba@481
    21
#include "test_tools.h"
deba@481
    22
#include <lemon/tolerance.h>
deba@481
    23
deba@481
    24
#include <lemon/config.h>
deba@481
    25
ladanyi@674
    26
#ifdef LEMON_HAVE_GLPK
alpar@484
    27
#include <lemon/glpk.h>
deba@481
    28
#endif
deba@481
    29
ladanyi@674
    30
#ifdef LEMON_HAVE_CPLEX
alpar@484
    31
#include <lemon/cplex.h>
deba@481
    32
#endif
deba@481
    33
ladanyi@674
    34
#ifdef LEMON_HAVE_SOPLEX
alpar@484
    35
#include <lemon/soplex.h>
deba@481
    36
#endif
deba@481
    37
ladanyi@674
    38
#ifdef LEMON_HAVE_CLP
alpar@484
    39
#include <lemon/clp.h>
deba@482
    40
#endif
deba@482
    41
alpar@1300
    42
#ifdef LEMON_HAVE_LP
alpar@1300
    43
#include <lemon/lp.h>
alpar@1300
    44
#endif
deba@481
    45
using namespace lemon;
deba@481
    46
alpar@1140
    47
int countCols(LpBase & lp) {
alpar@1140
    48
  int count=0;
alpar@1140
    49
  for (LpBase::ColIt c(lp); c!=INVALID; ++c) ++count;
alpar@1140
    50
  return count;
alpar@1140
    51
}
alpar@1140
    52
alpar@1140
    53
int countRows(LpBase & lp) {
alpar@1140
    54
  int count=0;
alpar@1140
    55
  for (LpBase::RowIt r(lp); r!=INVALID; ++r) ++count;
alpar@1140
    56
  return count;
alpar@1140
    57
}
alpar@1140
    58
alpar@1140
    59
deba@482
    60
void lpTest(LpSolver& lp)
deba@481
    61
{
deba@481
    62
deba@482
    63
  typedef LpSolver LP;
deba@481
    64
alpar@1140
    65
  // Test LpBase::clear()
alpar@1140
    66
  check(countRows(lp)==0, "Wrong number of rows");
alpar@1140
    67
  check(countCols(lp)==0, "Wrong number of cols");
alpar@1140
    68
  lp.addCol(); lp.addRow(); lp.addRow();
alpar@1140
    69
  check(countRows(lp)==2, "Wrong number of rows");
alpar@1140
    70
  check(countCols(lp)==1, "Wrong number of cols");
alpar@1140
    71
  lp.clear();
alpar@1140
    72
  check(countRows(lp)==0, "Wrong number of rows");
alpar@1140
    73
  check(countCols(lp)==0, "Wrong number of cols");
alpar@1140
    74
  lp.addCol(); lp.addCol(); lp.addCol(); lp.addRow();
alpar@1140
    75
  check(countRows(lp)==1, "Wrong number of rows");
alpar@1140
    76
  check(countCols(lp)==3, "Wrong number of cols");
alpar@1140
    77
  lp.clear();
alpar@1140
    78
deba@481
    79
  std::vector<LP::Col> x(10);
deba@481
    80
  //  for(int i=0;i<10;i++) x.push_back(lp.addCol());
deba@481
    81
  lp.addColSet(x);
deba@481
    82
  lp.colLowerBound(x,1);
deba@481
    83
  lp.colUpperBound(x,1);
deba@481
    84
  lp.colBounds(x,1,2);
deba@481
    85
deba@481
    86
  std::vector<LP::Col> y(10);
deba@481
    87
  lp.addColSet(y);
deba@481
    88
deba@481
    89
  lp.colLowerBound(y,1);
deba@481
    90
  lp.colUpperBound(y,1);
deba@481
    91
  lp.colBounds(y,1,2);
deba@481
    92
deba@481
    93
  std::map<int,LP::Col> z;
deba@481
    94
deba@481
    95
  z.insert(std::make_pair(12,INVALID));
deba@481
    96
  z.insert(std::make_pair(2,INVALID));
deba@481
    97
  z.insert(std::make_pair(7,INVALID));
deba@481
    98
  z.insert(std::make_pair(5,INVALID));
deba@481
    99
deba@481
   100
  lp.addColSet(z);
deba@481
   101
deba@481
   102
  lp.colLowerBound(z,1);
deba@481
   103
  lp.colUpperBound(z,1);
deba@481
   104
  lp.colBounds(z,1,2);
deba@481
   105
deba@481
   106
  {
deba@481
   107
    LP::Expr e,f,g;
deba@481
   108
    LP::Col p1,p2,p3,p4,p5;
deba@481
   109
    LP::Constr c;
deba@481
   110
deba@481
   111
    p1=lp.addCol();
deba@481
   112
    p2=lp.addCol();
deba@481
   113
    p3=lp.addCol();
deba@481
   114
    p4=lp.addCol();
deba@481
   115
    p5=lp.addCol();
deba@481
   116
deba@481
   117
    e[p1]=2;
deba@482
   118
    *e=12;
deba@481
   119
    e[p1]+=2;
deba@482
   120
    *e+=12;
deba@481
   121
    e[p1]-=2;
deba@482
   122
    *e-=12;
deba@481
   123
deba@481
   124
    e=2;
deba@481
   125
    e=2.2;
deba@481
   126
    e=p1;
deba@481
   127
    e=f;
deba@481
   128
deba@481
   129
    e+=2;
deba@481
   130
    e+=2.2;
deba@481
   131
    e+=p1;
deba@481
   132
    e+=f;
deba@481
   133
deba@481
   134
    e-=2;
deba@481
   135
    e-=2.2;
deba@481
   136
    e-=p1;
deba@481
   137
    e-=f;
deba@481
   138
deba@481
   139
    e*=2;
deba@481
   140
    e*=2.2;
deba@481
   141
    e/=2;
deba@481
   142
    e/=2.2;
deba@481
   143
deba@481
   144
    e=((p1+p2)+(p1-p2)+(p1+12)+(12+p1)+(p1-12)+(12-p1)+
deba@481
   145
       (f+12)+(12+f)+(p1+f)+(f+p1)+(f+g)+
deba@481
   146
       (f-12)+(12-f)+(p1-f)+(f-p1)+(f-g)+
deba@481
   147
       2.2*f+f*2.2+f/2.2+
deba@481
   148
       2*f+f*2+f/2+
deba@481
   149
       2.2*p1+p1*2.2+p1/2.2+
deba@481
   150
       2*p1+p1*2+p1/2
deba@481
   151
       );
deba@481
   152
deba@481
   153
deba@481
   154
    c = (e  <= f  );
deba@481
   155
    c = (e  <= 2.2);
deba@481
   156
    c = (e  <= 2  );
deba@481
   157
    c = (e  <= p1 );
deba@481
   158
    c = (2.2<= f  );
deba@481
   159
    c = (2  <= f  );
deba@481
   160
    c = (p1 <= f  );
deba@481
   161
    c = (p1 <= p2 );
deba@481
   162
    c = (p1 <= 2.2);
deba@481
   163
    c = (p1 <= 2  );
deba@481
   164
    c = (2.2<= p2 );
deba@481
   165
    c = (2  <= p2 );
deba@481
   166
deba@481
   167
    c = (e  >= f  );
deba@481
   168
    c = (e  >= 2.2);
deba@481
   169
    c = (e  >= 2  );
deba@481
   170
    c = (e  >= p1 );
deba@481
   171
    c = (2.2>= f  );
deba@481
   172
    c = (2  >= f  );
deba@481
   173
    c = (p1 >= f  );
deba@481
   174
    c = (p1 >= p2 );
deba@481
   175
    c = (p1 >= 2.2);
deba@481
   176
    c = (p1 >= 2  );
deba@481
   177
    c = (2.2>= p2 );
deba@481
   178
    c = (2  >= p2 );
deba@481
   179
deba@481
   180
    c = (e  == f  );
deba@481
   181
    c = (e  == 2.2);
deba@481
   182
    c = (e  == 2  );
deba@481
   183
    c = (e  == p1 );
deba@481
   184
    c = (2.2== f  );
deba@481
   185
    c = (2  == f  );
deba@481
   186
    c = (p1 == f  );
deba@481
   187
    //c = (p1 == p2 );
deba@481
   188
    c = (p1 == 2.2);
deba@481
   189
    c = (p1 == 2  );
deba@481
   190
    c = (2.2== p2 );
deba@481
   191
    c = (2  == p2 );
deba@481
   192
alpar@483
   193
    c = ((2 <= e) <= 3);
alpar@483
   194
    c = ((2 <= p1) <= 3);
deba@481
   195
alpar@483
   196
    c = ((2 >= e) >= 3);
alpar@483
   197
    c = ((2 >= p1) >= 3);
deba@481
   198
retvari@1092
   199
    { //Tests for #430
retvari@1092
   200
      LP::Col v=lp.addCol();
retvari@1092
   201
      LP::Constr c = v >= -3;
retvari@1092
   202
      c = c <= 4;
retvari@1092
   203
      LP::Constr c2;
alpar@1247
   204
#if ( __GNUC__ == 4 ) && ( __GNUC_MINOR__ == 3 )
alpar@1247
   205
      c2 = ( -3 <= v ) <= 4;
alpar@1247
   206
#else
retvari@1092
   207
      c2 = -3 <= v <= 4;
alpar@1247
   208
#endif
alpar@1247
   209
retvari@1092
   210
    }
retvari@1092
   211
deba@481
   212
    e[x[3]]=2;
deba@481
   213
    e[x[3]]=4;
deba@481
   214
    e[x[3]]=1;
deba@482
   215
    *e=12;
deba@481
   216
deba@482
   217
    lp.addRow(-LP::INF,e,23);
deba@482
   218
    lp.addRow(-LP::INF,3.0*(x[1]+x[2]/2)-x[3],23);
deba@482
   219
    lp.addRow(-LP::INF,3.0*(x[1]+x[2]*2-5*x[3]+12-x[4]/3)+2*x[4]-4,23);
deba@481
   220
deba@481
   221
    lp.addRow(x[1]+x[3]<=x[5]-3);
alpar@483
   222
    lp.addRow((-7<=x[1]+x[3]-12)<=3);
deba@481
   223
    lp.addRow(x[1]<=x[5]);
deba@481
   224
deba@481
   225
    std::ostringstream buf;
deba@481
   226
deba@481
   227
deba@481
   228
    e=((p1+p2)+(p1-0.99*p2));
deba@481
   229
    //e.prettyPrint(std::cout);
deba@481
   230
    //(e<=2).prettyPrint(std::cout);
deba@481
   231
    double tolerance=0.001;
deba@481
   232
    e.simplify(tolerance);
deba@481
   233
    buf << "Coeff. of p2 should be 0.01";
deba@481
   234
    check(e[p2]>0, buf.str());
deba@481
   235
deba@481
   236
    tolerance=0.02;
deba@481
   237
    e.simplify(tolerance);
deba@481
   238
    buf << "Coeff. of p2 should be 0";
deba@482
   239
    check(const_cast<const LpSolver::Expr&>(e)[p2]==0, buf.str());
deba@481
   240
alpar@587
   241
    //Test for clone/new
alpar@587
   242
    LP* lpnew = lp.newSolver();
alpar@587
   243
    LP* lpclone = lp.cloneSolver();
alpar@587
   244
    delete lpnew;
alpar@587
   245
    delete lpclone;
deba@481
   246
deba@481
   247
  }
deba@481
   248
deba@481
   249
  {
deba@481
   250
    LP::DualExpr e,f,g;
alpar@1232
   251
    LP::Row p1 = INVALID, p2 = INVALID;
deba@481
   252
deba@481
   253
    e[p1]=2;
deba@481
   254
    e[p1]+=2;
deba@481
   255
    e[p1]-=2;
deba@481
   256
deba@481
   257
    e=p1;
deba@481
   258
    e=f;
deba@481
   259
deba@481
   260
    e+=p1;
deba@481
   261
    e+=f;
deba@481
   262
deba@481
   263
    e-=p1;
deba@481
   264
    e-=f;
deba@481
   265
deba@481
   266
    e*=2;
deba@481
   267
    e*=2.2;
deba@481
   268
    e/=2;
deba@481
   269
    e/=2.2;
deba@481
   270
deba@481
   271
    e=((p1+p2)+(p1-p2)+
deba@481
   272
       (p1+f)+(f+p1)+(f+g)+
deba@481
   273
       (p1-f)+(f-p1)+(f-g)+
deba@481
   274
       2.2*f+f*2.2+f/2.2+
deba@481
   275
       2*f+f*2+f/2+
deba@481
   276
       2.2*p1+p1*2.2+p1/2.2+
deba@481
   277
       2*p1+p1*2+p1/2
deba@481
   278
       );
deba@481
   279
  }
deba@481
   280
deba@481
   281
}
deba@481
   282
deba@482
   283
void solveAndCheck(LpSolver& lp, LpSolver::ProblemType stat,
deba@481
   284
                   double exp_opt) {
deba@481
   285
  using std::string;
deba@481
   286
  lp.solve();
deba@482
   287
deba@481
   288
  std::ostringstream buf;
deba@482
   289
  buf << "PrimalType should be: " << int(stat) << int(lp.primalType());
deba@481
   290
deba@482
   291
  check(lp.primalType()==stat, buf.str());
deba@481
   292
deba@482
   293
  if (stat ==  LpSolver::OPTIMAL) {
deba@481
   294
    std::ostringstream sbuf;
alpar@587
   295
    sbuf << "Wrong optimal value (" << lp.primal() <<") with "
alpar@587
   296
         << lp.solverName() <<"\n     the right optimum is " << exp_opt;
deba@482
   297
    check(std::abs(lp.primal()-exp_opt) < 1e-3, sbuf.str());
deba@481
   298
  }
deba@481
   299
}
deba@481
   300
deba@482
   301
void aTest(LpSolver & lp)
deba@481
   302
{
deba@482
   303
  typedef LpSolver LP;
deba@481
   304
deba@481
   305
 //The following example is very simple
deba@481
   306
deba@482
   307
  typedef LpSolver::Row Row;
deba@482
   308
  typedef LpSolver::Col Col;
deba@481
   309
deba@481
   310
deba@481
   311
  Col x1 = lp.addCol();
deba@481
   312
  Col x2 = lp.addCol();
deba@481
   313
deba@481
   314
deba@481
   315
  //Constraints
deba@482
   316
  Row upright=lp.addRow(x1+2*x2 <=1);
deba@481
   317
  lp.addRow(x1+x2 >=-1);
deba@481
   318
  lp.addRow(x1-x2 <=1);
deba@481
   319
  lp.addRow(x1-x2 >=-1);
deba@481
   320
  //Nonnegativity of the variables
deba@481
   321
  lp.colLowerBound(x1, 0);
deba@481
   322
  lp.colLowerBound(x2, 0);
deba@481
   323
  //Objective function
deba@481
   324
  lp.obj(x1+x2);
deba@481
   325
deba@482
   326
  lp.sense(lp.MAX);
deba@481
   327
deba@481
   328
  //Testing the problem retrieving routines
deba@481
   329
  check(lp.objCoeff(x1)==1,"First term should be 1 in the obj function!");
deba@482
   330
  check(lp.sense() == lp.MAX,"This is a maximization!");
deba@481
   331
  check(lp.coeff(upright,x1)==1,"The coefficient in question is 1!");
deba@482
   332
  check(lp.colLowerBound(x1)==0,
deba@482
   333
        "The lower bound for variable x1 should be 0.");
deba@482
   334
  check(lp.colUpperBound(x1)==LpSolver::INF,
deba@482
   335
        "The upper bound for variable x1 should be infty.");
deba@482
   336
  check(lp.rowLowerBound(upright) == -LpSolver::INF,
deba@482
   337
        "The lower bound for the first row should be -infty.");
deba@482
   338
  check(lp.rowUpperBound(upright)==1,
deba@482
   339
        "The upper bound for the first row should be 1.");
deba@482
   340
  LpSolver::Expr e = lp.row(upright);
deba@482
   341
  check(e[x1] == 1, "The first coefficient should 1.");
deba@482
   342
  check(e[x2] == 2, "The second coefficient should 1.");
deba@481
   343
deba@482
   344
  lp.row(upright, x1+x2 <=1);
deba@482
   345
  e = lp.row(upright);
deba@482
   346
  check(e[x1] == 1, "The first coefficient should 1.");
deba@482
   347
  check(e[x2] == 1, "The second coefficient should 1.");
deba@482
   348
deba@482
   349
  LpSolver::DualExpr de = lp.col(x1);
deba@481
   350
  check(  de[upright] == 1, "The first coefficient should 1.");
deba@481
   351
deba@482
   352
  LpSolver* clp = lp.cloneSolver();
deba@481
   353
deba@481
   354
  //Testing the problem retrieving routines
deba@481
   355
  check(clp->objCoeff(x1)==1,"First term should be 1 in the obj function!");
deba@482
   356
  check(clp->sense() == clp->MAX,"This is a maximization!");
deba@481
   357
  check(clp->coeff(upright,x1)==1,"The coefficient in question is 1!");
deba@481
   358
  //  std::cout<<lp.colLowerBound(x1)<<std::endl;
deba@482
   359
  check(clp->colLowerBound(x1)==0,
deba@482
   360
        "The lower bound for variable x1 should be 0.");
deba@482
   361
  check(clp->colUpperBound(x1)==LpSolver::INF,
deba@482
   362
        "The upper bound for variable x1 should be infty.");
deba@481
   363
deba@482
   364
  check(lp.rowLowerBound(upright)==-LpSolver::INF,
deba@482
   365
        "The lower bound for the first row should be -infty.");
deba@482
   366
  check(lp.rowUpperBound(upright)==1,
deba@482
   367
        "The upper bound for the first row should be 1.");
deba@481
   368
  e = clp->row(upright);
deba@482
   369
  check(e[x1] == 1, "The first coefficient should 1.");
deba@482
   370
  check(e[x2] == 1, "The second coefficient should 1.");
deba@481
   371
deba@481
   372
  de = clp->col(x1);
deba@482
   373
  check(de[upright] == 1, "The first coefficient should 1.");
deba@481
   374
deba@481
   375
  delete clp;
deba@481
   376
deba@481
   377
  //Maximization of x1+x2
deba@481
   378
  //over the triangle with vertices (0,0) (0,1) (1,0)
deba@481
   379
  double expected_opt=1;
deba@482
   380
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
deba@481
   381
deba@481
   382
  //Minimization
deba@482
   383
  lp.sense(lp.MIN);
deba@481
   384
  expected_opt=0;
deba@482
   385
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
deba@481
   386
deba@481
   387
  //Vertex (-1,0) instead of (0,0)
deba@482
   388
  lp.colLowerBound(x1, -LpSolver::INF);
deba@481
   389
  expected_opt=-1;
deba@482
   390
  solveAndCheck(lp, LpSolver::OPTIMAL, expected_opt);
deba@481
   391
deba@481
   392
  //Erase one constraint and return to maximization
deba@482
   393
  lp.erase(upright);
deba@482
   394
  lp.sense(lp.MAX);
deba@482
   395
  expected_opt=LpSolver::INF;
deba@482
   396
  solveAndCheck(lp, LpSolver::UNBOUNDED, expected_opt);
deba@481
   397
deba@481
   398
  //Infeasibilty
deba@481
   399
  lp.addRow(x1+x2 <=-2);
deba@482
   400
  solveAndCheck(lp, LpSolver::INFEASIBLE, expected_opt);
deba@481
   401
deba@481
   402
}
deba@481
   403
alpar@587
   404
template<class LP>
alpar@587
   405
void cloneTest()
alpar@587
   406
{
alpar@587
   407
  //Test for clone/new
deba@598
   408
alpar@587
   409
  LP* lp = new LP();
alpar@587
   410
  LP* lpnew = lp->newSolver();
alpar@587
   411
  LP* lpclone = lp->cloneSolver();
alpar@587
   412
  delete lp;
alpar@587
   413
  delete lpnew;
alpar@587
   414
  delete lpclone;
alpar@587
   415
}
alpar@587
   416
deba@481
   417
int main()
deba@481
   418
{
deba@481
   419
  LpSkeleton lp_skel;
deba@481
   420
  lpTest(lp_skel);
deba@481
   421
alpar@1300
   422
#ifdef LEMON_HAVE_LP
alpar@1300
   423
  {
alpar@1300
   424
    Lp lp,lp2;
alpar@1300
   425
    lpTest(lp);
alpar@1300
   426
    aTest(lp2);
alpar@1300
   427
    cloneTest<Lp>();
alpar@1300
   428
  }
alpar@1300
   429
#endif
alpar@1300
   430
ladanyi@674
   431
#ifdef LEMON_HAVE_GLPK
deba@482
   432
  {
alpar@485
   433
    GlpkLp lp_glpk1,lp_glpk2;
deba@482
   434
    lpTest(lp_glpk1);
deba@482
   435
    aTest(lp_glpk2);
alpar@587
   436
    cloneTest<GlpkLp>();
deba@482
   437
  }
deba@481
   438
#endif
deba@481
   439
ladanyi@674
   440
#ifdef LEMON_HAVE_CPLEX
deba@482
   441
  try {
alpar@485
   442
    CplexLp lp_cplex1,lp_cplex2;
deba@482
   443
    lpTest(lp_cplex1);
deba@482
   444
    aTest(lp_cplex2);
deba@598
   445
    cloneTest<CplexLp>();
deba@482
   446
  } catch (CplexEnv::LicenseError& error) {
deba@482
   447
    check(false, error.what());
deba@482
   448
  }
deba@481
   449
#endif
deba@481
   450
ladanyi@674
   451
#ifdef LEMON_HAVE_SOPLEX
deba@482
   452
  {
alpar@485
   453
    SoplexLp lp_soplex1,lp_soplex2;
deba@482
   454
    lpTest(lp_soplex1);
deba@482
   455
    aTest(lp_soplex2);
alpar@587
   456
    cloneTest<SoplexLp>();
deba@482
   457
  }
deba@482
   458
#endif
deba@482
   459
ladanyi@674
   460
#ifdef LEMON_HAVE_CLP
deba@482
   461
  {
alpar@485
   462
    ClpLp lp_clp1,lp_clp2;
deba@482
   463
    lpTest(lp_clp1);
deba@482
   464
    aTest(lp_clp2);
alpar@587
   465
    cloneTest<ClpLp>();
deba@482
   466
  }
deba@481
   467
#endif
deba@481
   468
deba@481
   469
  return 0;
deba@481
   470
}