/home/coin/SVN-release/OS-2.4.2/OS/src/OSSolverInterfaces/OSBonminSolver.cpp

Go to the documentation of this file.
00001 /* $Id: OSBonminSolver.cpp 4292 2011-09-21 05:47:18Z kmartin $ */
00018 //#define DEBUG
00019 
00020 #include <iostream>
00021 
00022 
00023 #include "OSBonminSolver.h"
00024 #include "OSGeneral.h"
00025 #include "OSParameters.h"
00026 #include "OSMathUtil.h"
00027 #include "CoinFinite.hpp"
00028 #include "CoinTime.hpp"
00029 #include "BonOsiTMINLPInterface.hpp"
00030 #include "BonTMINLP.hpp"
00031 
00032 using namespace Bonmin;
00033 using namespace Ipopt;
00034 
00035 using std::cout;
00036 using std::endl;
00037 using std::ostringstream;
00038 
00039 BonminSolver::BonminSolver()
00040 {
00041     osrlwriter = new OSrLWriter();
00042     osresult = new OSResult();
00043     m_osilreader = NULL;
00044     m_osolreader = NULL;
00045     bonminErrorMsg = "";
00046 
00047 }
00048 
00049 BonminSolver::~BonminSolver()
00050 {
00051 #ifdef DEBUG
00052     cout << "inside BonminSolver destructor" << endl;
00053 #endif
00054     if(m_osilreader != NULL) delete m_osilreader;
00055     m_osilreader = NULL;
00056     if(m_osolreader != NULL) delete m_osolreader;
00057     m_osolreader = NULL;
00058 
00059     if(osrlwriter != NULL )delete osrlwriter;
00060     osrlwriter = NULL;
00061     if(osresult != NULL )
00062     {
00063         delete osresult;
00064         osresult = NULL;
00065         //cout << "DELETING OS RESULT" << endl;
00066     }
00067 #ifdef DEBUG
00068     cout << "leaving BonminSolver destructor" << endl;
00069 #endif
00070 }
00071 
00072 
00073 
00074 
00075 
00076 
00077 bool BonminProblem::get_variables_types(Index n, VariableType* var_types)
00078 {
00079     int i = 0;
00080     char *varType;
00081     varType = osinstance->getVariableTypes();
00082     n = osinstance->getVariableNumber();
00083     for(i = 0; i < n;   i++)
00084     {
00085         if( varType[i] == 'B')
00086         {
00087             var_types[i] = BINARY;
00088         }
00089         else
00090         {
00091             if( varType[i] == 'I')
00092             {
00093                 var_types[i] = INTEGER;
00094             }
00095             else
00096             {
00097                 if( varType[i] == 'C')
00098                 {
00099                     var_types[i] = CONTINUOUS;
00100                 }
00101                 else
00102                 {
00103                     throw ErrorClass("variable type not yet implemented");
00104                 }
00105             }
00106         }
00107     }
00108     return true;
00109 }
00110 
00111 bool BonminProblem::get_variables_linearity(Index n, Ipopt::TNLP::LinearityType* var_types)
00112 {
00113 
00114 
00115 
00116     std::cout << "Initialize Nonlinear Structures" << std::endl;
00117     try
00118     {
00119         osinstance->initForAlgDiff( );
00120     }
00121     catch(const ErrorClass& eclass)
00122     {
00123         bonminErrorMsg = eclass.errormsg;
00124         throw;
00125     }
00131     std::map<int, int> varIndexMap;
00132     std::map<int, int>::iterator posVarIndexMap;
00133     varIndexMap = osinstance->getAllNonlinearVariablesIndexMap( );
00134     /* first make all continuous */
00135 
00136     int i;
00137 
00138     for(i = 0; i < n; i++)
00139     {
00140         var_types[ i] = Ipopt::TNLP::LINEAR;
00141     }
00146     for(posVarIndexMap = varIndexMap.begin(); posVarIndexMap != varIndexMap.end(); ++posVarIndexMap)
00147     {
00148         std::cout <<  "Variable Index = "   << posVarIndexMap->first  << std::endl ;
00149         var_types[ posVarIndexMap->first] = Ipopt::TNLP::NON_LINEAR;
00150     }
00151     std::cout << "Number of nonlinear variables =  " << varIndexMap.size() << std::endl;
00152 
00153     return true;
00154 }
00155 
00156 bool BonminProblem::get_constraints_linearity(Index m, Ipopt::TNLP::LinearityType* const_types)
00157 {
00158 
00159     int i;
00160 
00161     for(i  = 0; i < m; i++)
00162     {
00163         const_types[ i] = Ipopt::TNLP::LINEAR;
00164     }
00165 
00166 
00167     int mm = osinstance->getNumberOfNonlinearExpressionTreeModIndexes();
00168 
00169 
00170 
00171     for(i = 0; i < mm; i++)
00172     {
00173         if(osinstance->getNonlinearExpressionTreeModIndexes()[ i] >= 0)
00174         {
00175             std::cout << osinstance->getNonlinearExpressionTreeModIndexes()[ i] << std::endl;
00176             const_types[ osinstance->getNonlinearExpressionTreeModIndexes()[ i] ] = Ipopt::TNLP::NON_LINEAR;
00177 
00178         }
00179 
00180     }
00181 
00182     return true;
00183 }
00184 
00185 // returns the size of the problem
00186 bool BonminProblem::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00187                                  Index& nnz_h_lag, TNLP::IndexStyleEnum& index_style)
00188 {
00189     //if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Bonmin NEEDS AN OBJECTIVE FUNCTION");
00190     // number of variables
00191     n = osinstance->getVariableNumber();
00192     // number of constraints
00193     m = osinstance->getConstraintNumber();
00194 #ifdef DEBUG
00195     cout << "Bonmin number variables  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << n << endl;
00196     cout << "Bonmin number constraints  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << m << endl;
00197 #endif
00198     try
00199     {
00200         osinstance->initForAlgDiff( );
00201     }
00202     catch(const ErrorClass& eclass)
00203     {
00204         bonminErrorMsg = eclass.errormsg;
00205         throw;
00206     }
00207     // use the OS Expression tree for function evaluations instead of CppAD
00208     osinstance->bUseExpTreeForFunEval = true;
00209     //std::cout << "Call sparse jacobian" << std::endl;
00210     SparseJacobianMatrix *sparseJacobian = NULL;
00211     try
00212     {
00213         sparseJacobian = osinstance->getJacobianSparsityPattern();
00214         //cout << "Get sparse Jacobian pattern" << std::endl;
00215     }
00216     catch(const ErrorClass& eclass)
00217     {
00218         bonminErrorMsg = eclass.errormsg;
00219         throw;
00220     }
00221     //std::cout << "Done calling sparse jacobian" << std::endl;
00222     nnz_jac_g = sparseJacobian->valueSize;
00223 #ifdef DEBUG
00224     cout << "nnz_jac_g  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_jac_g << endl;
00225 #endif
00226     // nonzeros in upper hessian
00227 
00228     if( (osinstance->getNumberOfNonlinearExpressions() == 0) && (osinstance->getNumberOfQuadraticTerms() == 0) )
00229     {
00230         cout << "This is a linear program"  << endl;
00231         nnz_h_lag = 0;
00232     }
00233     else
00234     {
00235         //std::cout << "Get Lagrangain Hessian Sparsity Pattern " << std::endl;
00236         SparseHessianMatrix *sparseHessian = osinstance->getLagrangianHessianSparsityPattern();
00237         //std::cout << "Done Getting Lagrangain Hessian Sparsity Pattern " << std::endl;
00238         nnz_h_lag = sparseHessian->hessDimension;
00239     }
00240 #ifdef DEBUG
00241     cout << "nnz_h_lag  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_h_lag << endl;
00242 #endif
00243     // use the C style indexing (0-based)
00244     index_style = TNLP::C_STYLE;
00245 
00246     return true;
00247 }//get_nlp_info
00248 
00249 
00250 bool  BonminProblem::get_bounds_info(Index n, Number* x_l, Number* x_u,
00251                                      Index m, Number* g_l, Number* g_u)
00252 {
00253     int i;
00254 
00255     double * mdVarLB = osinstance->getVariableLowerBounds();
00256     //std::cout << "GET BOUNDS INFORMATION FOR BONMIN !!!!!!!!!!!!!!!!!" << std::endl;
00257     // variables upper bounds
00258     double * mdVarUB = osinstance->getVariableUpperBounds();
00259 
00260     for(i = 0; i < n; i++)
00261     {
00262         x_l[ i] = mdVarLB[ i];
00263         x_u[ i] = mdVarUB[ i];
00264 #ifdef DEBUG
00265         cout << "x_l !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_l[i] << endl;
00266         cout << "x_u !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_u[i] << endl;
00267 #endif
00268     }
00269     // Bonmin interprets any number greater than nlp_upper_bound_inf as
00270     // infinity. The default value of nlp_upper_bound_inf and nlp_lower_bound_inf
00271     // is 1e19 and can be changed through bonmin options.
00272     // e.g. g_u[0] = 2e19;
00273 
00274     //constraint lower bounds
00275     double * mdConLB = osinstance->getConstraintLowerBounds();
00276     //constraint upper bounds
00277     double * mdConUB = osinstance->getConstraintUpperBounds();
00278 
00279     for(int i = 0; i < m; i++)
00280     {
00281         g_l[ i] = mdConLB[ i];
00282         g_u[ i] = mdConUB[ i];
00283 #ifdef DEBUG
00284         cout << "lower !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_l[i] << endl;
00285         cout << "upper !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_u[i] << endl;
00286 #endif
00287     }
00288     return true;
00289 }//get_bounds_info
00290 
00291 
00292 // returns the initial point for the problem
00293 bool BonminProblem::get_starting_point(Index n, bool init_x, Number* x,
00294                                        bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda,
00295                                        Number* lambda)
00296 {
00297     // Here, we assume we only have starting values for x, if you code
00298     // your own NLP, you can provide starting values for the dual variables
00299     // if you wish
00300     assert(init_x == true);
00301     assert(init_z == false);
00302     assert(init_lambda == false);
00303     int i, m1, n1;
00304 
00305 
00306 #ifdef DEBUG
00307     cout << "get initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00308 #endif
00309 
00310     //now set initial values
00311 #ifdef DEBUG
00312     cout << "get number of initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00313 #endif
00314     int k;
00315     if (osoption != NULL)
00316         m1 = osoption->getNumberOfInitVarValues();
00317     else
00318         m1 = 0;
00319 #ifdef DEBUG
00320     cout << "number of variables initialed: " << m1 << endl;
00321 #endif
00322 
00323     n1 = osinstance->getVariableNumber();
00324     bool* initialed;
00325     initialed = new bool[n1];
00326 #ifdef DEBUG
00327     cout << "number of variables in total: " << n1 << endl;
00328 #endif
00329 
00330 
00331     for(k = 0; k < n1; k++)
00332         initialed[k] = false;
00333 
00334     if (m1 > 0)
00335     {
00336 #ifdef DEBUG
00337         cout << "get initial values " << endl;
00338 #endif
00339 
00340         InitVarValue**  initVarVector = osoption->getInitVarValuesSparse();
00341 #ifdef DEBUG
00342         cout << "done " << endl;
00343 #endif
00344 
00345         double initval;
00346         try
00347         {
00348             for(k = 0; k < m1; k++)
00349             {
00350                 i = initVarVector[k]->idx;
00351                 if (initVarVector[k]->idx > n1)
00352                     throw ErrorClass ("Illegal index value in variable initialization");
00353 
00354                 initval = initVarVector[k]->value;
00355                 if (osinstance->instanceData->variables->var[i]->ub == OSDBL_MAX)
00356                 {
00357                     if (osinstance->instanceData->variables->var[i]->lb > initval)
00358                         throw ErrorClass ("Initial value outside of bounds");
00359                 }
00360                 else if (osinstance->instanceData->variables->var[i]->lb == -OSDBL_MAX)
00361                 {
00362                     if (osinstance->instanceData->variables->var[i]->ub < initval)
00363                         throw ErrorClass ("Initial value outside of bounds");
00364                 }
00365                 else
00366                 {
00367                     if ((osinstance->instanceData->variables->var[i]->lb > initval) ||
00368                             (osinstance->instanceData->variables->var[i]->ub < initval))
00369                         throw ErrorClass ("Initial value outside of bounds");
00370                 }
00371 
00372                 x[initVarVector[k]->idx] = initval;
00373                 initialed[initVarVector[k]->idx] = true;
00374             }
00375         }
00376         catch(const ErrorClass& eclass)
00377         {
00378             cout << "Error in BonminProblem::get_starting_point (OSBonminSolver.cpp)";
00379             cout << endl << endl << endl;
00380         }
00381     }  //  end if (m1 > 0)
00382 
00383     double default_initval;
00384     default_initval = 1.7171;
00385 
00386 
00387     for(k = 0; k < n1; k++)
00388     {
00389         if (!initialed[k])
00390             if (osinstance->instanceData->variables->var[k]->ub == OSDBL_MAX)
00391                 if (osinstance->instanceData->variables->var[k]->lb <= default_initval)
00392                     x[k] = default_initval;
00393                 else
00394                     x[k] = osinstance->instanceData->variables->var[k]->lb;
00395             else if (osinstance->instanceData->variables->var[k]->lb == -OSDBL_MAX)
00396                 if (osinstance->instanceData->variables->var[k]->ub >= default_initval)
00397                     x[k] = default_initval;
00398                 else
00399                     x[k] = osinstance->instanceData->variables->var[k]->ub;
00400             else if ((osinstance->instanceData->variables->var[k]->lb <= default_initval) &&
00401                      (osinstance->instanceData->variables->var[k]->ub >= default_initval))
00402                 x[k] = default_initval;
00403             else if (osinstance->instanceData->variables->var[k]->lb > default_initval)
00404                 x[k] = osinstance->instanceData->variables->var[k]->lb;
00405             else
00406                 x[k] = osinstance->instanceData->variables->var[k]->ub;
00407     }
00408     for(i = 0; i < n1; i++)
00409     {
00410 #ifdef DEBUG
00411         std::cout << "INITIAL VALUE !!!!!!!!!!!!!!!!!!!!  " << x[ i] << std::endl;
00412 #endif
00413     }
00414 
00415     osinstance->calculateAllObjectiveFunctionValues( x, true);
00416     delete[] initialed;
00417     return true;
00418 }//get_starting_point
00419 
00420 
00421 // returns the value of the objective function
00422 bool BonminProblem::eval_f(Index n, const Number* x, bool new_x, Number& obj_value)
00423 {
00424 
00425     try
00426     {
00427         if(osinstance->getObjectiveNumber() > 0)
00428         {
00429             if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") == 0)
00430             {
00431                 obj_value  = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 )[ 0];
00432             }
00433             else  // we have a max
00434             {
00435                 obj_value  = -osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 )[ 0];
00436             }
00437         }
00438     }
00439     catch(const ErrorClass& eclass)
00440     {
00441         bonminErrorMsg = eclass.errormsg;
00442         throw;
00443     }
00444     if( CoinIsnan( (double)obj_value) ) return false;
00445     return true;
00446 }
00447 
00448 bool BonminProblem::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f)
00449 {
00450     int i;
00451     double *objGrad = NULL;
00452     if(osinstance->getObjectiveNumber() > 0)
00453     {
00454         try
00455         {
00456             //objGrad = osinstance->calculateAllObjectiveFunctionGradients( const_cast<double*>(x), NULL, NULL,  new_x, 1)[ 0];
00457             //std::cout << "Calculate Objective function gradient " << std::endl;
00458             // we assume we are doing the objective function indexed by -1
00459             objGrad = osinstance->calculateObjectiveFunctionGradient( const_cast<double*>(x), NULL, NULL, -1,  new_x, 1);
00460         }
00461         catch(const ErrorClass& eclass)
00462         {
00463             bonminErrorMsg = eclass.errormsg;
00464             throw;
00465         }
00466 
00467         for(i = 0; i < n; i++)
00468         {
00469             if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") == 0)
00470             {
00471                 grad_f[ i]  = objGrad[ i];
00472             }
00473             else
00474             {
00475                 grad_f[ i]  = -objGrad[ i];
00476             }
00477             //std::cout << grad_f[ i]  << std::endl;
00478         }
00479     }
00480 //std::cout << "DONE WITH Calculate Objective function gradient " << std::endl;
00481     return true;
00482 }//eval_grad_f
00483 
00484 // return the value of the constraints: g(x)
00485 bool BonminProblem::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g)
00486 {
00487     try
00488     {
00489         double *conVals = osinstance->calculateAllConstraintFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 );
00490         int i;
00491         for(i = 0; i < m; i++)
00492         {
00493             if( CoinIsnan( (double)conVals[ i] ) ) return false;
00494             g[i] = conVals[ i]  ;
00495         }
00496         return true;
00497     }
00498     catch(const ErrorClass& eclass)
00499     {
00500         bonminErrorMsg = eclass.errormsg;
00501         throw;
00502     }
00503 }//eval_g
00504 
00505 
00506 // return the structure or values of the jacobian
00507 bool BonminProblem::eval_jac_g(Index n, const Number* x, bool new_x,
00508                                Index m, Index nele_jac, Index* iRow, Index *jCol,
00509                                Number* values)
00510 {
00511     SparseJacobianMatrix *sparseJacobian;
00512     if (values == NULL)
00513     {
00514         // return the values of the jacobian of the constraints
00515         //cout << "n: " << n << endl;
00516         //cout << "m: " << m << endl;
00517         //cout << "nele_jac: " <<  nele_jac << endl;
00518         // return the structure of the jacobian
00519         try
00520         {
00521             sparseJacobian = osinstance->getJacobianSparsityPattern();
00522         }
00523         catch(const ErrorClass& eclass)
00524         {
00525             bonminErrorMsg =  eclass.errormsg;
00526             throw;
00527         }
00528         int i = 0;
00529         int k, idx;
00530         for(idx = 0; idx < m; idx++)
00531         {
00532             for(k = *(sparseJacobian->starts + idx); k < *(sparseJacobian->starts + idx + 1); k++)
00533             {
00534                 iRow[i] = idx;
00535                 jCol[i] = *(sparseJacobian->indexes + k);
00536                 //cout << "ROW IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << iRow[i] << endl;
00537                 //cout << "COL IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << jCol[i] << endl;
00538                 i++;
00539             }
00540         }
00541     }
00542     else
00543     {
00544         //std::cout << "EVALUATING JACOBIAN" << std::endl;
00545         try
00546         {
00547             sparseJacobian = osinstance->calculateAllConstraintFunctionGradients( const_cast<double*>(x), NULL, NULL,  new_x, 1);
00548         }
00549         catch(const ErrorClass& eclass)
00550         {
00551             bonminErrorMsg = eclass.errormsg;
00552             throw;
00553         }
00554         //osinstance->getIterateResults( (double*)x, 0.0, NULL, -1, new_x,  1);
00555         for(int i = 0; i < nele_jac; i++)
00556         {
00557             values[ i] = sparseJacobian->values[i];
00558             //values[ i] = osinstance->m_mdJacValue[ i];
00559             //cout << "values[i]:!!!!!!!!!!!!  " <<  values[ i] << endl;
00560             //cout << "m_mdJacValue[ i]:!!!!!!!!!!!!  " <<  osinstance->m_mdJacValue[ i] << endl;
00561         }
00562     }
00563     return true;
00564 }//eval_jac_g
00565 
00566 //return the structure or values of the hessian
00567 bool BonminProblem::eval_h(Index n, const Number* x, bool new_x,
00568                            Number obj_factor, Index m, const Number* lambda,
00569                            bool new_lambda, Index nele_hess, Index* iRow,
00570                            Index* jCol, Number* values)
00571 {
00572 
00574     SparseHessianMatrix *sparseHessian;
00575 
00576     int i;
00577     if (values == NULL)
00578     {
00579         // return the structure. This is a symmetric matrix, fill the lower left triangle only.
00580         //cout << "get structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! "  << endl;
00581         try
00582         {
00583             sparseHessian = osinstance->getLagrangianHessianSparsityPattern( );
00584         }
00585         catch(const ErrorClass& eclass)
00586         {
00587             bonminErrorMsg = eclass.errormsg;
00588             throw;
00589         }
00590         //cout << "got structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! "  << endl;
00591         for(i = 0; i < nele_hess; i++)
00592         {
00593             iRow[i] = *(sparseHessian->hessColIdx + i);
00594             jCol[i] = *(sparseHessian->hessRowIdx + i);
00595             //cout << "ROW HESS IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << iRow[i] << endl;
00596             //cout << "COL HESS IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << jCol[i] << endl;
00597         }
00598     }
00599     else
00600     {
00601         //std::cout << "EVALUATING HESSIAN" << std::endl;
00602         // return the values. This is a symmetric matrix, fill the lower left triangle only
00603         double* objMultipliers = new double[1];
00604         objMultipliers[0] = obj_factor;
00605         try
00606         {
00607             sparseHessian = osinstance->calculateLagrangianHessian( const_cast<double*>(x), objMultipliers, const_cast<double*>(lambda) ,  new_x, 2);
00608             delete[]  objMultipliers;
00609         }
00610         catch(const ErrorClass& eclass)
00611         {
00612             bonminErrorMsg = eclass.errormsg;
00613             delete[]  objMultipliers;
00614             throw;
00615         }
00616         for(i = 0; i < nele_hess; i++)
00617         {
00618             values[ i]  = *(sparseHessian->hessValues + i);
00619         }
00620     }
00622     return true;
00623 }//eval_h
00624 
00625 bool BonminProblem::get_scaling_parameters(Number& obj_scaling,
00626         bool& use_x_scaling, Index n,
00627         Number* x_scaling,
00628         bool& use_g_scaling, Index m,
00629         Number* g_scaling)
00630 {
00631     /*
00632      * Bonmin assumes problems cast as a min
00633      * if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
00634         obj_scaling = -1;
00635     }
00636     else obj_scaling = 1;
00637     use_x_scaling = false;
00638     use_g_scaling = false;
00639     */
00640     return true;
00641 }//get_scaling_parameters
00642 
00643 
00644 
00645 
00646 void BonminProblem::finalize_solution(TMINLP::SolverReturn status_,
00647                                       Index n, const Number* x, Number obj_value)
00648 {
00649 
00650     status = status_;
00651 #ifdef DEBUG
00652     std::cout << "FINALIZE OBJ SOLUTION VALUE = " << obj_value << std::endl;
00653 #endif
00654 
00655 }
00656 
00657 
00658 
00659 
00660 void BonminSolver::buildSolverInstance() throw (ErrorClass)
00661 {
00662     try
00663     {
00664 
00665         if(osil.length() == 0 && osinstance == NULL) throw ErrorClass("there is no instance");
00666         if(osinstance == NULL)
00667         {
00668             m_osilreader = new OSiLReader();
00669             osinstance = m_osilreader->readOSiL( osil);
00670         }
00671         // Create a new instance of your nlp
00672         tminlp = new BonminProblem( osinstance, osoption);
00673         this->bCallbuildSolverInstance = true;
00674         //Now initialize from tminlp
00675         //bonminSetup.initialize( GetRawPtr(tminlp) );
00676     }
00677     catch(const ErrorClass& eclass)
00678     {
00679         std::cout << "THERE IS AN ERROR" << std::endl;
00680         osresult->setGeneralMessage( eclass.errormsg);
00681         osresult->setGeneralStatusType( "error");
00682         osrl = osrlwriter->writeOSrL( osresult);
00683         throw ErrorClass( osrl) ;
00684     }
00685 }//end buildSolverInstance()
00686 
00687 
00688 
00689 void BonminSolver::setSolverOptions() throw (ErrorClass)
00690 {
00691     try
00692     {
00693         this->bSetSolverOptions = true;
00694         bonminSetup.initializeOptionsAndJournalist();
00695         //Register an additional option -- just an example
00696         bonminSetup.roptions()->AddStringOption2("print_solution","Do we print the solution or not?",
00697                 "yes",
00698                 "no", "No, we don't.",
00699                 "yes", "Yes, we do.",
00700                 "A longer comment can be put here");
00701 
00702         // Here we can change the default value of some Bonmin or Ipopt option
00703         bonminSetup.options()->SetNumericValue("bonmin.time_limit", 5000); //changes bonmin's time limit
00704         //bonminSetup.options()->SetIntegerValue("bonmin.num_resolve_at_node", 3);
00705         //bonminSetup.options()->SetIntegerValue("bonmin.num_resolve_at_root", 3);
00706         //bonminSetup.options()->SetIntegerValue("bonmin.num_retry_unsolved_random_point", 30);
00707         //bonminSetup.options()->SetIntegerValue("print_level", 12);
00708         //bonminSetup.options()->SetIntegerValue("bonmin.bb_log_level", 4);
00709         //bonminSetup.options()->SetStringValue("bonmin.nlp_failure_behavior", "fathom");
00710         //bonminSetup.options()->SetNumericValue("bonmin.allowable_gap", -1);
00711         //bonminSetup.options()->SetNumericValue("bonmin.allowable_fraction_gap", -.1);
00712         //bonminSetup.options()->SetNumericValue("bonmin.cutoff_decr", -1);
00713         //bonminSetup.options()->SetStringValue("mu_oracle","loqo");
00714 
00715 
00716         //Here we read several option files
00717         //bonminSetup.readOptionsFile("Mybonmin.opt");
00718         //bonminSetup.readOptionsFile();// This reads the default file "bonmin.opt"
00719 
00720         // Options can also be set by using a string with a format similar to the bonmin.opt file
00721         bonminSetup.readOptionsString("bonmin.algorithm B-BB\n");
00722 
00723         //turn off a lot of output -- this can be overridden by using OSOptions
00724         bonminSetup.options()->SetIntegerValue("bonmin.bb_log_level", 0 );
00725         bonminSetup.options()->SetIntegerValue("bonmin.nlp_log_level", 0 );
00726 
00727         // Now we can obtain the value of the new option
00728         int printSolution;
00729         bonminSetup.options()->GetEnumValue("print_solution", printSolution,"");
00730         if(printSolution == 1)
00731         {
00732             tminlp->printSolutionAtEndOfAlgorithm();
00733         }
00734         //
00735         if(osoption == NULL && osol.length() > 0)
00736         {
00737             m_osolreader = new OSoLReader();
00738             osoption = m_osolreader->readOSoL( osol);
00739         }
00740 
00741         if(osoption != NULL && osoption->getNumberOfSolverOptions() > 0 )
00742         {
00743             char *pEnd;
00744             int i;
00745             std::vector<SolverOption*> optionsVector;
00746             optionsVector = osoption->getSolverOptions( "bonmin",true);
00747             int num_bonmin_options = optionsVector.size();
00748             for(i = 0; i < num_bonmin_options; i++)
00749             {
00750                 if(optionsVector[ i]->type == "numeric" )
00751                 {
00752                     std::cout << "FOUND A  NUMERIC OPTION  "  <<  os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) << std::endl;
00753                     if(optionsVector[ i]->category == "ipopt")
00754                     {
00755                         bonminSetup.options()->SetNumericValue(optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
00756                     }
00757                     else
00758                     {
00759                         if(optionsVector[ i]->category == "cbc" )
00760                         {
00761                             bonminSetup.options()->SetNumericValue("milp_solver."+optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
00762                         }
00763                         else
00764                         {
00765                             bonminSetup.options()->SetNumericValue("bonmin."+optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
00766                         }
00767                     }
00768                 }
00769                 else if(optionsVector[ i]->type == "integer" )
00770                 {
00771                     std::cout << "FOUND AN INTEGER OPTION  "  <<optionsVector[ i]->name << std::endl;
00772                     if(optionsVector[ i]->category == "ipopt")
00773                     {
00774                         bonminSetup.options()->SetIntegerValue(optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ) );
00775                     }
00776                     else
00777                     {
00778                         if(optionsVector[ i]->category == "cbc" )
00779                         {
00780                             std::cout << "SETTING INTEGER CBC OPTION" << std::endl;
00781                             bonminSetup.options()->SetIntegerValue("milp_solver."+optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ));
00782                         }
00783                         else
00784                         {
00785                             bonminSetup.options()->SetIntegerValue("bonmin."+optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() )  );
00786                         }
00787                     }
00788                 }
00789                 else if(optionsVector[ i]->type == "string" )
00790                 {
00791                     std::cout << "FOUND A STRING OPTION  "  <<optionsVector[ i]->name << std::endl;
00792                     if(optionsVector[ i]->category == "ipopt")
00793                     {
00794                         bonminSetup.options()->SetStringValue(optionsVector[ i]->name, optionsVector[ i]->value );
00795                     }
00796                     else
00797                     {
00798                         if(optionsVector[ i]->category == "cbc" )
00799                         {
00800                             bonminSetup.options()->SetStringValue("milp_solver."+optionsVector[ i]->name, optionsVector[ i]->value);
00801                         }
00802                         else
00803                         {
00804                             bonminSetup.options()->SetStringValue("bonmin."+optionsVector[ i]->name, optionsVector[ i]->value);
00805                         }
00806                     }
00807 
00808                 }
00809             }
00810         }
00811 
00812     }
00813 
00814     catch(const ErrorClass& eclass)
00815     {
00816         std::cout << "THERE IS AN ERROR" << std::endl;
00817         osresult->setGeneralMessage( eclass.errormsg);
00818         osresult->setGeneralStatusType( "error");
00819         osrl = osrlwriter->writeOSrL( osresult);
00820         throw ErrorClass( osrl) ;
00821     }
00822 }//end setSolverOptions()
00823 
00824 
00825 //void BonminSolver::solve() throw (ErrorClass) {
00826 void BonminSolver::solve() throw (ErrorClass)
00827 {
00828     if( this->bCallbuildSolverInstance == false) buildSolverInstance();
00829     if( this->bSetSolverOptions == false) setSolverOptions();
00830     try
00831     {
00832 
00833         //if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Bonmin NEEDS AN OBJECTIVE FUNCTION");
00834         //double start = CoinCpuTime();
00835         //OSiLWriter osilwriter;
00836         //cout << osilwriter.writeOSiL( osinstance) << endl;
00837         //if(osinstance->getVariableNumber() <= 0)throw ErrorClass("Bonmin requires decision variables");
00838         //double duration = CoinCpuTime() - start;
00839 
00840         try
00841         {
00842             bonminSetup.initialize( GetRawPtr(tminlp) );
00843             // bb is a Bonmin BonCbc object;;
00844             bb(  bonminSetup);  //process parameter file using Ipopt and do branch and bound using Cbc
00845 
00846         }
00847         catch(TNLPSolver::UnsolvedError *E)
00848         {
00849             //There has been a failure to solve a problem with Ipopt.
00850             std::cerr<<"Ipopt has failed to solve a problem"<<std::endl;
00851         }
00852         catch(OsiTMINLPInterface::SimpleError &E)
00853         {
00854             std::cerr<<E.className()<<"::"<<E.methodName()
00855                      <<std::endl
00856                      <<E.message()<<std::endl;
00857         }
00858         catch(CoinError &E)
00859         {
00860             std::cerr<<E.className()<<"::"<<E.methodName()
00861                      <<std::endl
00862                      <<E.message()<<std::endl;
00863         }
00864 
00865         if(( bb.model().isContinuousUnbounded() == true) && (osinstance->getNumberOfIntegerVariables() + osinstance->getNumberOfBinaryVariables() <= 0) )
00866         {
00867             std::string solutionDescription = "";
00868             std::string message = "Success";
00869             int solIdx = 0;
00870 
00871             if(osresult->setServiceName( "Bonin solver service") != true)
00872                 throw ErrorClass("OSResult error: setServiceName");
00873             if(osresult->setInstanceName(  osinstance->getInstanceName()) != true)
00874                 throw ErrorClass("OSResult error: setInstanceName");
00875             if(osresult->setVariableNumber( osinstance->getVariableNumber()  ) != true)
00876                 throw ErrorClass("OSResult error: setVariableNumber");
00877             if(osresult->setObjectiveNumber( 1) != true)
00878                 throw ErrorClass("OSResult error: setObjectiveNumber");
00879             if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00880                 throw ErrorClass("OSResult error: setConstraintNumber");
00881             if(osresult->setSolutionNumber(  1) != true)
00882                 throw ErrorClass("OSResult error: setSolutionNumer");
00883             if(osresult->setGeneralMessage( message) != true)
00884                 throw ErrorClass("OSResult error: setGeneralMessage");
00885             solutionDescription = "The problem is unbounded";
00886             osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00887             osresult->setGeneralStatusType("normal");
00888             if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00889             osrl = osrlwriter->writeOSrL( osresult);
00890             return;
00891         }
00892 
00893 
00894         if(( bb.model().isProvenInfeasible() == true)  )
00895         {
00896             std::string solutionDescription = "";
00897             std::string message = "Success";
00898             int solIdx = 0;
00899 
00900 
00901             if(osresult->setServiceName( "Bonin solver service") != true)
00902                 throw ErrorClass("OSResult error: setServiceName");
00903             if(osresult->setInstanceName(  osinstance->getInstanceName()) != true)
00904                 throw ErrorClass("OSResult error: setInstanceName");
00905             if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00906                 throw ErrorClass("OSResult error: setVariableNumer");
00907             if(osresult->setObjectiveNumber( 1) != true)
00908                 throw ErrorClass("OSResult error: setObjectiveNumber");
00909             if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00910                 throw ErrorClass("OSResult error: setConstraintNumber");
00911             if(osresult->setSolutionNumber(  1) != true)
00912                 throw ErrorClass("OSResult error: setSolutionNumber");
00913             if(osresult->setGeneralMessage( message) != true)
00914                 throw ErrorClass("OSResult error: setGeneralMessage");
00915             solutionDescription = "The problem is infeasible";
00916             osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00917             osresult->setGeneralStatusType("normal");
00918             if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00919             osrl = osrlwriter->writeOSrL( osresult);
00920             return;
00921         }
00922         status = tminlp->status;
00923         writeResult();
00924 
00925     }
00926     catch(const ErrorClass& eclass)
00927     {
00928         osresult->setGeneralMessage( eclass.errormsg);
00929         osresult->setGeneralStatusType( "error");
00930         osrl = osrlwriter->writeOSrL( osresult);
00931         throw ErrorClass( osrl) ;
00932     }
00933 }//solve
00934 
00935 
00936 
00937 
00938 void BonminSolver::writeResult()
00939 {
00940     double *x = NULL;
00941     double *z = NULL;
00942     int i = 0;
00943     int solIdx = 0;
00944     std::string solutionDescription = "";
00945     std::string message = "Bonmin solver finishes to the end.";
00946     if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(0, "Warning: this problem has zero decision variables!");
00947 
00948 
00949     try
00950     {
00951         if(osinstance->getVariableNumber()  > 0) x = new double[osinstance->getVariableNumber() ];
00952         z = new double[1];
00953         // resultHeader information
00954         if(osresult->setServiceName( OSgetVersionInfo()) != true)
00955             throw ErrorClass("OSResult error: setServiceName");
00956         if(osresult->setSolverInvoked( "COIN-OR Bonmin") != true)
00957             throw ErrorClass("OSResult error: setSolverInvoked");
00958         if(osresult->setInstanceName(  osinstance->getInstanceName()) != true)
00959             throw ErrorClass("OSResult error: setInstanceName");
00960         //if(osresult->setJobID( osoption->jobID) != true)
00961         //      throw ErrorClass("OSResult error: setJobID");
00962         // set basic problem parameters
00963 
00964         if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00965             throw ErrorClass("OSResult error: setVariableNumer");
00966         if(osresult->setObjectiveNumber( 1) != true)
00967             throw ErrorClass("OSResult error: setObjectiveNumber");
00968         if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00969             throw ErrorClass("OSResult error: setConstraintNumber");
00970         if(osresult->setSolutionNumber(  1) != true)
00971             throw ErrorClass("OSResult error: setSolutionNumer");
00972         if(osresult->setGeneralMessage( message) != true)
00973             throw ErrorClass("OSResult error: setGeneralMessage");
00974 
00975         switch( status)
00976         {
00977         case  TMINLP::SUCCESS:
00978             solutionDescription = "SUCCESS[BONMIN]: Algorithm terminated normally at a locally optimal point, satisfying the convergence tolerances.";
00979             //std::cout << solutionDescription << std::endl;
00980             osresult->setSolutionStatus(solIdx,  "locallyOptimal", solutionDescription);
00981             /* Retrieve the solution */
00982             if(osinstance->getObjectiveNumber() > 0)
00983             {
00984                 *(z + 0)  =  osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(bb.bestSolution()), true)[ 0];
00985                 osresult->setObjectiveValuesDense(solIdx, z);
00986             }
00987             if( osinstance->getVariableNumber() > 0)
00988             {
00989                 for(i=0; i < osinstance->getVariableNumber(); i++)
00990                 {
00991                     *(x + i) = bb.bestSolution()[i];
00992                     //std::cout <<  *(x + i)  << std::endl;
00993                 }
00994                 osresult->setPrimalVariableValuesDense(solIdx, x);
00995             }
00996             break;
00997 
00998         case TMINLP::LIMIT_EXCEEDED:
00999             solutionDescription = "LIMIT_EXCEEDED[BONMIN]: A resource limit was exceeded, we provide the current solution.";
01000             //std::cout << solutionDescription << std::endl;
01001             osresult->setSolutionStatus(solIdx,  "other", solutionDescription);
01002             //osresult->setPrimalVariableValuesDense(solIdx, const_cast<double*>(x));
01003             //osresult->setDualVariableValuesDense(solIdx, const_cast<double*>( lambda));
01004             /* Retrieve the solution */
01005             if(osinstance->getObjectiveNumber() > 0)
01006             {
01007                 *(z + 0)  =  osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(bb.bestSolution()), true)[ 0];
01008                 osresult->setObjectiveValuesDense(solIdx, z);
01009             }
01010             if( osinstance->getVariableNumber() > 0)
01011             {
01012                 for(i=0; i < osinstance->getVariableNumber(); i++)
01013                 {
01014                     *(x + i) = bb.bestSolution()[i];
01015                     //std::cout <<  *(x + i)  << std::endl;
01016                 }
01017                 osresult->setPrimalVariableValuesDense(solIdx, x);
01018             }
01019             break;
01020 
01021         case TMINLP::MINLP_ERROR:
01022             solutionDescription = "MINLP_ERROR [BONMIN]: Algorithm stopped with unspecified error.";
01023             //std::cout << solutionDescription << std::endl;
01024             osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
01025 
01026             break;
01027 
01028         case TMINLP::CONTINUOUS_UNBOUNDED:
01029             solutionDescription = "CONTINUOUS_UNBOUNDED [BONMIN]: The continuous relaxation is unbounded, the MINLP may or may not be unbounded.";
01030             //std::cout << solutionDescription << std::endl;
01031             osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
01032 
01033             break;
01034 
01035 
01036         case TMINLP::INFEASIBLE:
01037             solutionDescription = "INFEASIBLE [BONMIN]: Problem may be infeasible.";
01038             //std::cout << solutionDescription << std::endl;
01039             osresult->setSolutionStatus(solIdx,  "infeasible", solutionDescription);
01040             break;
01041 
01042 
01043         default:
01044             solutionDescription = "OTHER[BONMIN]: other unknown solution status from Bonmin solver";
01045             //std::cout << solutionDescription << std::endl;
01046             osresult->setSolutionStatus(solIdx,  "other", solutionDescription);
01047         }//switch end
01048         osresult->setGeneralStatusType("normal");
01049         osrl = osrlwriter->writeOSrL( osresult);
01050         if(osinstance->getVariableNumber() > 0 ) delete[] x;
01051         x = NULL;
01052         delete[] z;
01053         z = NULL;
01054     }//end try
01055 
01056 
01057     catch(const ErrorClass& eclass)
01058     {
01059         if(osinstance->getVariableNumber()  > 0) delete[] x;
01060         x = NULL;
01061         delete[] z;
01062         z = NULL;
01063         osresult->setGeneralMessage( eclass.errormsg);
01064         osresult->setGeneralStatusType( "error");
01065         osrl = osrlwriter->writeOSrL( osresult);
01066         throw ErrorClass( osrl) ;
01067     }
01068 
01069 
01070 }// end writeResult()
01071 
01072 
01073 void BonminSolver::dataEchoCheck()
01074 {
01075 
01076     int i;
01077 
01078     // print out problem parameters
01079     cout << "This is problem:  " << osinstance->getInstanceName() << endl;
01080     cout << "The problem source is:  " << osinstance->getInstanceSource() << endl;
01081     cout << "The problem description is:  " << osinstance->getInstanceDescription() << endl;
01082     cout << "number of variables = " << osinstance->getVariableNumber() << endl;
01083     cout << "number of Rows = " << osinstance->getConstraintNumber() << endl;
01084 
01085     // print out the variable information
01086     if(osinstance->getVariableNumber() > 0)
01087     {
01088         for(i = 0; i < osinstance->getVariableNumber(); i++)
01089         {
01090             if(osinstance->getVariableNames() != NULL) cout << "variable Names  " << osinstance->getVariableNames()[ i]  << endl;
01091             if(osinstance->getVariableTypes() != NULL) cout << "variable Types  " << osinstance->getVariableTypes()[ i]  << endl;
01092             if(osinstance->getVariableLowerBounds() != NULL) cout << "variable Lower Bounds  " << osinstance->getVariableLowerBounds()[ i]  << endl;
01093             if(osinstance->getVariableUpperBounds() != NULL) cout << "variable Upper Bounds  " <<  osinstance->getVariableUpperBounds()[i] << endl;
01094         }
01095     }
01096 
01097     // print out objective function information
01098     if(osinstance->getVariableNumber() > 0 || osinstance->instanceData->objectives->obj != NULL || osinstance->instanceData->objectives->numberOfObjectives > 0)
01099     {
01100         if( osinstance->getObjectiveMaxOrMins()[0] == "min")  cout <<  "problem is a minimization" << endl;
01101         else cout <<  "problem is a maximization" << endl;
01102         for(i = 0; i < osinstance->getVariableNumber(); i++)
01103         {
01104             cout << "OBJ COEFFICIENT =  " <<  osinstance->getDenseObjectiveCoefficients()[0][i] << endl;
01105         }
01106     }
01107     // print out constraint information
01108     if(osinstance->getConstraintNumber() > 0)
01109     {
01110         for(i = 0; i < osinstance->getConstraintNumber(); i++)
01111         {
01112             if(osinstance->getConstraintNames() != NULL) cout << "row name = " << osinstance->getConstraintNames()[i] <<  endl;
01113             if(osinstance->getConstraintLowerBounds() != NULL) cout << "row lower bound = " << osinstance->getConstraintLowerBounds()[i] <<  endl;
01114             if(osinstance->getConstraintUpperBounds() != NULL) cout << "row upper bound = " << osinstance->getConstraintUpperBounds()[i] <<  endl;
01115         }
01116     }
01117 
01118     // print out linear constraint data
01119     cout << endl;
01120     cout << "number of nonzeros =  " << osinstance->getLinearConstraintCoefficientNumber() << endl;
01121     for(i = 0; i <= osinstance->getVariableNumber(); i++)
01122     {
01123         cout << "Start Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->starts[ i] << endl;
01124     }
01125     cout << endl;
01126     for(i = 0; i < osinstance->getLinearConstraintCoefficientNumber(); i++)
01127     {
01128         cout << "Index Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->indexes[i] << endl;
01129         cout << "Nonzero Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->values[i] << endl;
01130     }
01131 
01132     // print out quadratic data
01133     cout << "number of qterms =  " <<  osinstance->getNumberOfQuadraticTerms() << endl;
01134     for(int i = 0; i <  osinstance->getNumberOfQuadraticTerms(); i++)
01135     {
01136         cout << "Row Index = " <<  osinstance->getQuadraticTerms()->rowIndexes[i] << endl;
01137         cout << "Var Index 1 = " << osinstance->getQuadraticTerms()->varOneIndexes[ i] << endl;
01138         cout << "Var Index 2 = " << osinstance->getQuadraticTerms()->varTwoIndexes[ i] << endl;
01139         cout << "Coefficient = " << osinstance->getQuadraticTerms()->coefficients[ i] << endl;
01140     }
01141 } // end dataEchoCheck
01142 
01143 
01144 BonminProblem::BonminProblem(OSInstance *osinstance_,  OSOption *osoption_)
01145 {
01146     osinstance = osinstance_;
01147     osoption = osoption_;
01148     printSol_ = false;
01149 }
01150 
01151 BonminProblem::~BonminProblem()
01152 {
01153 
01154 }
01155 
01156 
01157 

Generated on Wed Nov 30 03:04:23 2011 by  doxygen 1.4.7