00001
00018
00019
00020 #include "OSIpoptSolver.h"
00021 #include "OSGeneral.h"
00022 #include "OSParameters.h"
00023 #include "OSMathUtil.h"
00024 #include "CoinFinite.hpp"
00025
00026
00027 using std::cout;
00028 using std::endl;
00029 using std::ostringstream;
00030
00031
00032
00033 IpoptSolver::IpoptSolver() {
00034 osrlwriter = new OSrLWriter();
00035 osresult = new OSResult();
00036 m_osilreader = NULL;
00037 m_osolreader = NULL;
00038 ipoptErrorMsg = new std::string("");
00039 }
00040
00041 IpoptSolver::~IpoptSolver() {
00042 #ifdef DEBUG
00043 cout << "inside IpoptSolver destructor" << endl;
00044 #endif
00045 if(m_osilreader != NULL) delete m_osilreader;
00046 m_osilreader = NULL;
00047 if(m_osolreader != NULL) delete m_osolreader;
00048 m_osolreader = NULL;
00049 delete osresult;
00050 osresult = NULL;
00051 delete osrlwriter;
00052 osrlwriter = NULL;
00053
00054
00055 delete ipoptErrorMsg;
00056 #ifdef DEBUG
00057 cout << "leaving IpoptSolver destructor" << endl;
00058 #endif
00059 }
00060
00061
00062 bool IpoptProblem::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00063 Index& nnz_h_lag, IndexStyleEnum& index_style)
00064 {
00065 try{
00066
00067 if( (osinstance->getNumberOfIntegerVariables() + osinstance->getNumberOfBinaryVariables()) > 0 )
00068 throw ErrorClass("Ipopt does not solve integer programs -- please try Bonmin or Couenne");
00069
00070 n = osinstance->getVariableNumber();
00071
00072 m = osinstance->getConstraintNumber();
00073 #ifdef DEBUG
00074 cout << "number variables !!!!!!!!!!!!!!!!!!!!!!!!!!!" << n << endl;
00075 cout << "number constraints !!!!!!!!!!!!!!!!!!!!!!!!!!!" << m << endl;
00076 #endif
00077 try{
00078 osinstance->initForAlgDiff( );
00079 }
00080 catch(const ErrorClass& eclass){
00081 #ifdef DEBUG
00082 cout << "error in OSIpoptSolver, line 78:\n" << eclass.errormsg << endl;
00083 #endif
00084 *ipoptErrorMsg = eclass.errormsg;
00085 throw;
00086 }
00087
00088 osinstance->bUseExpTreeForFunEval = true;
00089
00090 SparseJacobianMatrix *sparseJacobian = NULL;
00091 try{
00092 sparseJacobian = osinstance->getJacobianSparsityPattern();
00093 }
00094 catch(const ErrorClass& eclass){
00095 #ifdef DEBUG
00096 cout << "error in OSIpoptSolver, line 91:\n" << eclass.errormsg << endl;
00097 #endif
00098 *ipoptErrorMsg = eclass.errormsg;
00099 throw;
00100 }
00101
00102 if (sparseJacobian != NULL){
00103 nnz_jac_g = sparseJacobian->valueSize;
00104 }else{
00105 nnz_jac_g = 0;
00106 }
00107
00108 #ifdef DEBUG
00109 cout << "nnz_jac_g !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_jac_g << endl;
00110 #endif
00111
00112
00113 if( (osinstance->getNumberOfNonlinearExpressions() == 0) && (osinstance->getNumberOfQuadraticTerms() == 0) ) {
00114 #ifdef DEBUG
00115 cout << "This is a linear program" << endl;
00116 #endif
00117 nnz_h_lag = 0;
00118 }
00119 else{
00120
00121 SparseHessianMatrix *sparseHessian = osinstance->getLagrangianHessianSparsityPattern();
00122
00123 if(sparseHessian != NULL){
00124 nnz_h_lag = sparseHessian->hessDimension;
00125 }else{
00126 nnz_h_lag = 0;
00127 }
00128 }
00129 #ifdef DEBUG
00130 cout << "print nnz_h_lag (OSIpoptSolver.cpp)" << endl;
00131 cout << "nnz_h_lag !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_h_lag << endl;
00132 cout << "set index_style (OSIpoptSolver.cpp)" << endl;
00133 #endif
00134
00135 index_style = TNLP::C_STYLE;
00136 #ifdef DEBUG
00137 cout << "return from get_nlp_info (OSIpoptSolver.cpp)" << nnz_h_lag << endl;
00138 #endif
00139
00141
00142 return true;
00143 }
00144 catch(const ErrorClass& eclass){
00145
00146 *ipoptErrorMsg = eclass.errormsg;
00147 throw;
00148 }
00149
00150 }
00151
00152
00153 bool IpoptProblem::get_bounds_info(Index n, Number* x_l, Number* x_u,
00154 Index m, Number* g_l, Number* g_u){
00155 int i;
00156 double * mdVarLB = osinstance->getVariableLowerBounds();
00157
00158
00159 double * mdVarUB = osinstance->getVariableUpperBounds();
00160
00161 for(i = 0; i < n; i++){
00162 x_l[ i] = mdVarLB[ i];
00163 x_u[ i] = mdVarUB[ i];
00164
00165
00166 }
00167
00168
00169
00170
00171
00172
00173 double * mdConLB = osinstance->getConstraintLowerBounds();
00174
00175 double * mdConUB = osinstance->getConstraintUpperBounds();
00176
00177 for(int i = 0; i < m; i++){
00178 g_l[ i] = mdConLB[ i];
00179 g_u[ i] = mdConUB[ i];
00180
00181
00182 }
00183 return true;
00184 }
00185
00186
00187
00188 bool IpoptProblem::get_starting_point(Index n, bool init_x, Number* x,
00189 bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda,
00190 Number* lambda) {
00191
00192
00193
00194 assert(init_x == true);
00195 assert(init_z == false);
00196 assert(init_lambda == false);
00197 int i, m1, n1;
00198
00199 #ifdef DEBUG
00200 cout << "get initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00201 #endif
00202
00203
00204 #ifdef DEBUG
00205 cout << "get number of initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00206 #endif
00207 int k;
00208 if (osoption != NULL)
00209 m1 = osoption->getNumberOfInitVarValues();
00210 else
00211 m1 = 0;
00212 #ifdef DEBUG
00213 cout << "number of variables initialed: " << m1 << endl;
00214 #endif
00215
00216 n1 = osinstance->getVariableNumber();
00217 bool* initialed;
00218 initialed = new bool[n1];
00219 #ifdef DEBUG
00220 cout << "number of variables in total: " << n1 << endl;
00221 #endif
00222
00223 for(k = 0; k < n1; k++)
00224 initialed[k] = false;
00225
00226 if (m1 > 0)
00227 {
00228 #ifdef DEBUG
00229 cout << "get initial values " << endl;
00230 #endif
00231
00232 InitVarValue** initVarVector = osoption->getInitVarValuesSparse();
00233 #ifdef DEBUG
00234 cout << "done " << endl;
00235 #endif
00236 try
00237 {
00238 double initval;
00239 for(k = 0; k < m1; k++)
00240 { i = initVarVector[k]->idx;
00241 if (initVarVector[k]->idx > n1)
00242 throw ErrorClass ("Illegal index value in variable initialization");
00243
00244 initval = initVarVector[k]->value;
00245 if (osinstance->instanceData->variables->var[k]->ub == OSDBL_MAX)
00246 { if (osinstance->instanceData->variables->var[k]->lb > initval)
00247 throw ErrorClass ("Initial value outside of bounds");
00248 }
00249 else
00250 if (osinstance->instanceData->variables->var[k]->lb == -OSDBL_MAX)
00251 { if (osinstance->instanceData->variables->var[k]->ub < initval)
00252 throw ErrorClass ("Initial value outside of bounds");
00253 }
00254 else
00255 { if ((osinstance->instanceData->variables->var[k]->lb > initval) ||
00256 (osinstance->instanceData->variables->var[k]->ub < initval))
00257 throw ErrorClass ("Initial value outside of bounds");
00258 }
00259
00260 x[initVarVector[k]->idx] = initval;
00261 initialed[initVarVector[k]->idx] = true;
00262 }
00263 }
00264 catch(const ErrorClass& eclass)
00265 { cout << "Error in IpoptProblem::get_starting_point (OSIpoptSolver.cpp, line 247)";
00266 cout << endl << endl << endl;
00267 }
00268 }
00269
00270 double default_initval;
00271 default_initval = 1.7171;
00272
00273
00274 for(k = 0; k < n1; k++)
00275 { if (!initialed[k])
00276 if (osinstance->instanceData->variables->var[k]->ub == OSDBL_MAX)
00277 if (osinstance->instanceData->variables->var[k]->lb <= default_initval)
00278 x[k] = default_initval;
00279 else
00280 x[k] = osinstance->instanceData->variables->var[k]->lb;
00281 else
00282 if (osinstance->instanceData->variables->var[k]->lb == -OSDBL_MAX)
00283 if (osinstance->instanceData->variables->var[k]->ub >= default_initval)
00284 x[k] = default_initval;
00285 else
00286 x[k] = osinstance->instanceData->variables->var[k]->ub;
00287 else
00288 if ((osinstance->instanceData->variables->var[k]->lb <= default_initval) &&
00289 (osinstance->instanceData->variables->var[k]->ub >= default_initval))
00290 x[k] = default_initval;
00291 else
00292 if (osinstance->instanceData->variables->var[k]->lb > default_initval)
00293 x[k] = osinstance->instanceData->variables->var[k]->lb;
00294 else
00295 x[k] = osinstance->instanceData->variables->var[k]->ub;
00296 }
00297
00298 #ifdef DEBUG
00299 for(i = 0; i < n1; i++){
00300 std::cout << "INITIAL VALUE !!!!!!!!!!!!!!!!!!!! " << x[ i] << std::endl;
00301 }
00302 #endif
00303
00304 osinstance->calculateAllObjectiveFunctionValues( x, true);
00305 delete[] initialed;
00306
00307 return true;
00308 }
00309
00310
00311 bool IpoptProblem::eval_f(Index n, const Number* x, bool new_x, Number& obj_value){
00312 try{
00313 if(osinstance->getObjectiveNumber() > 0){
00314
00315
00316 if(new_x == false) obj_value = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), false)[ 0];
00317 else obj_value = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, true, 0 )[ 0];
00318
00319
00320 }
00321
00322 }
00323 catch(const ErrorClass& eclass){
00324
00325 *ipoptErrorMsg = eclass.errormsg;
00326 throw;
00327 }
00328
00329 return true;
00330 }
00331
00332 bool IpoptProblem::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f){
00333 int i;
00334 double *objGrad = NULL;
00335 if(osinstance->getObjectiveNumber() > 0){
00336 try{
00337
00338 objGrad = osinstance->calculateObjectiveFunctionGradient( const_cast<double*>(x), NULL, NULL, -1, new_x, 1);
00339 }
00340 catch(const ErrorClass& eclass){
00341 #ifdef DEBUG
00342 cout << "error in OSIpoptSolver, line 314:\n" << eclass.errormsg << endl;
00343 #endif
00344 *ipoptErrorMsg = eclass.errormsg;
00345 throw;
00346 }
00347 for(i = 0; i < n; i++){
00348 grad_f[ i] = objGrad[ i];
00349 }
00350 }
00351 return true;
00352 }
00353
00354
00355 bool IpoptProblem::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) {
00356 try{
00357 double *conVals = osinstance->calculateAllConstraintFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 );
00358 int i;
00359 for(i = 0; i < m; i++){
00360 if( CoinIsnan( (double)conVals[ i] ) ) return false;
00361 g[i] = conVals[ i] ;
00362 }
00363 return true;
00364 }
00365 catch(const ErrorClass& eclass){
00366 #ifdef DEBUG
00367 cout << "error in OSIpoptSolver, line 338:\n" << eclass.errormsg << endl;
00368 #endif
00369 *ipoptErrorMsg = eclass.errormsg;
00370 throw;
00371 }
00372 }
00373
00374
00375
00376 bool IpoptProblem::eval_jac_g(Index n, const Number* x, bool new_x,
00377 Index m, Index nele_jac, Index* iRow, Index *jCol,
00378 Number* values){
00379 SparseJacobianMatrix *sparseJacobian;
00380 if (values == NULL) {
00381
00382
00383
00384
00385
00386 try{
00387 sparseJacobian = osinstance->getJacobianSparsityPattern();
00388 }
00389 catch(const ErrorClass& eclass){
00390 #ifdef DEBUG
00391 cout << "error in OSIpoptSolver, line 362:\n" << eclass.errormsg << endl;
00392 #endif
00393 *ipoptErrorMsg = eclass.errormsg;
00394 throw;
00395 }
00396 int i = 0;
00397 int k, idx;
00398 for(idx = 0; idx < m; idx++){
00399 for(k = *(sparseJacobian->starts + idx); k < *(sparseJacobian->starts + idx + 1); k++){
00400 iRow[i] = idx;
00401 jCol[i] = *(sparseJacobian->indexes + k);
00402
00403
00404 i++;
00405 }
00406 }
00407 }
00408 else {
00409
00410 try{
00411 sparseJacobian = osinstance->calculateAllConstraintFunctionGradients( const_cast<double*>(x), NULL, NULL, new_x, 1);
00412 }
00413 catch(const ErrorClass& eclass){
00414 #ifdef DEBUG
00415 cout << "error in OSIpoptSolver, line 386:\n" << eclass.errormsg << endl;
00416 #endif
00417 *ipoptErrorMsg = eclass.errormsg;
00418 throw;
00419 }
00420
00421 for(int i = 0; i < nele_jac; i++){
00422 values[ i] = sparseJacobian->values[i];
00423
00424
00425
00426 }
00427 }
00428 return true;
00429 }
00430
00431
00432 bool IpoptProblem::eval_h(Index n, const Number* x, bool new_x,
00433 Number obj_factor, Index m, const Number* lambda,
00434 bool new_lambda, Index nele_hess, Index* iRow,
00435 Index* jCol, Number* values){
00436
00438 SparseHessianMatrix *sparseHessian;
00439
00440 int i;
00441 if (values == NULL) {
00442
00443
00444 try{
00445 sparseHessian = osinstance->getLagrangianHessianSparsityPattern( );
00446 }
00447 catch(const ErrorClass& eclass){
00448
00449 *ipoptErrorMsg = eclass.errormsg;
00450 throw;
00451 }
00452
00453 for(i = 0; i < nele_hess; i++){
00454 iRow[i] = *(sparseHessian->hessColIdx + i);
00455 jCol[i] = *(sparseHessian->hessRowIdx + i);
00456
00457
00458 }
00459 }
00460 else {
00461
00462
00463 double* objMultipliers = new double[1];
00464 objMultipliers[0] = obj_factor;
00465 try{
00466 sparseHessian = osinstance->calculateLagrangianHessian( const_cast<double*>(x), objMultipliers, const_cast<double*>(lambda) , new_x, 2);
00467 delete[] objMultipliers;
00468 }
00469 catch(const ErrorClass& eclass){
00470 #ifdef DEBUG
00471 cout << "error in OSIpoptSolver, line 444:\n" << eclass.errormsg << endl;
00472 #endif
00473 *ipoptErrorMsg = eclass.errormsg;
00474 delete[] objMultipliers;
00475 throw;
00476 }
00477 for(i = 0; i < nele_hess; i++){
00478 values[ i] = *(sparseHessian->hessValues + i);
00479 }
00480 }
00482 return true;
00483 }
00484
00485 bool IpoptProblem::get_scaling_parameters(Number& obj_scaling,
00486 bool& use_x_scaling, Index n,
00487 Number* x_scaling,
00488 bool& use_g_scaling, Index m,
00489 Number* g_scaling){
00490 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
00491 obj_scaling = -1;
00492 }
00493 else obj_scaling = 1;
00494 use_x_scaling = false;
00495 use_g_scaling = false;
00496 return true;
00497 }
00498
00499 void IpoptProblem::finalize_solution(SolverReturn status,
00500 Index n, const Number* x, const Number* z_L, const Number* z_U,
00501 Index m, const Number* g, const Number* lambda, Number obj_value,
00502 const IpoptData* ip_data, IpoptCalculatedQuantities* ip_cq){
00503
00504
00505
00506 if(osinstance->getObjectiveNumber() > 0)
00507 obj_value = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), true)[ 0];
00508 OSrLWriter *osrlwriter ;
00509 osrlwriter = new OSrLWriter();
00510 #ifdef DEBUG
00511 printf("\n\nSolution of the primal variables, x\n");
00512 for (Index i=0; i<n; i++) {
00513 printf("x[%d] = %e\n", i, x[i]);
00514 }
00515
00516 printf("\n\nSolution of the bound multipliers, z_L and z_U\n");
00517 for (Index i=0; i<n; i++) {
00518 printf("z_L[%d] = %e\n", i, z_L[i]);
00519 }
00520 for (Index i=0; i<n; i++) {
00521 printf("z_U[%d] = %e\n", i, z_U[i]);
00522 }
00523 #endif
00524 if(osinstance->getObjectiveNumber() > 0) printf("\nObjective value f(x*) = %e\n", obj_value );
00525
00526 int solIdx = 0;
00527 int numberOfOtherVariableResults;
00528 int otherIdx;
00529 int numCon = osinstance->getConstraintNumber();
00530 ostringstream outStr;
00531
00532
00533
00534 double *dualValue = NULL;
00535 std::string *rcost = NULL;
00536 int* idx = NULL;
00537 double* mdObjValues = NULL;
00538 if(osinstance->getObjectiveNumber() > 0){
00539 mdObjValues = new double[1];
00540 }
00541
00542 std::string message = "Ipopt solver finishes to the end.";
00543 std::string solutionDescription = "";
00544
00545 try{
00546
00547 if(osresult->setSolverInvoked( "COIN-OR Ipopt") != true)
00548 throw ErrorClass("OSResult error: setSolverInvoked");
00549 if(osresult->setServiceName( getVersionInfo()) != true)
00550 throw ErrorClass("OSResult error: setServiceName");
00551 if(osresult->setInstanceName( osinstance->getInstanceName()) != true)
00552 throw ErrorClass("OSResult error: setInstanceName");
00553
00554
00555
00556
00557
00558 if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00559 throw ErrorClass("OSResult error: setVariableNumer");
00560 if(osresult->setObjectiveNumber( 1) != true)
00561 throw ErrorClass("OSResult error: setObjectiveNumber");
00562 if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00563 throw ErrorClass("OSResult error: setConstraintNumber");
00564 if(osresult->setSolutionNumber( 1) != true)
00565 throw ErrorClass("OSResult error: setSolutionNumer");
00566
00567
00568 if(osresult->setGeneralMessage( message) != true)
00569 throw ErrorClass("OSResult error: setGeneralMessage");
00570
00571 switch( status){
00572 case SUCCESS:
00573 solutionDescription = "SUCCESS[IPOPT]: Algorithm terminated normally at a locally optimal point, satisfying the convergence tolerances.";
00574 osresult->setSolutionStatus(solIdx, "locallyOptimal", solutionDescription);
00575
00576
00577 if(osinstance->getVariableNumber() > 0) osresult->setPrimalVariableValuesDense(solIdx, const_cast<double*>(x));
00578
00579 if(numCon > 0){
00580 dualValue = new double[ numCon];
00581 for (Index i=0; i < numCon; i++) {
00582 dualValue[ i] = -lambda[ i];
00583 }
00584
00585 osresult->setDualVariableValuesDense(solIdx, dualValue);
00586 }
00587
00588 if(osinstance->getObjectiveNumber() > 0){
00589 mdObjValues[0] = obj_value ;
00590 osresult->setObjectiveValuesDense(solIdx, mdObjValues);
00591 }
00592
00593
00594 if(osinstance->getVariableNumber() > 0){
00595 numberOfOtherVariableResults = 1;
00596 osresult->setNumberOfOtherVariableResults(solIdx, numberOfOtherVariableResults);
00597 rcost = new std::string[ osinstance->getVariableNumber()];
00598 idx = new int[ osinstance->getVariableNumber()];
00599 for (Index i = 0; i < n; i++) {
00600 rcost[ i] = os_dtoa_format( z_L[i] - z_U[i]);
00601 idx[ i] = i;
00602 }
00603 otherIdx = 0;
00604 osresult->setAnOtherVariableResultSparse(solIdx, otherIdx, "reduced costs", "", "the variable reduced costs",
00605 idx, rcost, osinstance->getVariableNumber());
00606 }
00607
00608
00609
00610
00611
00612
00613
00614
00615
00616
00617
00618
00619
00620
00621
00622
00623 if(osinstance->getVariableNumber() > 0) delete[] rcost;
00624 if(osinstance->getVariableNumber() > 0) delete[] idx;
00625 if(osinstance->getConstraintNumber() > 0) delete[] dualValue;
00626
00627
00628
00629
00630 break;
00631 case MAXITER_EXCEEDED:
00632 solutionDescription = "MAXITER_EXCEEDED[IPOPT]: Maximum number of iterations exceeded.";
00633 osresult->setSolutionStatus(solIdx, "stoppedByLimit", solutionDescription);
00634 if(x != NULL) osresult->setPrimalVariableValuesDense(solIdx, const_cast<double*>(x));
00635 if(lambda != NULL) osresult->setDualVariableValuesDense(solIdx, const_cast<double*>( lambda));
00636 if(osinstance->getObjectiveNumber() > 0){
00637 mdObjValues[0] = obj_value ;
00638 osresult->setObjectiveValuesDense(solIdx, mdObjValues);
00639 }
00640 break;
00641 case STOP_AT_TINY_STEP:
00642 solutionDescription = "STOP_AT_TINY_STEP[IPOPT]: Algorithm proceeds with very little progress.";
00643 osresult->setSolutionStatus(solIdx, "stoppedByLimit", solutionDescription);
00644 if(x != NULL) osresult->setPrimalVariableValuesDense(solIdx, const_cast<double*>( x));
00645 if(lambda != NULL) osresult->setDualVariableValuesDense(solIdx, const_cast<double*>( lambda));
00646 if(osinstance->getObjectiveNumber() > 0){
00647 mdObjValues[0] = obj_value ;
00648 osresult->setObjectiveValuesDense(solIdx, mdObjValues);
00649 }
00650 break;
00651 case STOP_AT_ACCEPTABLE_POINT:
00652 solutionDescription = "STOP_AT_ACCEPTABLE_POINT[IPOPT]: Algorithm stopped at a point that was converged, not to _desired_ tolerances, but to _acceptable_ tolerances";
00653 osresult->setSolutionStatus(solIdx, "IpoptAccetable", solutionDescription);
00654 if(lambda != NULL) osresult->setDualVariableValuesDense(solIdx, const_cast<double*>( lambda));
00655 if(x != NULL)osresult->setPrimalVariableValuesDense(solIdx, const_cast<double*>(x));
00656 if(osinstance->getObjectiveNumber() > 0){
00657 mdObjValues[0] = obj_value ;
00658 osresult->setObjectiveValuesDense(solIdx, mdObjValues);
00659 }
00660 break;
00661 case LOCAL_INFEASIBILITY:
00662 solutionDescription = "LOCAL_INFEASIBILITY[IPOPT]: Algorithm converged to a point of local infeasibility. Problem may be infeasible.";
00663 osresult->setSolutionStatus(solIdx, "infeasible", solutionDescription);
00664 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00665 break;
00666 case USER_REQUESTED_STOP:
00667 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.";
00668 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00669 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00670 break;
00671 case DIVERGING_ITERATES:
00672 solutionDescription = "DIVERGING_ITERATES[IPOPT]: It seems that the iterates diverge.";
00673 osresult->setSolutionStatus(solIdx, "unbounded", solutionDescription);
00674 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00675 break;
00676 case RESTORATION_FAILURE:
00677 solutionDescription = "RESTORATION_FAILURE[IPOPT]: Restoration phase failed, algorithm doesn't know how to proceed.";
00678 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00679 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00680 break;
00681 case ERROR_IN_STEP_COMPUTATION:
00682 solutionDescription = "ERROR_IN_STEP_COMPUTATION[IPOPT]: An unrecoverable error occurred while IPOPT tried to compute the search direction.";
00683 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00684 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00685 break;
00686 case INVALID_NUMBER_DETECTED:
00687 solutionDescription = "INVALID_NUMBER_DETECTED[IPOPT]: Algorithm received an invalid number (such as NaN or Inf) from the NLP; see also option check_derivatives_for_naninf.";
00688 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00689 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00690 break;
00691 case INTERNAL_ERROR:
00692 solutionDescription = "INTERNAL_ERROR[IPOPT]: An unknown internal error occurred. Please contact the IPOPT authors through the mailing list.";
00693 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00694 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00695 break;
00696 default:
00697 solutionDescription = "OTHER[IPOPT]: other unknown solution status from Ipopt solver";
00698 osresult->setSolutionStatus(solIdx, "other", solutionDescription);
00699 if( osinstance->getVariableNumber() == 0) osresult->setSolutionMessage(solIdx, "Warning: this problem has zero decision variables!");
00700 }
00701 osresult->setGeneralStatusType("normal");
00702 delete osrlwriter;
00703 if(osinstance->getObjectiveNumber() > 0){
00704 delete[] mdObjValues;
00705 }
00706 osrlwriter = NULL;
00707
00708 }
00709 catch(const ErrorClass& eclass){
00710 #ifdef DEBUG
00711 cout << "error in OSIpoptSolver, line 636:\n" << eclass.errormsg << endl;
00712 #endif
00713 osresult->setGeneralMessage( eclass.errormsg);
00714 osresult->setGeneralStatusType( "error");
00715 std::string osrl = osrlwriter->writeOSrL( osresult);
00716 delete osrlwriter;
00717 osrlwriter = NULL;
00718 throw ErrorClass( osrl) ;
00719 if(osinstance->getObjectiveNumber() > 0){
00720 delete[] mdObjValues;
00721 }
00722 mdObjValues = NULL;
00723 }
00725 }
00726
00727
00728 void IpoptSolver::setSolverOptions() throw (ErrorClass) {
00729 try{
00730
00731 this->bSetSolverOptions = true;
00732
00733
00734 app->Options()->SetIntegerValue("print_level", 0);
00735 app->Options()->SetIntegerValue("max_iter", 20000);
00736 app->Options()->SetNumericValue("bound_relax_factor", 0, true, true);
00737 app->Options()->SetStringValue("mu_strategy", "adaptive", true, true);
00738
00739 app->Options()->SetStringValue("check_derivatives_for_naninf", "yes");
00740
00741 if( (osinstance-> getNumberOfNonlinearExpressions() <= 0) && (osinstance->getNumberOfQuadraticTerms() <= 0) ){
00742 app->Options()->SetStringValue("hessian_constant", "yes", true, true);
00743 }
00744 if(osinstance->getObjectiveNumber() > 0){
00745 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
00746 app->Options()->SetStringValue("nlp_scaling_method", "user-scaling");
00747 }
00748 }
00749
00750
00751
00752 if(osoption == NULL && osol.length() > 0)
00753 {
00754 m_osolreader = new OSoLReader();
00755 osoption = m_osolreader->readOSoL( osol);
00756 }
00757
00758 if( osoption != NULL && osoption->getNumberOfSolverOptions() > 0 ){
00759
00760 std::vector<SolverOption*> optionsVector;
00761 optionsVector = osoption->getSolverOptions( "ipopt");
00762 char *pEnd;
00763 int i;
00764 int num_ipopt_options = optionsVector.size();
00765 for(i = 0; i < num_ipopt_options; i++){
00766
00767 if(optionsVector[ i]->type == "numeric" ){
00768
00769 app->Options()->SetNumericValue(optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
00770 }
00771 else if(optionsVector[ i]->type == "integer" ){
00772
00773 app->Options()->SetIntegerValue(optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ) );
00774 }
00775 else if(optionsVector[ i]->type == "string" ){
00776
00777 app->Options()->SetStringValue(optionsVector[ i]->name, optionsVector[ i]->value);
00778 }
00779 }
00780 }
00781 }
00782 catch(const ErrorClass& eclass){
00783 std::cout << "THERE IS AN ERROR" << std::endl;
00784 osresult->setGeneralMessage( eclass.errormsg);
00785 osresult->setGeneralStatusType( "error");
00786 osrl = osrlwriter->writeOSrL( osresult);
00787 throw ErrorClass( osrl) ;
00788 }
00789 }
00790
00791
00792
00793 void IpoptSolver::buildSolverInstance() throw (ErrorClass) {
00794 try{
00795
00796 if(osil.length() == 0 && osinstance == NULL) throw ErrorClass("there is no instance");
00797 if(osinstance == NULL){
00798 m_osilreader = new OSiLReader();
00799 osinstance = m_osilreader->readOSiL( osil);
00800 }
00801
00802 nlp = new IpoptProblem( osinstance, osoption, osresult, ipoptErrorMsg);
00803 app = new IpoptApplication();
00804 this->bCallbuildSolverInstance = true;
00805 }
00806 catch(const ErrorClass& eclass){
00807 #ifdef DEBUG
00808 cout << "error in OSIpoptSolver, line 722:\n" << eclass.errormsg << endl;
00809 #endif
00810 std::cout << "THERE IS AN ERROR" << std::endl;
00811 osresult->setGeneralMessage( eclass.errormsg);
00812 osresult->setGeneralStatusType( "error");
00813 osrl = osrlwriter->writeOSrL( osresult);
00814 throw ErrorClass( osrl) ;
00815 }
00816 }
00817
00818
00819
00820 void IpoptSolver::solve() throw (ErrorClass) {
00821 if( this->bCallbuildSolverInstance == false) buildSolverInstance();
00822 if( this->bSetSolverOptions == false) setSolverOptions();
00823 try{
00824 clock_t start, finish;
00825 double duration;
00826 start = clock();
00827
00828
00829
00830 finish = clock();
00831 duration = (double) (finish - start) / CLOCKS_PER_SEC;
00832
00833
00834
00835
00836
00837
00838
00839 app->Initialize();
00840
00841
00842
00843
00844 ApplicationReturnStatus status = app->OptimizeTNLP( nlp);
00845
00846 osrl = osrlwriter->writeOSrL( osresult);
00847
00848
00849 if (status < -2) {
00850 throw ErrorClass("Ipopt FAILED TO SOLVE THE PROBLEM: " + *ipoptErrorMsg);
00851 }
00852 }
00853 catch(const ErrorClass& eclass){
00854 #ifdef DEBUG
00855 cout << "error in OSIpoptSolver, line 775:\n" << eclass.errormsg << endl;
00856 #endif
00857 osresult->setGeneralMessage( eclass.errormsg);
00858 osresult->setGeneralStatusType( "error");
00859 osrl = osrlwriter->writeOSrL( osresult);
00860 throw ErrorClass( osrl) ;
00861 }
00862 }
00863
00864
00865 void IpoptSolver::dataEchoCheck(){
00866
00867 int i;
00868
00869
00870 cout << "This is problem: " << osinstance->getInstanceName() << endl;
00871 cout << "The problem source is: " << osinstance->getInstanceSource() << endl;
00872 cout << "The problem description is: " << osinstance->getInstanceDescription() << endl;
00873 cout << "number of variables = " << osinstance->getVariableNumber() << endl;
00874 cout << "number of Rows = " << osinstance->getConstraintNumber() << endl;
00875
00876
00877 if(osinstance->getVariableNumber() > 0){
00878 for(i = 0; i < osinstance->getVariableNumber(); i++){
00879 if(osinstance->getVariableNames() != NULL) cout << "variable Names " << osinstance->getVariableNames()[ i] << endl;
00880 if(osinstance->getVariableTypes() != NULL) cout << "variable Types " << osinstance->getVariableTypes()[ i] << endl;
00881 if(osinstance->getVariableLowerBounds() != NULL) cout << "variable Lower Bounds " << osinstance->getVariableLowerBounds()[ i] << endl;
00882 if(osinstance->getVariableUpperBounds() != NULL) cout << "variable Upper Bounds " << osinstance->getVariableUpperBounds()[i] << endl;
00883 }
00884 }
00885
00886
00887 if(osinstance->getVariableNumber() > 0 || osinstance->instanceData->objectives->obj != NULL || osinstance->instanceData->objectives->numberOfObjectives > 0){
00888 if( osinstance->getObjectiveMaxOrMins()[0] == "min") cout << "problem is a minimization" << endl;
00889 else cout << "problem is a maximization" << endl;
00890 for(i = 0; i < osinstance->getVariableNumber(); i++){
00891 cout << "OBJ COEFFICIENT = " << osinstance->getDenseObjectiveCoefficients()[0][i] << endl;
00892 }
00893 }
00894
00895 if(osinstance->getConstraintNumber() > 0){
00896 for(i = 0; i < osinstance->getConstraintNumber(); i++){
00897 if(osinstance->getConstraintNames() != NULL) cout << "row name = " << osinstance->getConstraintNames()[i] << endl;
00898 if(osinstance->getConstraintLowerBounds() != NULL) cout << "row lower bound = " << osinstance->getConstraintLowerBounds()[i] << endl;
00899 if(osinstance->getConstraintUpperBounds() != NULL) cout << "row upper bound = " << osinstance->getConstraintUpperBounds()[i] << endl;
00900 }
00901 }
00902
00903
00904 cout << endl;
00905 cout << "number of nonzeros = " << osinstance->getLinearConstraintCoefficientNumber() << endl;
00906 for(i = 0; i <= osinstance->getVariableNumber(); i++){
00907 cout << "Start Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->starts[ i] << endl;
00908 }
00909 cout << endl;
00910 for(i = 0; i < osinstance->getLinearConstraintCoefficientNumber(); i++){
00911 cout << "Index Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->indexes[i] << endl;
00912 cout << "Nonzero Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->values[i] << endl;
00913 }
00914
00915
00916 cout << "number of qterms = " << osinstance->getNumberOfQuadraticTerms() << endl;
00917 for(int i = 0; i < osinstance->getNumberOfQuadraticTerms(); i++){
00918 cout << "Row Index = " << osinstance->getQuadraticTerms()->rowIndexes[i] << endl;
00919 cout << "Var Index 1 = " << osinstance->getQuadraticTerms()->varOneIndexes[ i] << endl;
00920 cout << "Var Index 2 = " << osinstance->getQuadraticTerms()->varTwoIndexes[ i] << endl;
00921 cout << "Coefficient = " << osinstance->getQuadraticTerms()->coefficients[ i] << endl;
00922 }
00923 }
00924
00925
00926 IpoptProblem::IpoptProblem(OSInstance *osinstance_, OSOption *osoption_, OSResult *osresult_, std::string* ipoptErrorMsg_) {
00927 osinstance = osinstance_;
00928 osoption = osoption_;
00929 osresult = osresult_;
00930 ipoptErrorMsg = ipoptErrorMsg_;
00931 }
00932
00933 IpoptProblem::~IpoptProblem() {
00934
00935 }
00936
00937
00938