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

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

Generated on Tue Sep 30 03:01:25 2008 by  doxygen 1.4.7