/home/coin/SVN-release/OS-2.4.1/Bonmin/test/InterfaceTest.cpp

Go to the documentation of this file.
00001 // (C) Copyright Carnegie Mellon University 2005, 2007
00002 // All Rights Reserved.
00003 // This code is published under the Common Public License.
00004 //
00005 // Authors :
00006 // P. Bonami, Carnegie Mellon University
00007 //
00008 // Date :  07/01/2005
00009 
00010 #include "BonminConfig.h"
00011 #include "BonAmplInterface.hpp"
00012 #include "OsiClpSolverInterface.hpp"
00013 #include "BonTMINLP.hpp"
00014 #include "IpIpoptApplication.hpp"
00015 #ifdef COIN_HAS_ASL
00016 #include "BonAmplTMINLP.hpp"
00017 #include "BonAmplSetup.hpp"
00018 #endif
00019 
00020 #include "BonIpoptSolver.hpp"
00021 #include "BonminConfig.h"
00022 
00023 #ifdef COIN_HAS_FILTERSQP
00024 #include "BonFilterSolver.hpp"
00025 #endif
00026 
00027 #include "CoinError.hpp"
00028 
00029 #include <string>
00030 #include <cmath>
00031 using namespace Bonmin;
00032 
00033 void MyAssertFunc(bool c, const std::string &s, const std::string&  file, unsigned int line){
00034    if(c != true){
00035       fprintf(stderr, "Failed MyAssertion: %s in %s line %i.\n", s.c_str(), file.c_str(), line);
00036       exit(1);
00037    }
00038 }
00039 
00040 void DblEqAssertFunc(const double& a, const std::string &a_s, const double&b, const std::string& b_s,
00041                   const std::string&  file, unsigned int line){
00042    CoinRelFltEq eq;
00043    if(!eq(a,b)){
00044       fprintf(stderr, "Failed comparison: %s = %f != %s =%f in  %s line %i.\n", a_s.c_str(),
00045               a, b_s.c_str(), b, file.c_str(), line);
00046       exit(1);
00047    }
00048 }
00049 
00050 #define MAKE_STRING(exp) std::string(#exp)
00051 #define MyAssert(exp)  MyAssertFunc(exp, MAKE_STRING(exp), __FILE__, __LINE__);
00052 #define DblEqAssert(a,b)  DblEqAssertFunc(a,MAKE_STRING(a),b,MAKE_STRING(b), __FILE__, __LINE__);
00053 
00054 
00055 
00059 void testGetMethods(OsiTMINLPInterface &si)
00060 {
00061     CoinRelFltEq eq;// to test equality of doubles    
00062     std::cout<<"Checking get functions"<<std::endl;
00063       // Problem size
00064       MyAssert(si.getNumCols()==4);
00065       MyAssert(si.getNumRows()==3);
00066       
00067       //Check bounds on columns
00068       const double * colLow = si.getColLower();
00069       DblEqAssert(colLow[0],0.);
00070       DblEqAssert(colLow[1],0.);
00071       DblEqAssert(colLow[2],0.);
00072       DblEqAssert(colLow[3],0.);
00073       
00074       const double * colUp = si.getColUpper();
00075       MyAssert(colUp[0]>si.getInfinity());
00076       MyAssert(colUp[1]>si.getInfinity());
00077       DblEqAssert(colUp[2],1.);
00078       DblEqAssert(colUp[3],5.);
00079       //Check bounds on rows
00080       const double * rowLow = si.getRowLower();
00081       MyAssert(rowLow[0]<= -si.getInfinity());
00082       MyAssert(rowLow[1]<= -si.getInfinity());
00083       MyAssert(rowLow[2]<= -si.getInfinity());
00084                   
00085       const double * rowUp = si.getRowUpper();
00086       DblEqAssert(rowUp[0], 1./4.);
00087       DblEqAssert(rowUp[1], 0.);
00088       DblEqAssert(rowUp[2], 2.);
00089 
00090       //check objective sense
00091       MyAssert(si.getObjSense()==1);
00092       
00093       // check variables types
00094       MyAssert(si.isInteger(0)==0);
00095       MyAssert(si.isInteger(1)==0);
00096       MyAssert(si.isInteger(2)==1);
00097       MyAssert(si.isInteger(3)==1);
00098       
00099       MyAssert(si.isContinuous(0)==1);
00100       MyAssert(si.isContinuous(1)==1);
00101       MyAssert(si.isContinuous(2)==0);
00102       MyAssert(si.isContinuous(3)==0);
00103       
00104       MyAssert(si.isBinary(0)==0);
00105       MyAssert(si.isBinary(1)==0);
00106       MyAssert(si.isBinary(2)==1);
00107       MyAssert(si.isBinary(3)==0);
00108       
00109       MyAssert(si.isIntegerNonBinary(0)==0);
00110       MyAssert(si.isIntegerNonBinary(1)==0);
00111       MyAssert(si.isIntegerNonBinary(2)==0);
00112       MyAssert(si.isIntegerNonBinary(3)==1);
00113       
00114       MyAssert(si.isFreeBinary(2)==1);
00115       si.setColLower(2,1.);
00116       MyAssert(si.isFreeBinary(2)==0);
00117       si.setColLower(2,0.);
00118       
00119       //MyAssert(si.getInfinity()>1e50);
00120       std::cout<<"Test passed"<<std::endl;                  
00121 }
00122 void testOptimAndSolutionQuery(OsiTMINLPInterface & si)
00123 {
00124     CoinRelFltEq eq(1e-07);// to test equality of doubles    
00125     std::cout<<"Testing optimization methods and solution query"<<std::endl;
00126     si.initialSolve();
00127     
00128     MyAssert(si.isProvenOptimal());
00129 //    MyAssert(si.nCallOptimizeTNLP()==1);
00130     MyAssert(si.getIterationCount()>0);
00131     // Optimum of the problem is -( 3/2 + sqrt(5)/2)
00132     // with x = (1/2 + sqrt(5) y[1]=x and y[2] = 1/2 + sqrt(5)/2
00133     // (can easily be computed since constraint x-y[1]<=0 imply x = y[1] and the resulting problem has dimension 2
00134     if(!eq(si.getObjValue(),-( (3./2.) + sqrt(5.)/2.)))
00135         std::cout<<"Error in objective : "<<fabs(si.getObjValue()+( (3./2.) + sqrt(5.0)/2.))<<std::endl;
00136     
00137     //Test validity of primal solution
00138     const double * colsol = si.getColSolution();
00139     if(!eq(colsol[0],( (1./2.) + 1/sqrt(5.0))))
00140         std::cout<<"Error for y[1]  : "<<fabs(colsol[0]-( (1./2.) + 1/sqrt(5.0)))<<std::endl;
00141     if(!eq(colsol[1],( (1./2.) + 1/(2.*sqrt(5.0)))))
00142         std::cout<<"Error for y[2]  : "<<fabs(colsol[1]-( (1./2.) + 1/(2*sqrt(5.0))))<<std::endl;
00143     if(!eq(colsol[2],( (1./2.) + 1/sqrt(5.))))
00144         std::cout<<"Error for x  : "<<fabs(colsol[2]-( (1./2.) + 1/sqrt(5.0)))<<std::endl;
00145     //value of z is not tested
00146 
00147     //Test for row activity
00148     const double * rowAct = si.getRowActivity();
00149     if(!eq(rowAct[0],1./4.))
00150         std::cout<<"Error for row activity of c1 : "<<fabs(rowAct[0]-1./4.)<<std::endl;
00151     if(!eq(rowAct[1],0.))
00152         std::cout<<"Error for row activity of c2 : "<<fabs(rowAct[1])<<std::endl;
00153         
00154      //Check dual values dual for c1 = sqrt(5) c2=1 c3 not tested
00155      const double * duals = si.getRowPrice();
00156      if(!eq(duals[0],sqrt(5.0)))
00157              std::cout<<"Error dual of c1 : "<<fabs(duals[0]-sqrt(5.))<<std::endl;
00158      if(!eq(duals[1],1.))
00159              std::cout<<"Error dual of c2 : "<<fabs(duals[0]-1.)<<std::endl;
00160              
00161      std::cout<<"Test passed successfully"<<std::endl;
00162 }
00163 
00165 void testSetMethods(OsiTMINLPInterface &si)
00166 {
00167     CoinRelFltEq eq(1e-07);// to test equality of doubles    
00168     si.setColLower(2,1.);
00169     DblEqAssert(si.getColLower()[2],1.);
00170     si.initialSolve();    
00171     MyAssert(si.isProvenOptimal());
00172     DblEqAssert(si.getColSolution()[2],1.);
00173 
00174     CoinWarmStart * ws = si.getWarmStart();
00175     
00176     
00177     si.setColLower(2,0.);
00178     
00179     si.setColUpper(2,0.);
00180     DblEqAssert(si.getColUpper()[2],0.);
00181     si.setWarmStart(ws);
00182 
00183     si.resolve();
00184     MyAssert(si.isProvenOptimal());
00185     DblEqAssert(si.getColSolution()[2],0.);
00186     
00187     si.setColUpper(2,1.);
00188     delete ws;
00189 }
00190 
00191 void testOa(Bonmin::OsiTMINLPInterface &si)
00192 {
00193         CoinRelFltEq eq(1e-07);// to test equality of doubles    
00194     OsiClpSolverInterface lp;
00195     si.extractLinearRelaxation(lp);
00196 //    lp.writeMps("toy");
00197      MyAssert(lp.getNumCols()==4);
00198       MyAssert(lp.getNumRows()==3);
00199       //Check bounds on columns
00200       const double * colLow = lp.getColLower();
00201       DblEqAssert(colLow[0],0.);
00202       DblEqAssert(colLow[1],0.);
00203       DblEqAssert(colLow[2],0.);
00204       DblEqAssert(colLow[3],0.);
00205       
00206       const double * colUp = lp.getColUpper();
00207       MyAssert(colUp[0]>=lp.getInfinity());
00208       MyAssert(colUp[1]>=lp.getInfinity());
00209       DblEqAssert(colUp[2],1.);
00210       DblEqAssert(colUp[3],5.);
00211       //Check bounds on rows
00212       const double * rowLow = lp.getRowLower();
00213       std::cout<<rowLow[0]<<"\t"<<lp.getInfinity()<<std::endl;
00214       MyAssert(rowLow[0]<= -lp.getInfinity());
00215       MyAssert(rowLow[1]<= -lp.getInfinity());
00216       MyAssert(rowLow[2]<= -lp.getInfinity());
00217                   
00218       const double * rowUp = lp.getRowUpper();
00219       double sqrt5 = sqrt(5.);
00220       if(!eq(rowUp[0], 1./2. + 3./(2 * sqrt5))){
00221         double error = fabs(rowUp[0] - 1./2. - 3./(2 * sqrt5));
00222         std::cout<<"Error in OA for rowUp[0]: "
00223                  <<error<<std::endl;
00224       }
00225       DblEqAssert(rowUp[1], 0.);
00226       DblEqAssert(rowUp[2], 2.);
00227       //DblEqAssert(rowUp[3], 0.);
00228       
00229 
00230       //check objective sense
00231       MyAssert(si.getObjSense()==1);
00232       
00233       // check variables types
00234       MyAssert(si.isInteger(0)==0);
00235       MyAssert(si.isInteger(1)==0);
00236       MyAssert(si.isInteger(2)==1);
00237       MyAssert(si.isInteger(3)==1);
00238     
00239        //Now check the full matrix
00240        const CoinPackedMatrix * mat = lp.getMatrixByCol();
00241        int  inds[7] = {0, 1, 0, 2, 1, 2,  2};
00242        double vals[7] = {2. / sqrt(5.) , -1., 1./sqrt(5.), 1. ,  1. , 1., 1.};
00243        MyAssert(mat->getNumElements()==7);
00244        int k=0;
00245        for(int i = 0 ; i < si.getNumCols() ; i++)
00246        {
00247         for(int j = mat->getVectorStarts()[i] ; j < mat->getVectorStarts()[i] + mat->getVectorLengths()[i] ; j++)
00248         {
00249         MyAssert(inds[k]==mat->getIndices()[j]);
00250         if(!eq(vals[k],mat->getElements()[j])){
00251         double error = fabs(vals[k] - mat->getElements()[j]);
00252         std::cout<<"Error in OA for element of constraint matrix "<<k<<": "
00253                  <<error<<std::endl;
00254         if(error > 1e-06) throw -1;
00255       }
00256           k++;
00257         }
00258        }
00259 }
00260 
00261 void testFp(Bonmin::AmplInterface &si)
00262 {
00263         CoinRelFltEq eq(1e-07);// to test equality of doubles
00264         double x[1] = {0.};
00265         int ind[1]={1};
00266         si.solveFeasibilityProblem(1,x,ind, 0, 1);
00267         std::cout<<si.getColSolution()[0]<<std::endl;
00268          std::cout<<si.getColSolution()[1]<<std::endl;
00269        DblEqAssert(si.getColSolution()[1],(1./2.));
00270 }
00271 void interfaceTest(Ipopt::SmartPtr<TNLPSolver> solver)
00272 {
00273   /**********************************************************************************/
00274   /*   Test constructors                                                                                                              */
00275   /**********************************************************************************/
00276   std::cout<<"Test OsiTMINLPInterface with "
00277            <<solver->solverName()<<" solver"<<std::endl;
00278   // Test usefull constructor
00279 #ifdef COIN_HAS_ASL
00280   {
00281         //read a toy problem and do various tests
00282 //        var x binary;
00283 //        var z integer >= 0 <= 5;
00284 //        var y{1..2} >=0;
00285 //        
00286 //        
00287 //        minimize cost:
00288 //            - x - y[1] - y[2] ;
00289 //            
00290 //        subject to
00291 //            c1: ( y[1] - 1/2 )^2 + (y[2] - 1/2)^2 <= 1/4 ;
00292 //            c2: x - y[1] <= 0 ;
00293 //            c3: x + y[2] + z <= 2;
00294         
00295         //Setup Ipopt should be replaced if solver is changed
00296        const char * args[3] ={"name","mytoy",NULL}; //Ugly, but I don't know how to do differently
00297        const char ** argv = args;
00298        AmplInterface amplSi;
00299        amplSi.setSolver(solver);
00300        BonminSetup::registerAllOptions(solver->roptions());
00301        BonminAmplSetup bonmin;
00302        bonmin.initialize(amplSi,const_cast<char **&>(argv));
00303        OsiTMINLPInterface& si = *bonmin.nonlinearSolver();
00304        si.setWarmStartMode(2);
00305     std::cout<<"---------------------------------------------------------------------------------------------------------------------------------------------------------"
00306     <<std::endl<<"Testing useful constructor"<<std::endl
00307     <<"---------------------------------------------------------------------------------------------------------------------------------------------------------"<<std::endl;
00308       //Start of real tests
00309       testGetMethods(si);
00310       testOptimAndSolutionQuery(si);
00311       testSetMethods(si);
00312 
00313   }
00314   // Test copy constructor
00315   {
00316         //read a toy problem and do various tests
00317 //        var x binary;
00318 //        var z integer >= 0 <= 5;
00319 //        var y{1..2} >=0;
00320 //        
00321 //        
00322 //        minimize cost:
00323 //            - x - y[1] - y[2] ;
00324 //            
00325 //        subject to
00326 //            c1: ( y[1] - 1/2 )^2 + (y[2] - 1/2)^2 <= 1/4 ;
00327 //            c2: x - y[1] <= 0 ;
00328 //            c3: x + y[2] + z <= 2;
00329         
00330         //Setup should be replaced if solver is changed
00331       const char * args[3] ={"name","mytoy",NULL}; //Ugly, but I don't know how to do differently
00332       const char ** argv = args;
00333       AmplInterface amplSi;
00334       amplSi.setSolver(solver);
00335       BonminAmplSetup bonmin;
00336       bonmin.initialize(amplSi,const_cast<char **&>(argv));     
00337       OsiTMINLPInterface& si1 = *bonmin.nonlinearSolver();
00338       si1.setWarmStartMode(2);
00339       
00340       OsiTMINLPInterface si(si1);
00341     std::cout<<"---------------------------------------------------------------------------------------------------------------------------------------------------------"
00342     <<std::endl<<"Testing copy constructor"<<std::endl
00343     <<"---------------------------------------------------------------------------------------------------------------------------------------------------------"<<std::endl;
00344       //Start of real tests
00345       testGetMethods(si);
00346       testOptimAndSolutionQuery(si);
00347       testSetMethods(si);
00348   }
00349   
00350     // Test outer approximation methods
00351   {
00352         //Setup Ipopt should be replaced if solver is changed
00353         const char * args[3] ={"name","mytoy",NULL}; //Ugly, but I don't know how to do differently
00354         const char ** argv = args;
00355         AmplInterface amplSi;
00356         amplSi.setSolver(solver);
00357         BonminAmplSetup bonmin;
00358         bonmin.initialize(amplSi,const_cast<char **&>(argv));     
00359         OsiTMINLPInterface& si = *bonmin.nonlinearSolver();
00360         si.setWarmStartMode(2);
00361         std::cout<<"---------------------------------------------------------------------------------------------------------------------------------------------------------"
00362           <<std::endl<<"Testing outer approximations related methods"<<std::endl
00363           <<"---------------------------------------------------------------------------------------------------------------------------------------------------------"<<std::endl;
00364         testOa(si);
00365   }
00366   
00367   // Test Feasibility Pump methods
00368 //  {
00369 //    //Setup Ipopt should be replaced if solver is changed
00370 //    using namespace Ipopt;
00371 //    SmartPtr<Ipopt::IpoptApplication> app = new Ipopt::IpoptApplication();
00372 //    char * args[3] ={"name","toy3",NULL}; //Ugly, but I don't know how to do differently
00373 //    char ** argv = args;
00374 //    SmartPtr<TMINLP> ampl_tminlp = new AmplTMINLP(ConstPtr(app->Jnlst()),  app->Options(), argv);
00375 //    Bonmin::AmplInterface si(ampl_tminlp);
00376 //    std::cout<<"---------------------------------------------------------------------------------------------------------------------------------------------------------"
00377 //           <<std::endl<<"Testing optimization of some distance over feasible set"<<std::endl
00378 //           <<"---------------------------------------------------------------------------------------------------------------------------------------------------------"<<std::endl;
00379  //   testFp(si);
00380 //  }
00381 #endif // COIN_HAS_ASL
00382   std::cout<<"All test passed successfully"<<std::endl;
00383 } 
00384 
00385 int main()
00386 {
00387   WindowsErrorPopupBlocker();
00388 
00389   Ipopt::SmartPtr<IpoptSolver> ipopt_solver = new IpoptSolver;
00390   interfaceTest(GetRawPtr(ipopt_solver));
00391 
00392 #ifdef COIN_HAS_FILTERSQP
00393   Ipopt::SmartPtr<FilterSolver> filter_solver = new FilterSolver;
00394   interfaceTest(GetRawPtr(filter_solver));
00395 #endif
00396   return 0;
00397 }

Generated on Thu Nov 10 03:05:42 2011 by  doxygen 1.4.7