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