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

Go to the documentation of this file.
00001 
00020 #include "OSIpoptSolver.h"
00021 #include "OSCommonUtil.h"
00022 
00023 #include "OSDataStructures.h"
00024 #include "OSParameters.h" 
00025 #include "OSCommonUtil.h"
00026 #include "OSMathUtil.h"
00027 
00028 
00029 using std::cout; 
00030 using std::endl; 
00031 using std::ostringstream;
00032 using namespace Ipopt;
00033 
00034 
00035 IpoptSolver::IpoptSolver() {
00036         osrlwriter = new OSrLWriter();
00037         osresult = new OSResult();
00038         m_osilreader = NULL;
00039         ipoptErrorMsg = "";
00040 
00041 } 
00042 
00043 IpoptSolver::~IpoptSolver() {
00044         #ifdef DEBUG
00045         cout << "inside IpoptSolver destructor" << endl;
00046         #endif
00047         if(m_osilreader != NULL) delete m_osilreader;
00048         m_osilreader = NULL;
00049         delete osresult;
00050         osresult = NULL;
00051         delete osrlwriter;
00052         osrlwriter = NULL;
00053         //delete osinstance;
00054         //osinstance = NULL;
00055         #ifdef DEBUG
00056         cout << "leaving IpoptSolver destructor" << endl;
00057         #endif
00058 }
00059 
00060 // returns the size of the problem
00061 bool IpoptProblem::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00062                              Index& nnz_h_lag, IndexStyleEnum& index_style)
00063 {
00064         if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Ipopt NEEDS AN OBJECTIVE FUNCTION");                        
00065         // number of variables
00066         n = osinstance->getVariableNumber();
00067         // number of constraints
00068         m = osinstance->getConstraintNumber();
00069         cout << "number variables  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << n << endl;
00070         cout << "number constraints  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << m << endl;
00071         try{
00072                 osinstance->initForAlgDiff( );
00073         }
00074         catch(const ErrorClass& eclass){
00075                 ipoptErrorMsg = eclass.errormsg;
00076                 throw;  
00077         }       
00078         // use the OS Expression tree for function evaluations instead of CppAD
00079         osinstance->bUseExpTreeForFunEval = true;
00080         //std::cout << "Call sparse jacobian" << std::endl;
00081         SparseJacobianMatrix *sparseJacobian = NULL;
00082         try{
00083                 sparseJacobian = osinstance->getJacobianSparsityPattern();
00084         }
00085         catch(const ErrorClass& eclass){
00086                 ipoptErrorMsg = eclass.errormsg;
00087                 throw;  
00088         }
00089         //std::cout << "Done calling sparse jacobian" << std::endl;
00090         nnz_jac_g = sparseJacobian->valueSize;
00091         cout << "nnz_jac_g  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_jac_g << endl;  
00092         // nonzeros in upper hessian
00093         
00094         if( (osinstance->getNumberOfNonlinearExpressions() == 0) && (osinstance->getNumberOfQuadraticTerms() == 0) ) {
00095                 cout << "This is a linear program"  << endl;
00096                 nnz_h_lag = 0;
00097         }
00098         else{
00099                 //std::cout << "Get Lagrangain Hessian Sparsity Pattern " << std::endl;
00100                 SparseHessianMatrix *sparseHessian = osinstance->getLagrangianHessianSparsityPattern();
00101                 //std::cout << "Done Getting Lagrangain Hessian Sparsity Pattern " << std::endl;
00102                 nnz_h_lag = sparseHessian->hessDimension;
00103         }
00104         cout << "nnz_h_lag  !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_h_lag << endl;  
00105         // use the C style indexing (0-based)
00106         index_style = TNLP::C_STYLE;
00107   
00109 
00110   return true;
00111 }//get_nlp_info
00112 
00113 
00114 bool  IpoptProblem::get_bounds_info(Index n, Number* x_l, Number* x_u,
00115                                 Index m, Number* g_l, Number* g_u){
00116         int i; 
00117         double * mdVarLB = osinstance->getVariableLowerBounds();
00118         //std::cout << "GET BOUNDS INFORMATION FOR IPOPT !!!!!!!!!!!!!!!!!" << std::endl;
00119         // variables upper bounds
00120         double * mdVarUB = osinstance->getVariableUpperBounds();
00121 
00122         for(i = 0; i < n; i++){
00123                 x_l[ i] = mdVarLB[ i];
00124                 x_u[ i] = mdVarUB[ i];
00125                 //cout << "x_l !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_l[i] << endl;
00126                 //cout << "x_u !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_u[i] << endl;
00127         }
00128         // Ipopt interprets any number greater than nlp_upper_bound_inf as
00129         // infinity. The default value of nlp_upper_bound_inf and nlp_lower_bound_inf
00130         // is 1e19 and can be changed through ipopt options.
00131         // e.g. g_u[0] = 2e19;
00132 
00133         //constraint lower bounds
00134         double * mdConLB = osinstance->getConstraintLowerBounds();
00135         //constraint upper bounds
00136         double * mdConUB = osinstance->getConstraintUpperBounds();
00137         
00138         for(int i = 0; i < m; i++){
00139                 g_l[ i] = mdConLB[ i];
00140                 g_u[ i] = mdConUB[ i];
00141                 //cout << "lower !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_l[i] << endl;
00142                 //cout << "upper !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_u[i] << endl;
00143         }  
00144         return true;
00145 }//get_bounds_info
00146 
00147 
00148 // returns the initial point for the problem
00149 bool IpoptProblem::get_starting_point(Index n, bool init_x, Number* x,
00150      bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda,
00151      Number* lambda) {
00152         // Here, we assume we only have starting values for x, if you code
00153         // your own NLP, you can provide starting values for the dual variables
00154         // if you wish
00155         assert(init_x == true);
00156         assert(init_z == false);
00157         assert(init_lambda == false);
00158         int i;
00159         //cout << "get initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00160         double *mdXInit = osinstance->getVariableInitialValues(); 
00161         if( mdXInit != NULL) {
00162                 for(i = 0; i < n; i++){
00163                         if( CommonUtil::ISOSNAN( mdXInit[ i]) == true){ 
00164                                 x[ i] = 1.7171; 
00165                                 //std::cout << "INITIAL VALUE !!!!!!!!!!!!!!!!!!!!  " << x[ i] << std::endl;
00166                         }
00167                         else x[ i] = mdXInit[ i];
00168                         //std::cout << "INITIAL VALUE !!!!!!!!!!!!!!!!!!!!  " << x[ i] << std::endl;    
00169                 }       
00170         }
00171         else{
00172                 for(i = 0; i < n; i++){
00173                         x[ i] = 1.7171;
00174                 }
00175         }
00176         //cout << "got initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00177         return true;
00178 }//get_starting_point
00179 
00180 // returns the value of the objective function
00181 bool IpoptProblem::eval_f(Index n, const Number* x, bool new_x, Number& obj_value){
00182         try{
00183                 obj_value  = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 )[ 0];
00184         }
00185         catch(const ErrorClass& eclass){
00186                 ipoptErrorMsg = eclass.errormsg;
00187                 throw;  
00188         }
00189         if( CommonUtil::ISOSNAN( (double)obj_value) ) return false;
00190         return true;
00191 }
00192 
00193 bool IpoptProblem::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f){
00194         int i;
00195         double *objGrad;
00196         try{
00197                 //objGrad = osinstance->calculateAllObjectiveFunctionGradients( const_cast<double*>(x), NULL, NULL,  new_x, 1)[ 0];
00198                 objGrad = osinstance->calculateObjectiveFunctionGradient( const_cast<double*>(x), NULL, NULL, -1,  new_x, 1);
00199         }
00200         catch(const ErrorClass& eclass){
00201                 ipoptErrorMsg = eclass.errormsg;
00202                 throw;  
00203         }
00204         for(i = 0; i < n; i++){
00205                 grad_f[ i]  = objGrad[ i];
00206         }
00207         return true;
00208 }//eval_grad_f
00209 
00210 // return the value of the constraints: g(x)
00211 bool IpoptProblem::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) {
00212         try{
00213                 double *conVals = osinstance->calculateAllConstraintFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 );
00214                 int i;
00215                 for(i = 0; i < m; i++){
00216                         if( CommonUtil::ISOSNAN( (double)conVals[ i] ) ) return false;
00217                         g[i] = conVals[ i]  ;   
00218                 } 
00219                 return true;
00220         }
00221         catch(const ErrorClass& eclass){
00222                 ipoptErrorMsg = eclass.errormsg;
00223                 throw;  
00224         }
00225 }//eval_g
00226 
00227 
00228 // return the structure or values of the jacobian
00229 bool IpoptProblem::eval_jac_g(Index n, const Number* x, bool new_x,
00230                            Index m, Index nele_jac, Index* iRow, Index *jCol,
00231                            Number* values){
00232         SparseJacobianMatrix *sparseJacobian;
00233         if (values == NULL) {
00234                 // return the values of the jacobian of the constraints
00235                 //cout << "n: " << n << endl;
00236                 //cout << "m: " << m << endl;
00237                 //cout << "nele_jac: " <<  nele_jac << endl;
00238                 // return the structure of the jacobian
00239                 try{
00240                         sparseJacobian = osinstance->getJacobianSparsityPattern();
00241                 }
00242                 catch(const ErrorClass& eclass){
00243                         ipoptErrorMsg =  eclass.errormsg; 
00244                         throw; 
00245                 }
00246                 int i = 0;
00247                 int k, idx;
00248                 for(idx = 0; idx < m; idx++){
00249                         for(k = *(sparseJacobian->starts + idx); k < *(sparseJacobian->starts + idx + 1); k++){
00250                                 iRow[i] = idx;
00251                                 jCol[i] = *(sparseJacobian->indexes + k);
00252                                 //cout << "ROW IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << iRow[i] << endl;
00253                                 //cout << "COL IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << jCol[i] << endl;
00254                                 i++;
00255                         }
00256                 }
00257         }
00258         else {
00259                 //std::cout << "EVALUATING JACOBIAN" << std::endl; 
00260                 try{
00261                         sparseJacobian = osinstance->calculateAllConstraintFunctionGradients( const_cast<double*>(x), NULL, NULL,  new_x, 1);
00262                 }
00263                 catch(const ErrorClass& eclass){
00264                         ipoptErrorMsg = eclass.errormsg;
00265                         throw;  
00266                 }
00267                 //osinstance->getIterateResults( (double*)x, 0.0, NULL, -1, new_x,  1);
00268                 for(int i = 0; i < nele_jac; i++){
00269                         values[ i] = sparseJacobian->values[i];
00270                         //values[ i] = osinstance->m_mdJacValue[ i];
00271                         //cout << "values[i]:!!!!!!!!!!!!  " <<  values[ i] << endl;
00272                         //cout << "m_mdJacValue[ i]:!!!!!!!!!!!!  " <<  osinstance->m_mdJacValue[ i] << endl;
00273                 }
00274         }
00275   return true;
00276 }//eval_jac_g
00277 
00278 //return the structure or values of the hessian
00279 bool IpoptProblem::eval_h(Index n, const Number* x, bool new_x,
00280                        Number obj_factor, Index m, const Number* lambda,
00281                        bool new_lambda, Index nele_hess, Index* iRow,
00282                        Index* jCol, Number* values){
00283 
00285         SparseHessianMatrix *sparseHessian;
00286         
00287         int i;
00288         if (values == NULL) {
00289                 // return the structure. This is a symmetric matrix, fill the lower left triangle only.
00290                 //cout << "get structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! "  << endl;
00291                 try{
00292                         sparseHessian = osinstance->getLagrangianHessianSparsityPattern( );
00293                 }
00294                 catch(const ErrorClass& eclass){
00295                         ipoptErrorMsg = eclass.errormsg;
00296                         throw;  
00297                 }
00298                 //cout << "got structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! "  << endl;
00299                 for(i = 0; i < nele_hess; i++){
00300                         iRow[i] = *(sparseHessian->hessColIdx + i);
00301                         jCol[i] = *(sparseHessian->hessRowIdx + i);
00302                                 //cout << "ROW HESS IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << iRow[i] << endl;
00303                                 //cout << "COL HESS IDX  !!!!!!!!!!!!!!!!!!!!!!!!!!!"  << jCol[i] << endl;
00304                 }
00305         }
00306         else {
00307                 //std::cout << "EVALUATING HESSIAN" << std::endl; 
00308                 // return the values. This is a symmetric matrix, fill the lower left triangle only
00309                 double* objMultipliers = new double[1];
00310                 objMultipliers[0] = obj_factor;
00311                 try{
00312                         sparseHessian = osinstance->calculateLagrangianHessian( const_cast<double*>(x), objMultipliers, (double*)lambda ,  new_x, 2);
00313                 delete[]  objMultipliers;
00314                 }
00315                 catch(const ErrorClass& eclass){
00316                         ipoptErrorMsg = eclass.errormsg;
00317                         delete[]  objMultipliers;
00318                         throw;  
00319                 }
00320                 for(i = 0; i < nele_hess; i++){
00321                         values[ i]  = *(sparseHessian->hessValues + i);
00322                 }
00323         }
00325         return true;
00326 }//eval_h
00327 
00328 bool IpoptProblem::get_scaling_parameters(Number& obj_scaling,
00329                         bool& use_x_scaling, Index n,
00330                         Number* x_scaling,
00331                     bool& use_g_scaling, Index m,
00332                     Number* g_scaling){
00333         if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
00334                 obj_scaling = -1;
00335         }
00336     else obj_scaling = 1;
00337     use_x_scaling = false;
00338     use_g_scaling = false;
00339         return true;
00340 }//get_scaling_parameters
00341 
00342 void IpoptProblem::finalize_solution(SolverReturn status,
00343                                Index n, const Number* x, const Number* z_L, const Number* z_U,
00344                                   Index m, const Number* g, const Number* lambda,
00345                                   Number obj_value,
00346                                    const IpoptData* ip_data,
00347                                    IpoptCalculatedQuantities* ip_cq)
00348 {
00349           // here is where we would store the solution to variables, or write to a file, etc
00350           // so we could use the solution.
00351           // For this example, we write the solution to the console
00352         OSrLWriter *osrlwriter ;
00353         osrlwriter = new OSrLWriter();
00354 #ifdef DEBUG
00355           printf("\n\nSolution of the primal variables, x\n");
00356           for (Index i=0; i<n; i++) {
00357             printf("x[%d] = %e\n", i, x[i]);
00358           }
00359         
00360           printf("\n\nSolution of the bound multipliers, z_L and z_U\n");
00361           for (Index i=0; i<n; i++) {
00362             printf("z_L[%d] = %e\n", i, z_L[i]);
00363           }
00364           for (Index i=0; i<n; i++) {
00365             printf("z_U[%d] = %e\n", i, z_U[i]);
00366           }
00367 #endif
00368           printf("\n\nObjective value\n");
00369           printf("f(x*) = %e\n", obj_value);
00370           
00371         int solIdx = 0;
00372         int numberOfOtherVariableResult;
00373         int otherIdx;
00374         ostringstream outStr;
00375         
00376         std::string *rcost = NULL;
00377         double* mdObjValues = new double[1];
00378         std::string message = "Ipopt solver finishes to the end.";
00379         std::string solutionDescription = "";   
00380 
00381         try{
00382                 // resultHeader infomration
00383                 if(osresult->setServiceName( "Ipopt solver service") != true)
00384                         throw ErrorClass("OSResult error: setServiceName");
00385                 if(osresult->setInstanceName(  osinstance->getInstanceName()) != true)
00386                         throw ErrorClass("OSResult error: setInstanceName");
00387 
00388                 //if(osresult->setJobID( osoption->jobID) != true)
00389                 //      throw ErrorClass("OSResult error: setJobID");
00390 
00391                 // set basic problem parameters
00392                 if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00393                         throw ErrorClass("OSResult error: setVariableNumer");
00394                 if(osresult->setObjectiveNumber( 1) != true)
00395                         throw ErrorClass("OSResult error: setObjectiveNumber");
00396                 if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00397                         throw ErrorClass("OSResult error: setConstraintNumber");
00398                 if(osresult->setSolutionNumber(  1) != true)
00399                         throw ErrorClass("OSResult error: setSolutionNumer");   
00400 
00401 
00402                 if(osresult->setGeneralMessage( message) != true)
00403                         throw ErrorClass("OSResult error: setGeneralMessage");
00404 
00405                 switch( status){
00406                         case SUCCESS:
00407                                 solutionDescription = "SUCCESS[IPOPT]: Algorithm terminated successfully at a locally optimal point, satisfying the convergence tolerances.";
00408                                 osresult->setSolutionStatus(solIdx,  "locallyOptimal", solutionDescription);
00409                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x));
00410                                 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00411                                 mdObjValues[0] = obj_value;
00412                                 osresult->setObjectiveValues(solIdx, mdObjValues);
00413                                 
00414                                 
00415                                 // set other
00416                                 
00417                                 numberOfOtherVariableResult = 2;
00418                                 osresult->setNumberOfOtherVariableResult(solIdx, numberOfOtherVariableResult);
00419                                 
00420                         
00421                                 rcost = new std::string[ osinstance->getVariableNumber()];
00422                                 
00423 
00424                                 for (Index i = 0; i < n; i++) {
00425                                     rcost[ i] =  os_dtoa_format( z_L[i]); 
00426                                  }
00427                                 otherIdx = 0;
00428                                 osresult->setAnOtherVariableResult(solIdx, otherIdx, "varL", "Lagrange Multiplier on the Variable Lower Bound", rcost);
00429                                 
00430                                 
00431                                 for (Index i = 0; i < n; i++) {
00432                                     rcost[ i] =  os_dtoa_format( z_U[i]); 
00433                                  }
00434                                 otherIdx = 1;
00435                                 osresult->setAnOtherVariableResult(solIdx, otherIdx, "varU", "Lagrange Multiplier on the Variable Upper Bound", rcost);
00436                                 
00437                                 delete[] rcost;
00438                                 
00439                                 // done with other
00440                                 
00441                                 
00442                         break;
00443                         case MAXITER_EXCEEDED:
00444                                 solutionDescription = "MAXITER_EXCEEDED[IPOPT]: Maximum number of iterations exceeded.";
00445                                 osresult->setSolutionStatus(solIdx,  "stoppedByLimit", solutionDescription);
00446                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x));
00447                                 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00448                                 mdObjValues[0] = obj_value;
00449                                 osresult->setObjectiveValues(solIdx, mdObjValues);
00450                         break;
00451                         case STOP_AT_TINY_STEP:
00452                                 solutionDescription = "STOP_AT_TINY_STEP[IPOPT]: Algorithm proceeds with very little progress.";
00453                                 osresult->setSolutionStatus(solIdx,  "stoppedByLimit", solutionDescription);
00454                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>( x));
00455                                 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00456                                 mdObjValues[0] = obj_value;
00457                                 osresult->setObjectiveValues(solIdx, mdObjValues);
00458                         break;
00459                         case STOP_AT_ACCEPTABLE_POINT:
00460                                 solutionDescription = "STOP_AT_ACCEPTABLE_POINT[IPOPT]: Algorithm stopped at a point that was converged, not to _desired_ tolerances, but to _acceptable_ tolerances";
00461                                 osresult->setSolutionStatus(solIdx,  "IpoptAccetable", solutionDescription);
00462                                 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x));
00463                                 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda));
00464                                 mdObjValues[0] = obj_value;
00465                                 osresult->setObjectiveValues(solIdx, mdObjValues);
00466                         break;
00467                         case LOCAL_INFEASIBILITY:
00468                                 solutionDescription = "LOCAL_INFEASIBILITY[IPOPT]: Algorithm converged to a point of local infeasibility. Problem may be infeasible.";
00469                                 osresult->setSolutionStatus(solIdx,  "infeasible", solutionDescription);
00470                         break;
00471                         case USER_REQUESTED_STOP:
00472                                 solutionDescription = "USER_REQUESTED_STOP[IPOPT]: The user call-back function  intermediate_callback returned false, i.e., the user code requested a premature termination of the optimization.";
00473                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00474                         break;
00475                         case DIVERGING_ITERATES:
00476                                 solutionDescription = "DIVERGING_ITERATES[IPOPT]: It seems that the iterates diverge.";
00477                                 osresult->setSolutionStatus(solIdx,  "unbounded", solutionDescription);
00478                         break;
00479                         case RESTORATION_FAILURE:
00480                                 solutionDescription = "RESTORATION_FAILURE[IPOPT]: Restoration phase failed, algorithm doesn't know how to proceed.";
00481                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00482                         break;
00483                         case ERROR_IN_STEP_COMPUTATION:
00484                                 solutionDescription = "ERROR_IN_STEP_COMPUTATION[IPOPT]: An unrecoverable error occurred while IPOPT tried to compute the search direction.";
00485                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00486                         break;
00487                         case INVALID_NUMBER_DETECTED:
00488                                 solutionDescription = "INVALID_NUMcatBER_DETECTED[IPOPT]: Algorithm received an invalid number (such as NaN or Inf) from the NLP; see also option check_derivatives_for_naninf.";
00489                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00490                         break;
00491                         case INTERNAL_ERROR:
00492                                 solutionDescription = "INTERNAL_ERROR[IPOPT]: An unknown internal error occurred. Please contact the IPOPT authors through the mailing list.";
00493                                 osresult->setSolutionStatus(solIdx,  "error", solutionDescription);
00494                         break;
00495                         default:
00496                                 solutionDescription = "OTHER[IPOPT]: other unknown solution status from Ipopt solver";
00497                                 osresult->setSolutionStatus(solIdx,  "other", solutionDescription);
00498                 }
00499                 osresult->setGeneralStatusType("success");
00500                 delete osrlwriter;
00501                 delete[] mdObjValues;
00502                 osrlwriter = NULL;
00503 
00504         }
00505         catch(const ErrorClass& eclass){
00506                 osresult->setGeneralMessage( eclass.errormsg);
00507                 osresult->setGeneralStatusType( "error");
00508                 std::string osrl = osrlwriter->writeOSrL( osresult);
00509                 delete osrlwriter;
00510                 osrlwriter = NULL;
00511                 throw ErrorClass(  osrl) ;
00512                 delete[] mdObjValues;
00513                 mdObjValues = NULL;
00514         }
00516 }
00517 
00518 
00519 void IpoptSolver::buildSolverInstance() throw (ErrorClass) {
00520         try{
00521                 
00522                 if(osil.length() == 0 && osinstance == NULL) throw ErrorClass("there is no instance");
00523                 if(osinstance == NULL){
00524                         m_osilreader = new OSiLReader();
00525                         osinstance = m_osilreader->readOSiL( osil);
00526                 }
00527                 // Create a new instance of your nlp 
00528                 nlp = new IpoptProblem( osinstance, osresult);
00529                 app = new IpoptApplication();
00530                 this->bCallbuildSolverInstance = true;
00531         }
00532         catch(const ErrorClass& eclass){
00533                 std::cout << "THERE IS AN ERROR" << std::endl;
00534                 osresult->setGeneralMessage( eclass.errormsg);
00535                 osresult->setGeneralStatusType( "error");
00536                 osrl = osrlwriter->writeOSrL( osresult);
00537                 throw ErrorClass( osrl) ;
00538         }                               
00539 }//end buildSolverInstance() 
00540 
00541 
00542 //void IpoptSolver::solve() throw (ErrorClass) {
00543 void IpoptSolver::solve() throw (ErrorClass) {
00544         if( this->bCallbuildSolverInstance == false) buildSolverInstance();
00545         try{
00546                 clock_t start, finish;
00547                 double duration;
00548                 start = clock();
00549                 //OSiLWriter osilwriter;
00550                 //cout << osilwriter.writeOSiL( osinstance) << endl;
00551                 if(osinstance->getVariableNumber() <= 0)throw ErrorClass("Ipopt requires decision variables");
00552                 finish = clock();
00553                 duration = (double) (finish - start) / CLOCKS_PER_SEC;
00554                 cout << "Parsing took (seconds): "<< duration << endl;
00555                 //dataEchoCheck();
00556                 /***************now the ipopt invokation*********************/
00557                 app->Options()->SetNumericValue("tol", 1e-9);
00558                 app->Options()->SetIntegerValue("print_level", 0);
00559                 app->Options()->SetIntegerValue("max_iter", 20000);
00560                 app->Options()->SetStringValue("mu_strategy", "adaptive");
00561                 app->Options()->SetStringValue("output_file", "ipopt.out");
00562                 app->Options()->SetStringValue("check_derivatives_for_naninf", "yes");
00563                 // see if we have a linear program
00564                 if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Ipopt NEEDS AN OBJECTIVE FUNCTION");
00565                 if( (osinstance->getNumberOfNonlinearExpressions() == 0) && (osinstance->getNumberOfQuadraticTerms() == 0) ) 
00566                         app->Options()->SetStringValue("hessian_approximation", "limited-memory");
00567                 // if it is a max problem call scaling and set to -1
00568                 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
00569                         app->Options()->SetStringValue("nlp_scaling_method", "user-scaling");
00570                 }
00571                 // Intialize the IpoptApplication and process the options
00572                 std::cout << "Call Ipopt Initialize" << std::endl;
00573                 app->Initialize();
00574                 std::cout << "Finished Ipopt Initialize" << std::endl;
00575                 //nlp->osinstance = this->osinstance;
00576                 // Ask Ipopt to solve the problem
00577                 std::cout << "Call Ipopt Optimize" << std::endl;
00578                 ApplicationReturnStatus status = app->OptimizeTNLP( nlp);
00579                 std::cout << "Finish Ipopt Optimize" << std::endl;
00580                 osrl = osrlwriter->writeOSrL( osresult);
00581                 std::cout << "Finish writing the osrl" << std::endl;
00582                 //if (status != Solve_Succeeded) {
00583                 if (status < -2) {
00584                         throw ErrorClass("Ipopt FAILED TO SOLVE THE PROBLEM: " + ipoptErrorMsg);
00585                 }       
00586         }
00587         catch(const ErrorClass& eclass){
00588                 osresult->setGeneralMessage( eclass.errormsg);
00589                 osresult->setGeneralStatusType( "error");
00590                 osrl = osrlwriter->writeOSrL( osresult);
00591                 throw ErrorClass( osrl) ;
00592         }
00593 }//solve
00594 
00595 
00596 void IpoptSolver::dataEchoCheck(){
00597         
00598         int i;
00599         
00600         // print out problem parameters
00601         cout << "This is problem:  " << osinstance->getInstanceName() << endl;
00602         cout << "The problem source is:  " << osinstance->getInstanceSource() << endl;
00603         cout << "The problem description is:  " << osinstance->getInstanceDescription() << endl;
00604         cout << "number of variables = " << osinstance->getVariableNumber() << endl;
00605         cout << "number of Rows = " << osinstance->getConstraintNumber() << endl;
00606 
00607         // print out the variable information
00608         if(osinstance->getVariableNumber() > 0){
00609                 for(i = 0; i < osinstance->getVariableNumber(); i++){
00610                         if(osinstance->getVariableNames() != NULL) cout << "variable Names  " << osinstance->getVariableNames()[ i]  << endl;
00611                         if(osinstance->getVariableTypes() != NULL) cout << "variable Types  " << osinstance->getVariableTypes()[ i]  << endl;
00612                         if(osinstance->getVariableLowerBounds() != NULL) cout << "variable Lower Bounds  " << osinstance->getVariableLowerBounds()[ i]  << endl;
00613                         if(osinstance->getVariableUpperBounds() != NULL) cout << "variable Upper Bounds  " <<  osinstance->getVariableUpperBounds()[i] << endl;
00614                 }
00615         }
00616         
00617         // print out objective function information
00618         if(osinstance->getVariableNumber() > 0 || osinstance->instanceData->objectives->obj != NULL || osinstance->instanceData->objectives->numberOfObjectives > 0){
00619                 if( osinstance->getObjectiveMaxOrMins()[0] == "min")  cout <<  "problem is a minimization" << endl;
00620                 else cout <<  "problem is a maximization" << endl;
00621                 for(i = 0; i < osinstance->getVariableNumber(); i++){
00622                         cout << "OBJ COEFFICIENT =  " <<  osinstance->getDenseObjectiveCoefficients()[0][i] << endl;
00623                 }
00624         }
00625         // print out constraint information
00626         if(osinstance->getConstraintNumber() > 0){
00627                 for(i = 0; i < osinstance->getConstraintNumber(); i++){
00628                         if(osinstance->getConstraintNames() != NULL) cout << "row name = " << osinstance->getConstraintNames()[i] <<  endl;
00629                         if(osinstance->getConstraintLowerBounds() != NULL) cout << "row lower bound = " << osinstance->getConstraintLowerBounds()[i] <<  endl;
00630                         if(osinstance->getConstraintUpperBounds() != NULL) cout << "row upper bound = " << osinstance->getConstraintUpperBounds()[i] <<  endl; 
00631                 }
00632         }
00633         
00634         // print out linear constraint data
00635         cout << endl;
00636         cout << "number of nonzeros =  " << osinstance->getLinearConstraintCoefficientNumber() << endl;
00637         for(i = 0; i <= osinstance->getVariableNumber(); i++){
00638                 cout << "Start Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->starts[ i] << endl;
00639         }
00640         cout << endl;
00641         for(i = 0; i < osinstance->getLinearConstraintCoefficientNumber(); i++){
00642                 cout << "Index Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->indexes[i] << endl;
00643                 cout << "Nonzero Value =  " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->values[i] << endl;
00644         }
00645 
00646         // print out quadratic data
00647         cout << "number of qterms =  " <<  osinstance->getNumberOfQuadraticTerms() << endl;
00648         for(int i = 0; i <  osinstance->getNumberOfQuadraticTerms(); i++){
00649                 cout << "Row Index = " <<  osinstance->getQuadraticTerms()->rowIndexes[i] << endl;
00650                 cout << "Var Index 1 = " << osinstance->getQuadraticTerms()->varOneIndexes[ i] << endl;
00651                 cout << "Var Index 2 = " << osinstance->getQuadraticTerms()->varTwoIndexes[ i] << endl;
00652                 cout << "Coefficient = " << osinstance->getQuadraticTerms()->coefficients[ i] << endl;
00653         }
00654 } // end dataEchoCheck
00655 
00656 
00657 IpoptProblem::IpoptProblem(OSInstance *osinstance_,  OSResult *osresult_) {
00658         osinstance = osinstance_;
00659         osresult = osresult_;
00660 }
00661 
00662 IpoptProblem::~IpoptProblem() {
00663 
00664 }
00665 
00666 
00667 

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