/home/coin/SVN-release/OS-2.3.2/OS/v2.0/OSBonminSolver.cpp

Go to the documentation of this file.
00001 
00018 #include <iostream>
00019 
00020 
00021 #include "OSBonminSolver.h"
00022 #include "OSDataStructures.h"
00023 #include "OSParameters.h" 
00024 #include "OSCommonUtil.h"
00025 #include "OSMathUtil.h"
00026 
00027 #include "BonOsiTMINLPInterface.hpp"
00028 #include "BonCbc.hpp"
00029 
00030 #include "CoinTime.hpp"
00031 
00032 using std::cout; 
00033 using std::endl; 
00034 using std::ostringstream;
00035 
00036 
00037 BonminSolver::BonminSolver() {
00038         osrlwriter = new OSrLWriter();
00039         osresult = new OSResult();
00040         m_osilreader = NULL;
00041         m_osolreader = NULL;
00042         bonminErrorMsg = "";
00043 
00044 } 
00045 
00046 BonminSolver::~BonminSolver() {
00047         #ifdef DEBUG
00048         cout << "inside BonminSolver destructor" << endl;
00049         #endif
00050         if(m_osilreader != NULL) delete m_osilreader;
00051         m_osilreader = NULL;
00052         if(m_osolreader != NULL) delete m_osolreader;
00053         m_osolreader = NULL;
00054         delete osresult;
00055         osresult = NULL;
00056         delete osrlwriter;
00057         osrlwriter = NULL;
00058         //delete osinstance;
00059         //osinstance = NULL;
00060         #ifdef DEBUG
00061         cout << "leaving BonminSolver destructor" << endl;
00062         #endif
00063 }
00064 
00065 
00066 
00067 
00068 
00069 
00070 bool BonminProblem::get_variables_types(Index n, VariableType* var_types){
00071         int i = 0;
00072         char *varType;
00073         varType = osinstance->getVariableTypes();
00074         n = osinstance->getVariableNumber();
00075         for(i = 0; i < n;   i++){
00076                 if( varType[i] == 'B') {
00077                         var_types[i] = BINARY;
00078                 }
00079                 else{
00080                         if( varType[i] == 'I') {
00081                                 var_types[i] = INTEGER;
00082                         }
00083                         else{
00084                                 var_types[i] = CONTINUOUS;
00085                         }       
00086                 }
00087         }
00088         return true;
00089 }
00090 
00091 bool BonminProblem::get_variables_linearity(Index n, Ipopt::TNLP::LinearityType* var_types){
00092         
00093         
00094         
00095         std::cout << "Initialize Nonlinear Structures" << std::endl;
00096         try{
00097                 osinstance->initForAlgDiff( );
00098         }
00099         catch(const ErrorClass& eclass){
00100                 bonminErrorMsg = eclass.errormsg;
00101                 throw;  
00102         }       
00108         std::map<int, int> varIndexMap;
00109         std::map<int, int>::iterator posVarIndexMap;
00110         varIndexMap = osinstance->getAllNonlinearVariablesIndexMap( );
00111         /* first make all continuous */
00112         
00113         int i;
00114         
00115         for(i = 0; i < n; i++){
00116                 var_types[ i] = Ipopt::TNLP::LINEAR;    
00117         }
00122         for(posVarIndexMap = varIndexMap.begin(); posVarIndexMap != varIndexMap.end(); ++posVarIndexMap){
00123                         std::cout <<  "Variable Index = "   << posVarIndexMap->first  << std::endl ;
00124                         var_types[ posVarIndexMap->first] = Ipopt::TNLP::NON_LINEAR;
00125         }
00126         std::cout << "Number of nonlinear variables =  " << varIndexMap.size() << std::endl;    
00127 
00128         return true;
00129 }
00130 
00131 bool BonminProblem::get_constraints_linearity(Index m, Ipopt::TNLP::LinearityType* const_types){
00132         
00133         int i;
00134         
00135         for(i  = 0; i < m; i++){
00136                 const_types[ i] = Ipopt::TNLP::LINEAR;
00137         }
00138         
00139         
00140         int mm = osinstance->getNumberOfNonlinearExpressionTreeModIndexes();
00141         
00142 
00143         
00144         for(i = 0; i < mm; i++){
00145                 if(osinstance->getNonlinearExpressionTreeModIndexes()[ i] >= 0){
00146                         std::cout << osinstance->getNonlinearExpressionTreeModIndexes()[ i] << std::endl;
00147                         const_types[ osinstance->getNonlinearExpressionTreeModIndexes()[ i] ] = Ipopt::TNLP::NON_LINEAR;
00148                         
00149                 }
00150                 
00151         }
00152                 
00153         return true;
00154 }
00155 
00156 // returns the size of the problem
00157 bool BonminProblem::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00158                              Index& nnz_h_lag, TNLP::IndexStyleEnum& index_style)
00159 {
00160         if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Bonmin NEEDS AN OBJECTIVE FUNCTION");                       
00161         // number of variables
00162         n = osinstance->getVariableNumber();
00163         // number of constraints
00164         m = osinstance->getConstraintNumber();
00165         cout << "Bonmin number variables  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << n << endl;
00166         cout << "Bonmin number constraints  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << m << endl;
00167         try{
00168                 osinstance->initForAlgDiff( );
00169         }
00170         catch(const ErrorClass& eclass){
00171                 bonminErrorMsg = eclass.errormsg;
00172                 throw;  
00173         }       
00174         // use the OS Expression tree for function evaluations instead of CppAD
00175         osinstance->bUseExpTreeForFunEval = true;
00176         //std::cout << "Call sparse jacobian" << std::endl;
00177         SparseJacobianMatrix *sparseJacobian = NULL;
00178         try{
00179                 sparseJacobian = osinstance->getJacobianSparsityPattern();
00180         }
00181         catch(const ErrorClass& eclass){
00182                 bonminErrorMsg = eclass.errormsg;
00183                 throw;  
00184         }
00185         //std::cout << "Done calling sparse jacobian" << std::endl;
00186         nnz_jac_g = sparseJacobian->valueSize;
00187         cout << "nnz_jac_g  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_jac_g << endl;  
00188         // nonzeros in upper hessian
00189         
00190         if( (osinstance->getNumberOfNonlinearExpressions() == 0) && (osinstance->getNumberOfQuadraticTerms() == 0) ) {
00191                 cout << "This is a linear program"  << endl;
00192                 nnz_h_lag = 0;
00193         }
00194         else{
00195                 //std::cout << "Get Lagrangian Hessian Sparsity Pattern " << std::endl;
00196                 SparseHessianMatrix *sparseHessian = osinstance->getLagrangianHessianSparsityPattern();
00197                 //std::cout << "Done Getting Lagrangian Hessian Sparsity Pattern " << std::endl;
00198                 nnz_h_lag = sparseHessian->hessDimension;
00199         }
00200         cout << "print nnz_h_lag (OSBonminSolver.cpp)" << endl; 
00201         cout << "nnz_h_lag  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_h_lag << endl;  
00202         cout << "set index_style (OSBonminSolver.cpp)" << endl; 
00203         // use the C style indexing (0-based)
00204         index_style = TNLP::C_STYLE;
00205         cout << "return from get_nlp_info (OSBonminSolver.cpp)" << endl;        
00206   
00208 
00209   return true;
00210 }//get_nlp_info
00211 
00212 
00213 bool  BonminProblem::get_bounds_info(Index n, Number* x_l, Number* x_u,
00214                                 Index m, Number* g_l, Number* g_u){
00215         int i; 
00216 
00217         double * mdVarLB = osinstance->getVariableLowerBounds();
00218         //std::cout << "GET BOUNDS INFORMATION FOR BONMIN !!!!!!!!!!!!!!!!!" << std::endl;
00219         // variables upper bounds
00220         double * mdVarUB = osinstance->getVariableUpperBounds();
00221 
00222         for(i = 0; i < n; i++){
00223                 x_l[ i] = mdVarLB[ i];
00224                 x_u[ i] = mdVarUB[ i];
00225                 cout << "x_l !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_l[i] << endl;
00226                 cout << "x_u !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_u[i] << endl;
00227         }
00228         // Bonmin interprets any number greater than nlp_upper_bound_inf as
00229         // infinity. The default value of nlp_upper_bound_inf and nlp_lower_bound_inf
00230         // is 1e19 and can be changed through bonmin options.
00231         // e.g. g_u[0] = 2e19;
00232 
00233         //constraint lower bounds
00234         double * mdConLB = osinstance->getConstraintLowerBounds();
00235         //constraint upper bounds
00236         double * mdConUB = osinstance->getConstraintUpperBounds();
00237         
00238         for(int i = 0; i < m; i++){
00239                 g_l[ i] = mdConLB[ i];
00240                 g_u[ i] = mdConUB[ i];
00241                 cout << "lower !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_l[i] << endl;
00242                 cout << "upper !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_u[i] << endl;
00243         }  
00244         return true;
00245 }//get_bounds_info
00246 
00247 
00248 // returns the initial point for the problem
00249 bool BonminProblem::get_starting_point(Index n, bool init_x, Number* x,
00250      bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda,
00251      Number* lambda) {
00252         // Here, we assume we only have starting values for x, if you code
00253         // your own NLP, you can provide starting values for the dual variables
00254         // if you wish
00255         assert(init_x == true);
00256         assert(init_z == false);
00257         assert(init_lambda == false);
00258         int i, m1, n1;
00259 
00260 #ifdef DEBUG
00261         cout << "get initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00262 #endif
00263         
00264         //now set initial values
00265 #ifdef DEBUG
00266         cout << "get number of initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00267         cout << "Is osoption = NULL? " << (osoption == NULL) << endl;
00268 #endif
00269         int k;
00270         if (osoption != NULL)
00271                 m1 = osoption->getNumberOfInitVarValues();
00272         else
00273                 m1 = 0;
00274 #ifdef DEBUG
00275         cout << "number of variables initialed: " << m1 << endl;
00276 #endif
00277 
00278         n1 = osinstance->getVariableNumber();
00279         bool* initialed;
00280         initialed = new bool[n1];
00281 #ifdef DEBUG
00282         cout << "number of variables in total: " << n1 << endl;
00283 #endif
00284 
00285 
00286         for(k = 0; k < n1; k++)
00287                 initialed[k] = false;
00288 
00289         if (m1 > 0)
00290         {
00291 #ifdef DEBUG
00292                         cout << "get initial values " << endl;
00293 #endif
00294 
00295                 InitVarValue**  initVarVector = osoption->getInitVarValuesSparse();
00296 #ifdef DEBUG
00297                 cout << "done " << endl;
00298 #endif
00299 
00300                 double initval;
00301                 try 
00302                 {
00303                         for(k = 0; k < m1; k++)
00304                         {       
00305                                 i = initVarVector[k]->idx;
00306                                 if (initVarVector[k]->idx > n1)
00307                                         throw ErrorClass ("Illegal index value in variable initialization");
00308 
00309                                 initval = initVarVector[k]->value;
00310                                 if (osinstance->instanceData->variables->var[i]->ub == OSDBL_MAX)
00311                                 {       if (osinstance->instanceData->variables->var[i]->lb > initval)
00312                                                 throw ErrorClass ("Initial value outside of bounds");
00313                                 }
00314                                 else
00315                                         if (osinstance->instanceData->variables->var[i]->lb == -OSDBL_MAX)
00316                                         {       if (osinstance->instanceData->variables->var[i]->ub < initval)
00317                                                         throw ErrorClass ("Initial value outside of bounds");
00318                                         }
00319                                         else
00320                                         {       if ((osinstance->instanceData->variables->var[i]->lb > initval) ||
00321                                                         (osinstance->instanceData->variables->var[i]->ub < initval))
00322                                                         throw ErrorClass ("Initial value outside of bounds");
00323                                         }
00324 
00325                                 x[initVarVector[k]->idx] = initval;
00326                                 initialed[initVarVector[k]->idx] = true;
00327                         }
00328                 }
00329                 catch(const ErrorClass& eclass)
00330                 {       cout << "Error in BonminProblem::get_starting_point (OSBonminSolver.cpp)";
00331                         cout << endl << endl << endl;
00332                 }       
00333         }  //  end if (m1 > 0)          
00334 
00335         double default_initval;
00336         default_initval = 1.7171;
00337 
00338         for(k = 0; k < n1; k++)
00339         {
00340                 if (!initialed[k])
00341                         if (osinstance->instanceData->variables->var[k]->ub == OSDBL_MAX)
00342                                 if (osinstance->instanceData->variables->var[k]->lb <= default_initval)
00343                                         x[k] = default_initval;
00344                                 else
00345                                         x[k] = osinstance->instanceData->variables->var[k]->lb;
00346                         else
00347                                 if (osinstance->instanceData->variables->var[k]->lb == -OSDBL_MAX)
00348                                         if (osinstance->instanceData->variables->var[k]->ub >= default_initval)
00349                                                 x[k] = default_initval;
00350                                         else
00351                                                 x[k] = osinstance->instanceData->variables->var[k]->ub;
00352                                 else
00353                                         if ((osinstance->instanceData->variables->var[k]->lb <= default_initval) && 
00354                                                 (osinstance->instanceData->variables->var[k]->ub >= default_initval))
00355                                                 x[k] = default_initval;
00356                                         else
00357                                                 if (osinstance->instanceData->variables->var[k]->lb > default_initval)
00358                                                         x[k] = osinstance->instanceData->variables->var[k]->lb;
00359                                                 else
00360                                                         x[k] = osinstance->instanceData->variables->var[k]->ub;
00361         }
00362         for(i = 0; i < n1; i++){
00363                 std::cout << "INITIAL VALUE !!!!!!!!!!!!!!!!!!!!  " << x[ i] << std::endl;
00364         }
00365  
00366 
00367         delete[] initialed;
00368         
00369         return true;
00370 }//get_starting_point
00371 
00372 
00373 // returns the value of the objective function
00374 bool BonminProblem::eval_f(Index n, const Number* x, bool new_x, Number& obj_value){
00375         try{
00376                 obj_value  = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 )[ 0];
00377         }
00378         catch(const ErrorClass& eclass){
00379                 bonminErrorMsg = eclass.errormsg;
00380                 throw;  
00381         }
00382         if( CommonUtil::ISOSNAN( (double)obj_value) ) return false;
00383         return true;
00384 }
00385 
00386 bool BonminProblem::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f){
00387         int i;
00388         double *objGrad;
00389         try{
00390                 //objGrad = osinstance->calculateAllObjectiveFunctionGradients( const_cast<double*>(x), NULL, NULL,  new_x, 1)[ 0];
00391                 //std::cout << "Calculate Objective function gradient " << std::endl;
00392                 // we assume we are doing the objective function indexed by -1
00393                 objGrad = osinstance->calculateObjectiveFunctionGradient( const_cast<double*>(x), NULL, NULL, -1,  new_x, 1);
00394         }
00395         catch(const ErrorClass& eclass){
00396                 bonminErrorMsg = eclass.errormsg;
00397                 throw;  
00398         }
00399         for(i = 0; i < n; i++){
00400                 grad_f[ i]  = objGrad[ i];
00401                 //std::cout << grad_f[ i]  << std::endl;
00402         }
00403         std::cout << "DONE WITH Calculate Objective function gradient " << std::endl;
00404         return true;
00405 }//eval_grad_f
00406 
00407 // return the value of the constraints: g(x)
00408 bool BonminProblem::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) {
00409         try{
00410                 double *conVals = osinstance->calculateAllConstraintFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 );
00411                 int i;
00412                 for(i = 0; i < m; i++){
00413                         if( CommonUtil::ISOSNAN( (double)conVals[ i] ) ) return false;
00414                         g[i] = conVals[ i]  ;   
00415                 } 
00416                 return true;
00417         }
00418         catch(const ErrorClass& eclass){
00419                 bonminErrorMsg = eclass.errormsg;
00420                 throw;  
00421         }
00422 }//eval_g
00423 
00424 
00425 // return the structure or values of the jacobian
00426 bool BonminProblem::eval_jac_g(Index n, const Number* x, bool new_x,
00427                            Index m, Index nele_jac, Index* iRow, Index *jCol,
00428                            Number* values){
00429         SparseJacobianMatrix *sparseJacobian;
00430         if (values == NULL) {
00431                 // return the values of the jacobian of the constraints
00432                 //cout << "n: " << n << endl;
00433                 //cout << "m: " << m << endl;
00434                 //cout << "nele_jac: " <<  nele_jac << endl;
00435                 // return the structure of the jacobian
00436                 try{
00437                         sparseJacobian = osinstance->getJacobianSparsityPattern();
00438                 }
00439                 catch(const ErrorClass& eclass){
00440                         bonminErrorMsg =  eclass.errormsg; 
00441                         throw; 
00442                 }
00443                 int i = 0;
00444                 int k, idx;
00445                 for(idx = 0; idx < m; idx++){
00446                         for(k = *(sparseJacobian->starts + idx); k < *(sparseJacobian->starts + idx + 1); k++){
00447                                 iRow[i] = idx;
00448                                 jCol[i] = *(sparseJacobian->indexes + k);
00449                                 //cout << "ROW IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << iRow[i] << endl;
00450                                 //cout << "COL IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << jCol[i] << endl;
00451                                 i++;
00452                         }
00453                 }
00454         }
00455         else {
00456                 //std::cout << "EVALUATING JACOBIAN" << std::endl; 
00457                 try{
00458                         sparseJacobian = osinstance->calculateAllConstraintFunctionGradients( const_cast<double*>(x), NULL, NULL,  new_x, 1);
00459                 }
00460                 catch(const ErrorClass& eclass){
00461                         bonminErrorMsg = eclass.errormsg;
00462                         throw;  
00463                 }
00464                 //osinstance->getIterateResults( (double*)x, 0.0, NULL, -1, new_x,  1);
00465                 for(int i = 0; i < nele_jac; i++){
00466                         values[ i] = sparseJacobian->values[i];
00467                         //values[ i] = osinstance->m_mdJacValue[ i];
00468                         //cout << "values[i]:!!!!!!!!!!!!  " <<  values[ i] << endl;
00469                         //cout << "m_mdJacValue[ i]:!!!!!!!!!!!!  " <<  osinstance->m_mdJacValue[ i] << endl;
00470                 }
00471         }
00472   return true;
00473 }//eval_jac_g
00474 
00475 //return the structure or values of the hessian
00476 bool BonminProblem::eval_h(Index n, const Number* x, bool new_x,
00477                        Number obj_factor, Index m, const Number* lambda,
00478                        bool new_lambda, Index nele_hess, Index* iRow,
00479                        Index* jCol, Number* values){
00480 
00482         SparseHessianMatrix *sparseHessian;
00483         
00484         int i;
00485         if (values == NULL) {
00486                 // return the structure. This is a symmetric matrix, fill the lower left triangle only.
00487                 //cout << "get structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! "  << endl;
00488                 try{
00489                         sparseHessian = osinstance->getLagrangianHessianSparsityPattern( );
00490                 }
00491                 catch(const ErrorClass& eclass){
00492                         bonminErrorMsg = eclass.errormsg;
00493                         throw;  
00494                 }
00495                 //cout << "got structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! "  << endl;
00496                 for(i = 0; i < nele_hess; i++){
00497                         iRow[i] = *(sparseHessian->hessColIdx + i);
00498                         jCol[i] = *(sparseHessian->hessRowIdx + i);
00499                                 //cout << "ROW HESS IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << iRow[i] << endl;
00500                                 //cout << "COL HESS IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << jCol[i] << endl;
00501                 }
00502         }
00503         else {
00504                 //std::cout << "EVALUATING HESSIAN" << std::endl; 
00505                 // return the values. This is a symmetric matrix, fill the lower left triangle only
00506                 double* objMultipliers = new double[1];
00507                 objMultipliers[0] = obj_factor;
00508                 try{
00509                         sparseHessian = osinstance->calculateLagrangianHessian( const_cast<double*>(x), objMultipliers, const_cast<double*>(lambda) ,  new_x, 2);
00510                 delete[]  objMultipliers;
00511                 }
00512                 catch(const ErrorClass& eclass){
00513                         bonminErrorMsg = eclass.errormsg;
00514                         delete[]  objMultipliers;
00515                         throw;  
00516                 }
00517                 for(i = 0; i < nele_hess; i++){
00518                         values[ i]  = *(sparseHessian->hessValues + i);
00519                 }
00520         }
00522         return true;
00523 }//eval_h
00524 
00525 bool BonminProblem::get_scaling_parameters(Number& obj_scaling,
00526                         bool& use_x_scaling, Index n,
00527                         Number* x_scaling,
00528                     bool& use_g_scaling, Index m,
00529                     Number* g_scaling){
00530         /*
00531          * Bonmin assumes problems cast as a min
00532          * if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
00533                 obj_scaling = -1;
00534         }
00535     else obj_scaling = 1;
00536     use_x_scaling = false;
00537     use_g_scaling = false;
00538     */
00539         return true;
00540 }//get_scaling_parameters
00541 
00542 
00543 
00544 
00545 void
00546 BonminProblem::finalize_solution(TMINLP::SolverReturn status,
00547                             Index n, const Number* x, Number obj_value)
00548 {
00549         OSrLWriter *osrlwriter ;
00550         osrlwriter = new OSrLWriter();
00551         std::cout<<"Problem status: "<<status<<std::endl;
00552         std::cout<<"Objective value: "<<obj_value<<std::endl;
00553         if(printSol_ && x != NULL){
00554                 std::cout<<"Solution:"<<std::endl;
00555                 for(int i = 0 ; i < n ; i++){
00556                         std::cout<<"x["<<i<<"] = "<<x[i];
00557                         if(i < n-1) std::cout<<", ";
00558                 }
00559                 std::cout<<std::endl;
00560         }
00561         
00562           printf("\n\nObjective value\n");
00563           printf("f(x*) = %e\n", obj_value);
00564         int solIdx = 0;
00565         ostringstream outStr;
00566         double* mdObjValues = new double[1];
00567         std::string message = "Bonmin solver finishes to the end.";
00568         std::string solutionDescription = "";   
00569 
00570         try{
00571                 
00572                 // resultHeader information
00573                 if(osresult->setServiceName( "Bonmin solver service") != true)
00574                         throw ErrorClass("OSResult error: setServiceName");
00575                 if(osresult->setInstanceName(  osinstance->getInstanceName()) != true)
00576                         throw ErrorClass("OSResult error: setInstanceName");
00577 
00578                 //if(osresult->setJobID( osoption->jobID) != true)
00579                 //      throw ErrorClass("OSResult error: setJobID");
00580 
00581                 // set basic problem parameters
00582                 if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00583                         throw ErrorClass("OSResult error: setVariableNumer");
00584                 if(osresult->setObjectiveNumber( 1) != true)
00585                         throw ErrorClass("OSResult error: setObjectiveNumber");
00586                 if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00587                         throw ErrorClass("OSResult error: setConstraintNumber");
00588                 if(osresult->setSolutionNumber(  1) != true)
00589                         throw ErrorClass("OSResult error: setSolutionNumer");   
00590 
00591 
00592                 if(osresult->setGeneralMessage( message) != true)
00593                         throw ErrorClass("OSResult error: setGeneralMessage");
00594 
00595                 switch( status){
00596                         case SUCCESS:
00597                                 solutionDescription = "SUCCESS[BONMIN]: Algorithm terminated successfully at a locally optimal point, satisfying the convergence tolerances.";
00598                                 osresult->setSolutionStatus(solIdx,  "locallyOptimal", solutionDescription);
00599                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x), osinstance->getVariableNumber());
00600                                 //osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00601                                 mdObjValues[0] = obj_value;
00602                                 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
00603                         break;
00604                         case MAXITER_EXCEEDED:
00605                                 solutionDescription = "MAXITER_EXCEEDED[BONMIN]: Maximum number of iterations exceeded.";
00606                                 osresult->setSolutionStatus(solIdx,  "stoppedByLimit", solutionDescription);
00607                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x), osinstance->getVariableNumber());
00608                                 //osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00609                                 mdObjValues[0] = obj_value;
00610                                 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
00611                         break;
00612                         case STOP_AT_TINY_STEP:
00613                                 solutionDescription = "STOP_AT_TINY_STEP[BONMIN]: Algorithm proceeds with very little progress.";
00614                                 osresult->setSolutionStatus(solIdx,  "stoppedByLimit", solutionDescription);
00615                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>( x), osinstance->getVariableNumber());
00616                                 //osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00617                                 mdObjValues[0] = obj_value;
00618                                 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
00619                         break;
00620                         case STOP_AT_ACCEPTABLE_POINT:
00621                                 solutionDescription = "STOP_AT_ACCEPTABLE_POINT[BONMIN]: Algorithm stopped at a point that was converged, not to _desired_ tolerances, but to _acceptable_ tolerances";
00622                                 osresult->setSolutionStatus(solIdx,  "BonminAccetable", solutionDescription);
00623                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x), osinstance->getVariableNumber());
00624                                 //osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00625                                 mdObjValues[0] = obj_value;
00626                                 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
00627                         break;
00628                         case LOCAL_INFEASIBILITY:
00629                                 solutionDescription = "LOCAL_INFEASIBILITY[BONMIN]: Algorithm converged to a point of local infeasibility. Problem may be infeasible.";
00630                                 osresult->setSolutionStatus(solIdx,  "infeasible", solutionDescription);
00631                         break;
00632                         case USER_REQUESTED_STOP:
00633                                 solutionDescription = "USER_REQUESTED_STOP[BONMIN]: The user call-back function  intermediate_callback returned false, i.e., the user code requested a premature termination of the optimization.";
00634                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00635                         break;
00636                         case DIVERGING_ITERATES:
00637                                 solutionDescription = "DIVERGING_ITERATES[BONMIN]: It seems that the iterates diverge.";
00638                                 osresult->setSolutionStatus(solIdx,  "unbounded", solutionDescription);
00639                         break;
00640                         case RESTORATION_FAILURE:
00641                                 solutionDescription = "RESTORATION_FAILURE[BONMIN]: Restoration phase failed, algorithm doesn't know how to proceed.";
00642                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00643                         break;
00644                         case ERROR_IN_STEP_COMPUTATION:
00645                                 solutionDescription = "ERROR_IN_STEP_COMPUTATION[BONMIN]: An unrecoverable error occurred while IPOPT tried to compute the search direction.";
00646                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00647                         break;
00648                         case INVALID_NUMBER_DETECTED:
00649                                 solutionDescription = "INVALID_NUMcatBER_DETECTED[BONMIN]: Algorithm received an invalid number (such as NaN or Inf) from the NLP; see also option check_derivatives_for_naninf.";
00650                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00651                         break;
00652                         case INTERNAL_ERROR:
00653                                 solutionDescription = "INTERNAL_ERROR[BONMIN]: An unknown internal error occurred. Please contact the IPOPT authors through the mailing list.";
00654                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00655                         break;
00656                         default:
00657                                 solutionDescription = "OTHER[BONMIN]: other unknown solution status from Bonmin solver";
00658                                 osresult->setSolutionStatus(solIdx,  "other", solutionDescription);
00659                 }
00660                 osresult->setGeneralStatusType("normal");
00661                 delete osrlwriter;
00662                 delete[] mdObjValues;
00663                 osrlwriter = NULL;
00664 
00665         }
00666         
00667         catch(const ErrorClass& eclass){
00668                 osresult->setGeneralMessage( eclass.errormsg);
00669                 osresult->setGeneralStatusType( "error");
00670                 std::string osrl = osrlwriter->writeOSrL( osresult);
00671                 delete osrlwriter;
00672                 osrlwriter = NULL;
00673                 throw ErrorClass(  osrl) ;
00674                 delete[] mdObjValues;
00675                 mdObjValues = NULL;
00676         }
00677 }
00678 
00679 
00680    
00681 
00682 void BonminSolver::buildSolverInstance() throw (ErrorClass) {
00683         try{            
00684                 if(osil.length() == 0 && osinstance == NULL) throw ErrorClass("there is no instance");
00685                 if(osinstance == NULL){
00686                         m_osilreader = new OSiLReader();
00687                         osinstance = m_osilreader->readOSiL( osil);
00688                 }
00689                 // Create a new instance of your nlp 
00690                 tminlp = new BonminProblem( osinstance, osoption, osresult);
00691                 this->bCallbuildSolverInstance = true;
00692                 //Now initialize from tminlp
00693                 bonmin.initialize( GetRawPtr(tminlp) );
00694         }
00695         catch(const ErrorClass& eclass){
00696                 std::cout << "THERE IS AN ERROR" << std::endl;
00697                 osresult->setGeneralMessage( eclass.errormsg);
00698                 osresult->setGeneralStatusType( "error");
00699                 osrl = osrlwriter->writeOSrL( osresult);
00700                 throw ErrorClass( osrl) ;
00701         }                               
00702 }//end buildSolverInstance() 
00703 
00704 
00705 
00706 void BonminSolver::setSolverOptions() throw (ErrorClass) {
00707         try{
00708 
00709                 this->bSetSolverOptions = true;
00710                 bonmin.initializeOptionsAndJournalist();
00711                 //Register an additional option
00712                 bonmin.roptions()->AddStringOption2("print_solution","Do we print the solution or not?",
00713                                                  "yes",
00714                                                  "no", "No, we don't.",
00715                                                  "yes", "Yes, we do.",
00716                                                  "A longer comment can be put here");
00717                   
00718                  // Here we can change the default value of some Bonmin or Ipopt option
00719                 bonmin.options()->SetNumericValue("bonmin.time_limit", 1000); //changes bonmin's time limit
00720                 bonmin.options()->SetStringValue("mu_oracle","loqo");
00721                   
00722                 //Here we read several option files
00723                 bonmin.readOptionsFile("Mybonmin.opt");
00724                 bonmin.readOptionsFile();// This reads the default file "bonmin.opt"
00725                   
00726                 // Options can also be set by using a string with a format similar to the bonmin.opt file
00727                 bonmin.readOptionsString("bonmin.algorithm B-BB\n");
00728                   
00729                 // Now we can obtain the value of the new option
00730                 int printSolution;
00731                 bonmin.options()->GetEnumValue("print_solution", printSolution,"");
00732                 if(printSolution == 1){
00733                         tminlp->printSolutionAtEndOfAlgorithm();
00734                 }       
00735                 
00736                 if(osoption == NULL && osol.length() > 0)
00737                 {
00738                         m_osolreader = new OSoLReader();
00739                         osoption = m_osolreader->readOSoL( osol);
00740                 }
00741 
00742                 if(osoption != NULL){
00743                         std::cout << "number of solver options "  <<  osoption->getNumberOfSolverOptions() << std::endl;
00744                         std::vector<SolverOption*> optionsVector;
00745                         optionsVector = osoption->getSolverOptions( "bonmin");
00746                         char *pEnd;
00747                         int i;
00748                         int num_bonmin_options = optionsVector.size();
00749                         for(i = 0; i < num_bonmin_options; i++){
00750                                 std::cout << "bonmin solver option  "  << optionsVector[ i]->name << std::endl;
00751                                 if(optionsVector[ i]->type == "numeric" ){
00752                                         std::cout << "FOUND A NUMERIC OPTION  "  <<  os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) << std::endl;
00753                                         bonmin.options()->SetNumericValue(optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );      
00754                                 }
00755                                 else if(optionsVector[ i]->type == "integer" ){
00756                                         std::cout << "FOUND AN INTEGER OPTION  "  << atoi( optionsVector[ i]->value.c_str() ) << std::endl;
00757                                         bonmin.options()->SetIntegerValue(optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ) );
00758                                 }
00759                                 else if(optionsVector[ i]->type == "string" ){
00760                                         bonmin.options()->SetStringValue(optionsVector[ i]->name, optionsVector[ i]->value);
00761                                 }
00762                         }       
00763                 }
00764                 
00765         }
00766         
00767         catch(const ErrorClass& eclass){
00768                 std::cout << "THERE IS AN ERROR" << std::endl;
00769                 osresult->setGeneralMessage( eclass.errormsg);
00770                 osresult->setGeneralStatusType( "error");
00771                 osrl = osrlwriter->writeOSrL( osresult);
00772                 throw ErrorClass( osrl) ;
00773         }                               
00774 }//end setSolverOptions() 
00775 
00776 
00777 //void BonminSolver::solve() throw (ErrorClass) {
00778 void BonminSolver::solve() throw (ErrorClass) {
00779         if( this->bCallbuildSolverInstance == false) buildSolverInstance();
00780         std::cout << "set Solver Options for Gus" << std::endl;
00781         if( this->bSetSolverOptions == false) setSolverOptions();
00782         std::cout << "done setting set Solver Options for Gus" << std::endl;
00783 
00784         try{
00785                 double start = CoinCpuTime();
00786                 //OSiLWriter osilwriter;
00787                 //cout << osilwriter.writeOSiL( osinstance) << endl;
00788                 if(osinstance->getVariableNumber() <= 0)throw ErrorClass("Bonmin requires decision variables");
00789 
00790 //              The next line should come out eventually...
00791                 if(osinstance->getConstraintNumber() <= 0)throw ErrorClass("Bonmin cannot handle unconstrained problems");
00792 
00793                 double duration = CoinCpuTime() - start;
00794                 cout << "Parsing took (seconds): "<< duration << endl;
00795 
00796                   try {
00797                     Bab bb;
00798                     bb(  bonmin);  //process parameter file using Ipopt and do branch and bound using Cbc
00799 
00800 
00801                   }
00802                   catch(TNLPSolver::UnsolvedError *E) {
00803                     //There has been a failure to solve a problem with Ipopt.
00804                     std::cerr<<"Ipopt has failed to solve a problem"<<std::endl;
00805                   }
00806                   catch(OsiTMINLPInterface::SimpleError &E) {
00807                     std::cerr<<E.className()<<"::"<<E.methodName()
00808                              <<std::endl
00809                              <<E.message()<<std::endl;
00810                   }
00811                   catch(CoinError &E) {
00812                     std::cerr<<E.className()<<"::"<<E.methodName()
00813                              <<std::endl
00814                              <<E.message()<<std::endl;
00815                   }
00816                 
00817 
00818                 // see if we have a linear program
00819                 if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Bonmin NEEDS AN OBJECTIVE FUNCTION");
00820                 if( (osinstance->getNumberOfNonlinearExpressions() == 0) && (osinstance->getNumberOfQuadraticTerms() == 0) ) 
00821                         //app->Options()->SetStringValue("hessian_approximation", "limited-memory");
00822                 // if it is a max problem call scaling and set to -1
00823                 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
00824                         //app->Options()->SetStringValue("nlp_scaling_method", "user-scaling");
00825                 }
00826 
00827                 std::cout << "Finish Bonmin Optimize" << std::endl;
00828                 osrl = osrlwriter->writeOSrL( osresult);
00829                 std::cout << "Finish writing the osrl" << std::endl;
00830 
00831 
00832         }
00833         catch(const ErrorClass& eclass){
00834                 osresult->setGeneralMessage( eclass.errormsg);
00835                 osresult->setGeneralStatusType( "error");
00836                 osrl = osrlwriter->writeOSrL( osresult);
00837                 throw ErrorClass( osrl) ;
00838         }
00839 }//solve
00840 
00841 
00842 void BonminSolver::dataEchoCheck(){
00843         
00844         int i;
00845         
00846         // print out problem parameters
00847         cout << "This is problem:  " << osinstance->getInstanceName() << endl;
00848         cout << "The problem source is:  " << osinstance->getInstanceSource() << endl;
00849         cout << "The problem description is:  " << osinstance->getInstanceDescription() << endl;
00850         cout << "number of variables = " << osinstance->getVariableNumber() << endl;
00851         cout << "number of Rows = " << osinstance->getConstraintNumber() << endl;
00852 
00853         // print out the variable information
00854         if(osinstance->getVariableNumber() > 0){
00855                 for(i = 0; i < osinstance->getVariableNumber(); i++){
00856                         if(osinstance->getVariableNames() != NULL) cout << "variable Names  " << osinstance->getVariableNames()[ i]  << endl;
00857                         if(osinstance->getVariableTypes() != NULL) cout << "variable Types  " << osinstance->getVariableTypes()[ i]  << endl;
00858                         if(osinstance->getVariableLowerBounds() != NULL) cout << "variable Lower Bounds  " << osinstance->getVariableLowerBounds()[ i]  << endl;
00859                         if(osinstance->getVariableUpperBounds() != NULL) cout << "variable Upper Bounds  " <<  osinstance->getVariableUpperBounds()[i] << endl;
00860                 }
00861         }
00862         
00863         // print out objective function information
00864         if(osinstance->getVariableNumber() > 0 || osinstance->instanceData->objectives->obj != NULL || osinstance->instanceData->objectives->numberOfObjectives > 0){
00865                 if( osinstance->getObjectiveMaxOrMins()[0] == "min")  cout <<  "problem is a minimization" << endl;
00866                 else cout <<  "problem is a maximization" << endl;
00867                 for(i = 0; i < osinstance->getVariableNumber(); i++){
00868                         cout << "OBJ COEFFICIENT =  " <<  osinstance->getDenseObjectiveCoefficients()[0][i] << endl;
00869                 }
00870         }
00871         // print out constraint information
00872         if(osinstance->getConstraintNumber() > 0){
00873                 for(i = 0; i < osinstance->getConstraintNumber(); i++){
00874                         if(osinstance->getConstraintNames() != NULL) cout << "row name = " << osinstance->getConstraintNames()[i] <<  endl;
00875                         if(osinstance->getConstraintLowerBounds() != NULL) cout << "row lower bound = " << osinstance->getConstraintLowerBounds()[i] <<  endl;
00876                         if(osinstance->getConstraintUpperBounds() != NULL) cout << "row upper bound = " << osinstance->getConstraintUpperBounds()[i] <<  endl; 
00877                 }
00878         }
00879         
00880         // print out linear constraint data
00881         cout << endl;
00882         cout << "number of nonzeros =  " << osinstance->getLinearConstraintCoefficientNumber() << endl;
00883         for(i = 0; i <= osinstance->getVariableNumber(); i++){
00884                 cout << "Start Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->starts[ i] << endl;
00885         }
00886         cout << endl;
00887         for(i = 0; i < osinstance->getLinearConstraintCoefficientNumber(); i++){
00888                 cout << "Index Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->indexes[i] << endl;
00889                 cout << "Nonzero Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->values[i] << endl;
00890         }
00891 
00892         // print out quadratic data
00893         cout << "number of qterms =  " <<  osinstance->getNumberOfQuadraticTerms() << endl;
00894         for(int i = 0; i <  osinstance->getNumberOfQuadraticTerms(); i++){
00895                 cout << "Row Index = " <<  osinstance->getQuadraticTerms()->rowIndexes[i] << endl;
00896                 cout << "Var Index 1 = " << osinstance->getQuadraticTerms()->varOneIndexes[ i] << endl;
00897                 cout << "Var Index 2 = " << osinstance->getQuadraticTerms()->varTwoIndexes[ i] << endl;
00898                 cout << "Coefficient = " << osinstance->getQuadraticTerms()->coefficients[ i] << endl;
00899         }
00900 } // end dataEchoCheck
00901 
00902 
00903 BonminProblem::BonminProblem(OSInstance *osinstance_,  OSOption *osoption_, OSResult *osresult_) {
00904         osinstance = osinstance_;
00905         osoption = osoption_;
00906         osresult = osresult_;
00907         printSol_ = false;
00908 }
00909 
00910 BonminProblem::~BonminProblem() {
00911 
00912 }
00913 
00914 
00915 

Generated on Fri Jan 7 03:24:42 2011 by  doxygen 1.4.7