/home/coin/SVN-release/OS-2.0.0/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 "BonAmplInterface.hpp"
00011 #include "OsiClpSolverInterface.hpp"
00012 #include "BonTMINLP.hpp"
00013 #include "IpIpoptApplication.hpp"
00014 #ifdef COIN_HAS_ASL
00015 #include "BonAmplTMINLP.hpp"
00016 #include "BonAmplSetup.hpp"
00017 #endif
00018 
00019 #include "BonIpoptSolver.hpp"
00020 #include "BonminConfig.h"
00021 
00022 #ifdef COIN_HAS_FILTERSQP
00023 #include "BonFilterSolver.hpp"
00024 #endif
00025 
00026 #include "CoinError.hpp"
00027 
00028 #include <string>
00029 #include <cmath>
00030 using namespace Bonmin;
00031 
00032 void MyAssertFunc(bool c, const std::string &s, const std::string&  file, unsigned int line){
00033    if(c != true){
00034       fprintf(stderr, "Failed MyAssertion: %s in %s line %i.\n", s.c_str(), file.c_str(), line);
00035       exit(1);
00036    }
00037 }
00038 
00039 void DblEqAssertFunc(const double& a, const std::string &a_s, const double&b, const std::string& b_s,
00040                   const std::string&  file, unsigned int line){
00041    CoinRelFltEq eq;
00042    if(!eq(a,b)){
00043       fprintf(stderr, "Failed comparison: %s = %f != %s =%f in  %s line %i.\n", a_s.c_str(),
00044               a, b_s.c_str(), b, file.c_str(), line);
00045       exit(1);
00046    }
00047 }
00048 
00049 #define MAKE_STRING(exp) std::string(#exp)
00050 #define MyAssert(exp)  MyAssertFunc(exp, MAKE_STRING(exp), __FILE__, __LINE__);
00051 #define DblEqAssert(a,b)  DblEqAssertFunc(a,MAKE_STRING(a),b,MAKE_STRING(b), __FILE__, __LINE__);
00052 
00053 
00054 
00058 void testGetMethods(OsiTMINLPInterface &si)
00059 {
00060     CoinRelFltEq eq;// to test equality of doubles    
00061     std::cout<<"Checking get functions"<<std::endl;
00062       // Problem size
00063       MyAssert(si.getNumCols()==4);
00064       MyAssert(si.getNumRows()==3);
00065       
00066       //Check bounds on columns
00067       const double * colLow = si.getColLower();
00068       DblEqAssert(colLow[0],0.);
00069       DblEqAssert(colLow[1],0.);
00070       DblEqAssert(colLow[2],0.);
00071       DblEqAssert(colLow[3],0.);
00072       
00073       const double * colUp = si.getColUpper();
00074       MyAssert(colUp[0]>si.getInfinity());
00075       MyAssert(colUp[1]>si.getInfinity());
00076       DblEqAssert(colUp[2],1.);
00077       DblEqAssert(colUp[3],5.);
00078       //Check bounds on rows
00079       const double * rowLow = si.getRowLower();
00080       MyAssert(rowLow[0]<= -si.getInfinity());
00081       MyAssert(rowLow[1]<= -si.getInfinity());
00082       MyAssert(rowLow[2]<= -si.getInfinity());
00083                   
00084       const double * rowUp = si.getRowUpper();
00085       DblEqAssert(rowUp[0], 1./4.);
00086       DblEqAssert(rowUp[1], 0.);
00087       DblEqAssert(rowUp[2], 2.);
00088 
00089       //check objective sense
00090       MyAssert(si.getObjSense()==1);
00091       
00092       // check variables types
00093       MyAssert(si.isInteger(0)==0);
00094       MyAssert(si.isInteger(1)==0);
00095       MyAssert(si.isInteger(2)==1);
00096       MyAssert(si.isInteger(3)==1);
00097       
00098       MyAssert(si.isContinuous(0)==1);
00099       MyAssert(si.isContinuous(1)==1);
00100       MyAssert(si.isContinuous(2)==0);
00101       MyAssert(si.isContinuous(3)==0);
00102       
00103       MyAssert(si.isBinary(0)==0);
00104       MyAssert(si.isBinary(1)==0);
00105       MyAssert(si.isBinary(2)==1);
00106       MyAssert(si.isBinary(3)==0);
00107       
00108       MyAssert(si.isIntegerNonBinary(0)==0);
00109       MyAssert(si.isIntegerNonBinary(1)==0);
00110       MyAssert(si.isIntegerNonBinary(2)==0);
00111       MyAssert(si.isIntegerNonBinary(3)==1);
00112       
00113       MyAssert(si.isFreeBinary(2)==1);
00114       si.setColLower(2,1.);
00115       MyAssert(si.isFreeBinary(2)==0);
00116       si.setColLower(2,0.);
00117       
00118       //MyAssert(si.getInfinity()>1e50);
00119       std::cout<<"Test passed"<<std::endl;                  
00120 }
00121 void testOptimAndSolutionQuery(OsiTMINLPInterface & si)
00122 {
00123     CoinRelFltEq eq(1e-07);// to test equality of doubles    
00124     std::cout<<"Testing optimization methods and solution query"<<std::endl;
00125     si.initialSolve();
00126     
00127     MyAssert(si.isProvenOptimal());
00128 //    MyAssert(si.nCallOptimizeTNLP()==1);
00129     MyAssert(si.getIterationCount()>0);
00130     // Optimum of the problem is -( 3/2 + sqrt(5)/2)
00131     // with x = (1/2 + sqrt(5) y[1]=x and y[2] = 1/2 + sqrt(5)/2
00132     // (can easily be computed since constraint x-y[1]<=0 imply x = y[1] and the resulting problem has dimension 2
00133     if(!eq(si.getObjValue(),-( (3./2.) + sqrt(5.)/2.)))
00134         std::cout<<"Error in objective : "<<fabs(si.getObjValue()+( (3./2.) + sqrt(5.0)/2.))<<std::endl;
00135     
00136     //Test validity of primal solution
00137     const double * colsol = si.getColSolution();
00138     if(!eq(colsol[0],( (1./2.) + 1/sqrt(5.0))))
00139         std::cout<<"Error for y[1]  : "<<fabs(colsol[0]-( (1./2.) + 1/sqrt(5.0)))<<std::endl;
00140     if(!eq(colsol[1],( (1./2.) + 1/(2.*sqrt(5.0)))))
00141         std::cout<<"Error for y[2]  : "<<fabs(colsol[1]-( (1./2.) + 1/(2*sqrt(5.0))))<<std::endl;
00142     if(!eq(colsol[2],( (1./2.) + 1/sqrt(5.))))
00143         std::cout<<"Error for x  : "<<fabs(colsol[2]-( (1./2.) + 1/sqrt(5.0)))<<std::endl;
00144     //value of z is not tested
00145 
00146     //Test for row activity
00147     const double * rowAct = si.getRowActivity();
00148     if(!eq(rowAct[0],1./4.))
00149         std::cout<<"Error for row activity of c1 : "<<fabs(rowAct[0]-1./4.)<<std::endl;
00150     if(!eq(rowAct[1],0.))
00151         std::cout<<"Error for row activity of c2 : "<<fabs(rowAct[1])<<std::endl;
00152         
00153      //Check dual values dual for c1 = sqrt(5) c2=1 c3 not tested
00154      const double * duals = si.getRowPrice();
00155      if(!eq(duals[0],sqrt(5.0)))
00156              std::cout<<"Error dual of c1 : "<<fabs(duals[0]-sqrt(5.))<<std::endl;
00157      if(!eq(duals[1],1.))
00158              std::cout<<"Error dual of c2 : "<<fabs(duals[0]-1.)<<std::endl;
00159              
00160      std::cout<<"Test passed successfully"<<std::endl;
00161 }
00162 
00164 void testSetMethods(OsiTMINLPInterface &si)
00165 {
00166     CoinRelFltEq eq(1e-07);// to test equality of doubles    
00167     si.setColLower(2,1.);
00168     DblEqAssert(si.getColLower()[2],1.);
00169     si.initialSolve();    
00170     MyAssert(si.isProvenOptimal());
00171     DblEqAssert(si.getColSolution()[2],1.);
00172 
00173     CoinWarmStart * ws = si.getWarmStart();
00174     
00175     
00176     si.setColLower(2,0.);
00177     
00178     si.setColUpper(2,0.);
00179     DblEqAssert(si.getColUpper()[2],0.);
00180     si.setWarmStart(ws);
00181 
00182     si.resolve();
00183     MyAssert(si.isProvenOptimal());
00184     DblEqAssert(si.getColSolution()[2],0.);
00185     
00186     si.setColUpper(2,1.);
00187     delete ws;
00188 }
00189 
00190 void testOa(Bonmin::OsiTMINLPInterface &si)
00191 {
00192         CoinRelFltEq eq(1e-07);// to test equality of doubles    
00193     OsiClpSolverInterface lp;
00194     si.extractLinearRelaxation(lp);
00195 //    lp.writeMps("toy");
00196      MyAssert(lp.getNumCols()==4);
00197       MyAssert(lp.getNumRows()==3);
00198       //Check bounds on columns
00199       const double * colLow = lp.getColLower();
00200       DblEqAssert(colLow[0],0.);
00201       DblEqAssert(colLow[1],0.);
00202       DblEqAssert(colLow[2],0.);
00203       DblEqAssert(colLow[3],0.);
00204       
00205       const double * colUp = lp.getColUpper();
00206       MyAssert(colUp[0]>=lp.getInfinity());
00207       MyAssert(colUp[1]>=lp.getInfinity());
00208       DblEqAssert(colUp[2],1.);
00209       DblEqAssert(colUp[3],5.);
00210       //Check bounds on rows
00211       const double * rowLow = lp.getRowLower();
00212       std::cout<<rowLow[0]<<"\t"<<lp.getInfinity()<<std::endl;
00213       MyAssert(rowLow[0]<= -lp.getInfinity());
00214       MyAssert(rowLow[1]<= -lp.getInfinity());
00215       MyAssert(rowLow[2]<= -lp.getInfinity());
00216                   
00217       const double * rowUp = lp.getRowUpper();
00218       double sqrt5 = sqrt(5.);
00219       if(!eq(rowUp[0], 1./2. + 3./(2 * sqrt5))){
00220         double error = fabs(rowUp[0] - 1./2. - 3./(2 * sqrt5));
00221         std::cout<<"Error in OA for rowUp[0]: "
00222                  <<error<<std::endl;
00223       }
00224       DblEqAssert(rowUp[1], 0.);
00225       DblEqAssert(rowUp[2], 2.);
00226       //DblEqAssert(rowUp[3], 0.);
00227       
00228 
00229       //check objective sense
00230       MyAssert(si.getObjSense()==1);
00231       
00232       // check variables types
00233       MyAssert(si.isInteger(0)==0);
00234       MyAssert(si.isInteger(1)==0);
00235       MyAssert(si.isInteger(2)==1);
00236       MyAssert(si.isInteger(3)==1);
00237     
00238        //Now check the full matrix
00239        const CoinPackedMatrix * mat = lp.getMatrixByCol();
00240        int  inds[7] = {0, 1, 0, 2, 1, 2,  2};
00241        double vals[7] = {2. / sqrt(5.) , -1., 1./sqrt(5.), 1. ,  1. , 1., 1.};
00242        MyAssert(mat->getNumElements()==7);
00243        int k=0;
00244        for(int i = 0 ; i < si.getNumCols() ; i++)
00245        {
00246         for(int j = mat->getVectorStarts()[i] ; j < mat->getVectorStarts()[i] + mat->getVectorLengths()[i] ; j++)
00247         {
00248         MyAssert(inds[k]==mat->getIndices()[j]);
00249         if(!eq(vals[k],mat->getElements()[j])){
00250         double error = fabs(vals[k] - mat->getElements()[j]);
00251         std::cout<<"Error in OA for element of constraint matrix "<<k<<": "
00252                  <<error<<std::endl;
00253         if(error > 1e-06) throw -1;
00254       }
00255           k++;
00256         }
00257        }
00258 }
00259 
00260 void testFp(Bonmin::AmplInterface &si)
00261 {
00262         CoinRelFltEq eq(1e-07);// to test equality of doubles
00263         OsiCuts cuts;
00264         double x[1] = {0.};
00265         int ind[1]={1};
00266         si.getFeasibilityOuterApproximation(1,x,ind,cuts, 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.setExposeWarmStart(true);
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.setExposeWarmStart(true);
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.setExposeWarmStart(true);
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_FSQP
00393   Ipopt::SmartPtr<FilterSolver> filter_solver = new FilterSolver;
00394   interfaceTest(GetRawPtr(filter_solver));
00395 #endif
00396   return 0;
00397 }

Generated on Mon Aug 3 03:02:19 2009 by  doxygen 1.4.7