/home/coin/SVN-release/OS-1.0.0/OS/src/OSSolverInterfaces/IpoptSolver.cpp

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

Generated on Thu May 15 22:15:05 2008 by  doxygen 1.4.7