My Project
OSIpoptSolver.cpp
Go to the documentation of this file.
1
19#define DEBUG
20
21#include "OSIpoptSolver.h"
22#include "OSDataStructures.h"
23#include "OSParameters.h"
24#include "OSCommonUtil.h"
25#include "OSMathUtil.h"
26
27
28using std::cout;
29using std::endl;
30using std::ostringstream;
31using namespace Ipopt;
32
33
35 osrlwriter = new OSrLWriter();
36 osresult = new OSResult();
37 m_osilreader = NULL;
38 m_osolreader = NULL;
39 ipoptErrorMsg = "";
40}
41
43 #ifdef DEBUG
44 cout << "inside IpoptSolver destructor" << endl;
45 #endif
46 if(m_osilreader != NULL) delete m_osilreader;
47 m_osilreader = NULL;
48 if(m_osolreader != NULL) delete m_osolreader;
49 m_osolreader = NULL;
50 delete osresult;
51 osresult = NULL;
52 delete osrlwriter;
53 osrlwriter = NULL;
54 //delete osinstance;
55 //osinstance = NULL;
56 #ifdef DEBUG
57 cout << "leaving IpoptSolver destructor" << endl;
58 #endif
59}
60
61// returns the size of the problem
62bool IpoptProblem::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
63 Index& nnz_h_lag, IndexStyleEnum& index_style)
64{
65 if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Ipopt NEEDS AN OBJECTIVE FUNCTION");
66 // number of variables
68 // number of constraints
70#ifdef DEBUG
71 cout << "number variables !!!!!!!!!!!!!!!!!!!!!!!!!!!" << n << endl;
72 cout << "number constraints !!!!!!!!!!!!!!!!!!!!!!!!!!!" << m << endl;
73#endif
74 try{
76 }
77 catch(const ErrorClass& eclass){
78#ifdef DEBUG
79 cout << "error in OSIpoptSolver, line 78:\n" << eclass.errormsg << endl;
80#endif
81 ipoptErrorMsg = eclass.errormsg;
82 throw;
83 }
84 // use the OS Expression tree for function evaluations instead of CppAD
86 //std::cout << "Call sparse jacobian" << std::endl;
87 SparseJacobianMatrix *sparseJacobian = NULL;
88 try{
89 sparseJacobian = osinstance->getJacobianSparsityPattern();
90 }
91 catch(const ErrorClass& eclass){
92#ifdef DEBUG
93 cout << "error in OSIpoptSolver, line 91:\n" << eclass.errormsg << endl;
94#endif
95 ipoptErrorMsg = eclass.errormsg;
96 throw;
97 }
98 //std::cout << "Done calling sparse jacobian" << std::endl;
99 nnz_jac_g = sparseJacobian->valueSize;
100#ifdef DEBUG
101 cout << "nnz_jac_g !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_jac_g << endl;
102#endif
103 // nonzeros in upper hessian
104
106#ifdef DEBUG
107 cout << "This is a linear program" << endl;
108#endif
109 nnz_h_lag = 0;
110 }
111 else{
112 //std::cout << "Get Lagrangain Hessian Sparsity Pattern " << std::endl;
114 //std::cout << "Done Getting Lagrangain Hessian Sparsity Pattern " << std::endl;
115 nnz_h_lag = sparseHessian->hessDimension;
116 }
117#ifdef DEBUG
118 cout << "print nnz_h_lag (OSIpoptSolver.cpp)" << endl;
119 cout << "nnz_h_lag !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_h_lag << endl;
120 cout << "set index_style (OSIpoptSolver.cpp)" << endl;
121#endif
122 // use the C style indexing (0-based)
123 index_style = TNLP::C_STYLE;
124#ifdef DEBUG
125 cout << "return from get_nlp_info (OSIpoptSolver.cpp)" << nnz_h_lag << endl;
126#endif
127
129
130 return true;
131}//get_nlp_info
132
133
134bool IpoptProblem::get_bounds_info(Index n, Number* x_l, Number* x_u,
135 Index m, Number* g_l, Number* g_u){
136 int i;
137 double * mdVarLB = osinstance->getVariableLowerBounds();
138 //std::cout << "GET BOUNDS INFORMATION FOR IPOPT !!!!!!!!!!!!!!!!!" << std::endl;
139 // variables upper bounds
140 double * mdVarUB = osinstance->getVariableUpperBounds();
141
142 for(i = 0; i < n; i++){
143 x_l[ i] = mdVarLB[ i];
144 x_u[ i] = mdVarUB[ i];
145 //cout << "x_l !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_l[i] << endl;
146 //cout << "x_u !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_u[i] << endl;
147 }
148 // Ipopt interprets any number greater than nlp_upper_bound_inf as
149 // infinity. The default value of nlp_upper_bound_inf and nlp_lower_bound_inf
150 // is 1e19 and can be changed through ipopt options.
151 // e.g. g_u[0] = 2e19;
152
153 //constraint lower bounds
154 double * mdConLB = osinstance->getConstraintLowerBounds();
155 //constraint upper bounds
156 double * mdConUB = osinstance->getConstraintUpperBounds();
157
158 for(int i = 0; i < m; i++){
159 g_l[ i] = mdConLB[ i];
160 g_u[ i] = mdConUB[ i];
161 //cout << "lower !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_l[i] << endl;
162 //cout << "upper !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_u[i] << endl;
163 }
164 return true;
165}//get_bounds_info
166
167
168// returns the initial point for the problem
169bool IpoptProblem::get_starting_point(Index n, bool init_x, Number* x,
170 bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda,
171 Number* lambda) {
172 // Here, we assume we only have starting values for x, if you code
173 // your own NLP, you can provide starting values for the dual variables
174 // if you wish
175 assert(init_x == true);
176 assert(init_z == false);
177 assert(init_lambda == false);
178 int i, m1, n1;
179
180#ifdef DEBUG
181 cout << "get initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
182#endif
183
184 //now set initial values
185#ifdef DEBUG
186 cout << "get number of initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
187 cout << "Is osoption = NULL? " << (osoption == NULL) << endl;
188#endif
189 int k;
190 if (osoption != NULL)
192 else
193 m1 = 0;
194#ifdef DEBUG
195 cout << "number of variables initialed: " << m1 << endl;
196#endif
197
199 bool* initialed;
200 initialed = new bool[n1];
201#ifdef DEBUG
202 cout << "number of variables in total: " << n1 << endl;
203#endif
204
205 for(k = 0; k < n1; k++)
206 initialed[k] = false;
207
208 if (m1 > 0)
209 {
210#ifdef DEBUG
211 cout << "get initial values " << endl;
212#endif
213
214 InitVarValue** initVarVector = osoption->getInitVarValuesSparse();
215#ifdef DEBUG
216 cout << "done " << endl;
217#endif
218 try
219 {
220 double initval;
221 for(k = 0; k < m1; k++)
222 { cout << "process component " << k << " -- index " << initVarVector[k]->idx << endl;
223 i = initVarVector[k]->idx;
224 if (initVarVector[k]->idx > n1)
225 throw ErrorClass ("Illegal index value in variable initialization");
226
227 initval = initVarVector[k]->value;
229 { if (osinstance->instanceData->variables->var[k]->lb > initval)
230 throw ErrorClass ("Initial value outside of bounds");
231 }
232 else
234 { if (osinstance->instanceData->variables->var[k]->ub < initval)
235 throw ErrorClass ("Initial value outside of bounds");
236 }
237 else
238 { if ((osinstance->instanceData->variables->var[k]->lb > initval) ||
239 (osinstance->instanceData->variables->var[k]->ub < initval))
240 throw ErrorClass ("Initial value outside of bounds");
241 }
242
243 x[initVarVector[k]->idx] = initval;
244 initialed[initVarVector[k]->idx] = true;
245 }
246 }
247 catch(const ErrorClass& eclass)
248 { cout << "Error in IpoptProblem::get_starting_point (OSIpoptSolver.cpp, line 247)";
249 cout << endl << endl << endl;
250 }
251 } // end if (m1 > 0)
252
253 double default_initval;
254 default_initval = 1.7171;
255
256 for(k = 0; k < n1; k++)
257 { cout << "verify component " << k << endl;
258 if (!initialed[k])
260 if (osinstance->instanceData->variables->var[k]->lb <= default_initval)
261 x[k] = default_initval;
262 else
263 x[k] = osinstance->instanceData->variables->var[k]->lb;
264 else
266 if (osinstance->instanceData->variables->var[k]->ub >= default_initval)
267 x[k] = default_initval;
268 else
269 x[k] = osinstance->instanceData->variables->var[k]->ub;
270 else
271 if ((osinstance->instanceData->variables->var[k]->lb <= default_initval) &&
272 (osinstance->instanceData->variables->var[k]->ub >= default_initval))
273 x[k] = default_initval;
274 else
275 if (osinstance->instanceData->variables->var[k]->lb > default_initval)
276 x[k] = osinstance->instanceData->variables->var[k]->lb;
277 else
278 x[k] = osinstance->instanceData->variables->var[k]->ub;
279 }
280
281 for(i = 0; i < n1; i++){
282 std::cout << "INITIAL VALUE !!!!!!!!!!!!!!!!!!!! " << x[ i] << std::endl;
283 }
284
285
286 delete[] initialed;
287
288 return true;
289}//get_starting_point
290
291// returns the value of the objective function
292bool IpoptProblem::eval_f(Index n, const Number* x, bool new_x, Number& obj_value){
293 try{
294 obj_value = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 )[ 0];
295 }
296 catch(const ErrorClass& eclass){
297#ifdef DEBUG
298 cout << "error in OSIpoptSolver, line 296:\n" << eclass.errormsg << endl;
299#endif
300 ipoptErrorMsg = eclass.errormsg;
301 throw;
302 }
303 if( CommonUtil::ISOSNAN( (double)obj_value) ) return false;
304 return true;
305}
306
307bool IpoptProblem::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f){
308 int i;
309 double *objGrad;
310 try{
311 //objGrad = osinstance->calculateAllObjectiveFunctionGradients( const_cast<double*>(x), NULL, NULL, new_x, 1)[ 0];
312 objGrad = osinstance->calculateObjectiveFunctionGradient( const_cast<double*>(x), NULL, NULL, -1, new_x, 1);
313 }
314 catch(const ErrorClass& eclass){
315#ifdef DEBUG
316 cout << "error in OSIpoptSolver, line 314:\n" << eclass.errormsg << endl;
317#endif
318 ipoptErrorMsg = eclass.errormsg;
319 throw;
320 }
321 for(i = 0; i < n; i++){
322 grad_f[ i] = objGrad[ i];
323 }
324 return true;
325}//eval_grad_f
326
327// return the value of the constraints: g(x)
328bool IpoptProblem::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) {
329 try{
330 double *conVals = osinstance->calculateAllConstraintFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 );
331 int i;
332 for(i = 0; i < m; i++){
333 if( CommonUtil::ISOSNAN( (double)conVals[ i] ) ) return false;
334 g[i] = conVals[ i] ;
335 }
336 return true;
337 }
338 catch(const ErrorClass& eclass){
339#ifdef DEBUG
340 cout << "error in OSIpoptSolver, line 338:\n" << eclass.errormsg << endl;
341#endif
342 ipoptErrorMsg = eclass.errormsg;
343 throw;
344 }
345}//eval_g
346
347
348// return the structure or values of the jacobian
349bool IpoptProblem::eval_jac_g(Index n, const Number* x, bool new_x,
350 Index m, Index nele_jac, Index* iRow, Index *jCol,
351 Number* values){
352 SparseJacobianMatrix *sparseJacobian;
353 if (values == NULL) {
354 // return the values of the jacobian of the constraints
355 //cout << "n: " << n << endl;
356 //cout << "m: " << m << endl;
357 //cout << "nele_jac: " << nele_jac << endl;
358 // return the structure of the jacobian
359 try{
360 sparseJacobian = osinstance->getJacobianSparsityPattern();
361 }
362 catch(const ErrorClass& eclass){
363#ifdef DEBUG
364 cout << "error in OSIpoptSolver, line 362:\n" << eclass.errormsg << endl;
365#endif
366 ipoptErrorMsg = eclass.errormsg;
367 throw;
368 }
369 int i = 0;
370 int k, idx;
371 for(idx = 0; idx < m; idx++){
372 for(k = *(sparseJacobian->starts + idx); k < *(sparseJacobian->starts + idx + 1); k++){
373 iRow[i] = idx;
374 jCol[i] = *(sparseJacobian->indexes + k);
375 //cout << "ROW IDX !!!!!!!!!!!!!!!!!!!!!!!!!!!" << iRow[i] << endl;
376 //cout << "COL IDX !!!!!!!!!!!!!!!!!!!!!!!!!!!" << jCol[i] << endl;
377 i++;
378 }
379 }
380 }
381 else {
382 //std::cout << "EVALUATING JACOBIAN" << std::endl;
383 try{
384 sparseJacobian = osinstance->calculateAllConstraintFunctionGradients( const_cast<double*>(x), NULL, NULL, new_x, 1);
385 }
386 catch(const ErrorClass& eclass){
387#ifdef DEBUG
388 cout << "error in OSIpoptSolver, line 386:\n" << eclass.errormsg << endl;
389#endif
390 ipoptErrorMsg = eclass.errormsg;
391 throw;
392 }
393 //osinstance->getIterateResults( (double*)x, 0.0, NULL, -1, new_x, 1);
394 for(int i = 0; i < nele_jac; i++){
395 values[ i] = sparseJacobian->values[i];
396 //values[ i] = osinstance->m_mdJacValue[ i];
397 //cout << "values[i]:!!!!!!!!!!!! " << values[ i] << endl;
398 //cout << "m_mdJacValue[ i]:!!!!!!!!!!!! " << osinstance->m_mdJacValue[ i] << endl;
399 }
400 }
401 return true;
402}//eval_jac_g
403
404//return the structure or values of the hessian
405bool IpoptProblem::eval_h(Index n, const Number* x, bool new_x,
406 Number obj_factor, Index m, const Number* lambda,
407 bool new_lambda, Index nele_hess, Index* iRow,
408 Index* jCol, Number* values){
409
411 SparseHessianMatrix *sparseHessian;
412
413 int i;
414 if (values == NULL) {
415 // return the structure. This is a symmetric matrix, fill the lower left triangle only.
416 //cout << "get structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
417 try{
419 }
420 catch(const ErrorClass& eclass){
421#ifdef DEBUG
422 cout << "error in OSIpoptSolver, line 420:\n" << eclass.errormsg << endl;
423#endif
424 ipoptErrorMsg = eclass.errormsg;
425 throw;
426 }
427 //cout << "got structure of HESSIAN !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
428 for(i = 0; i < nele_hess; i++){
429 iRow[i] = *(sparseHessian->hessColIdx + i);
430 jCol[i] = *(sparseHessian->hessRowIdx + i);
431 //cout << "ROW HESS IDX !!!!!!!!!!!!!!!!!!!!!!!!!!!" << iRow[i] << endl;
432 //cout << "COL HESS IDX !!!!!!!!!!!!!!!!!!!!!!!!!!!" << jCol[i] << endl;
433 }
434 }
435 else {
436 //std::cout << "EVALUATING HESSIAN" << std::endl;
437 // return the values. This is a symmetric matrix, fill the lower left triangle only
438 double* objMultipliers = new double[1];
439 objMultipliers[0] = obj_factor;
440 try{
441 sparseHessian = osinstance->calculateLagrangianHessian( const_cast<double*>(x), objMultipliers, (double*)lambda , new_x, 2);
442 delete[] objMultipliers;
443 }
444 catch(const ErrorClass& eclass){
445#ifdef DEBUG
446 cout << "error in OSIpoptSolver, line 444:\n" << eclass.errormsg << endl;
447#endif
448 ipoptErrorMsg = eclass.errormsg;
449 delete[] objMultipliers;
450 throw;
451 }
452 for(i = 0; i < nele_hess; i++){
453 values[ i] = *(sparseHessian->hessValues + i);
454 }
455 }
457 return true;
458}//eval_h
459
460bool IpoptProblem::get_scaling_parameters(Number& obj_scaling,
461 bool& use_x_scaling, Index n,
462 Number* x_scaling,
463 bool& use_g_scaling, Index m,
464 Number* g_scaling){
465 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
466 obj_scaling = -1;
467 }
468 else obj_scaling = 1;
469 use_x_scaling = false;
470 use_g_scaling = false;
471 return true;
472}//get_scaling_parameters
473
474void IpoptProblem::finalize_solution(SolverReturn status,
475 Index n, const Number* x, const Number* z_L, const Number* z_U,
476 Index m, const Number* g, const Number* lambda,
477 Number obj_value,
478 const IpoptData* ip_data,
479 IpoptCalculatedQuantities* ip_cq)
480{
481 // here is where we would store the solution to variables, or write to a file, etc
482 // so we could use the solution.
483 // For this example, we write the solution to the console
484 OSrLWriter *osrlwriter ;
485 osrlwriter = new OSrLWriter();
486#ifdef DEBUG
487 printf("\n\nSolution of the primal variables, x\n");
488 for (Index i=0; i<n; i++) {
489 printf("x[%d] = %e\n", i, x[i]);
490 }
491
492 printf("\n\nSolution of the bound multipliers, z_L and z_U\n");
493 for (Index i=0; i<n; i++) {
494 printf("z_L[%d] = %e\n", i, z_L[i]);
495 }
496 for (Index i=0; i<n; i++) {
497 printf("z_U[%d] = %e\n", i, z_U[i]);
498 }
499#endif
500 printf("\nObjective value f(x*) = %e\n", obj_value);
501
502 int solIdx = 0;
503 int numberOfOtherVariableResults;
504 int otherIdx;
505 ostringstream outStr;
506
507 std::string *rcost = NULL;
508 double* mdObjValues = new double[1];
509 std::string message = "Ipopt solver finishes to the end.";
510 std::string solutionDescription = "";
511
512 try{
513 // resultHeader infomration
514 if(osresult->setServiceName( "Ipopt solver service") != true)
515 throw ErrorClass("OSResult error: setServiceName");
517 throw ErrorClass("OSResult error: setInstanceName");
518
519 //if(osresult->setJobID( osoption->jobID) != true)
520 // throw ErrorClass("OSResult error: setJobID");
521
522 // set basic problem parameters
524 throw ErrorClass("OSResult error: setVariableNumer");
525 if(osresult->setObjectiveNumber( 1) != true)
526 throw ErrorClass("OSResult error: setObjectiveNumber");
528 throw ErrorClass("OSResult error: setConstraintNumber");
529 if(osresult->setSolutionNumber( 1) != true)
530 throw ErrorClass("OSResult error: setSolutionNumer");
531
532
533 if(osresult->setGeneralMessage( message) != true)
534 throw ErrorClass("OSResult error: setGeneralMessage");
535
536 switch( status){
537 case SUCCESS:
538 solutionDescription = "SUCCESS[IPOPT]: Algorithm terminated successfully at a locally optimal point, satisfying the convergence tolerances.";
539 osresult->setSolutionStatus(solIdx, "locallyOptimal", solutionDescription);
540 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x), osinstance->getVariableNumber());
541 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda), osinstance->getConstraintNumber());
542 mdObjValues[0] = obj_value;
543 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
544
545
546 // set other
547
548 numberOfOtherVariableResults = 2;
549 osresult->setNumberOfOtherVariableResults(solIdx, numberOfOtherVariableResults);
550
551
552 rcost = new std::string[ osinstance->getVariableNumber()];
553
554
555 for (Index i = 0; i < n; i++) {
556 rcost[ i] = os_dtoa_format( z_L[i]);
557 }
558 otherIdx = 0;
559 osresult->setAnOtherVariableResult(solIdx, otherIdx, "varL", "Lagrange Multiplier on the Variable Lower Bound", rcost, osinstance->getVariableNumber());
560
561
562 for (Index i = 0; i < n; i++) {
563 rcost[ i] = os_dtoa_format( z_U[i]);
564 }
565 otherIdx = 1;
566 osresult->setAnOtherVariableResult(solIdx, otherIdx, "varU", "Lagrange Multiplier on the Variable Upper Bound", rcost, osinstance->getVariableNumber());
567
568 delete[] rcost;
569
570 // done with other
571
572
573 break;
574 case MAXITER_EXCEEDED:
575 solutionDescription = "MAXITER_EXCEEDED[IPOPT]: Maximum number of iterations exceeded.";
576 osresult->setSolutionStatus(solIdx, "stoppedByLimit", solutionDescription);
577 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x), osinstance->getVariableNumber());
578 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda), osinstance->getConstraintNumber());
579 mdObjValues[0] = obj_value;
580 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
581 break;
582 case STOP_AT_TINY_STEP:
583 solutionDescription = "STOP_AT_TINY_STEP[IPOPT]: Algorithm proceeds with very little progress.";
584 osresult->setSolutionStatus(solIdx, "stoppedByLimit", solutionDescription);
585 osresult->setPrimalVariableValues(solIdx, const_cast<double*>( x), osinstance->getVariableNumber());
586 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda), osinstance->getConstraintNumber());
587 mdObjValues[0] = obj_value;
588 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
589 break;
590 case STOP_AT_ACCEPTABLE_POINT:
591 solutionDescription = "STOP_AT_ACCEPTABLE_POINT[IPOPT]: Algorithm stopped at a point that was converged, not to _desired_ tolerances, but to _acceptable_ tolerances";
592 osresult->setSolutionStatus(solIdx, "IpoptAccetable", solutionDescription);
593 osresult->setPrimalVariableValues(solIdx, const_cast<double*>(x), osinstance->getVariableNumber());
594 osresult->setDualVariableValues(solIdx, const_cast<double*>( lambda), osinstance->getConstraintNumber());
595 mdObjValues[0] = obj_value;
596 osresult->setObjectiveValues(solIdx, mdObjValues, 1);
597 break;
598 case LOCAL_INFEASIBILITY:
599 solutionDescription = "LOCAL_INFEASIBILITY[IPOPT]: Algorithm converged to a point of local infeasibility. Problem may be infeasible.";
600 osresult->setSolutionStatus(solIdx, "infeasible", solutionDescription);
601 break;
602 case USER_REQUESTED_STOP:
603 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.";
604 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
605 break;
606 case DIVERGING_ITERATES:
607 solutionDescription = "DIVERGING_ITERATES[IPOPT]: It seems that the iterates diverge.";
608 osresult->setSolutionStatus(solIdx, "unbounded", solutionDescription);
609 break;
610 case RESTORATION_FAILURE:
611 solutionDescription = "RESTORATION_FAILURE[IPOPT]: Restoration phase failed, algorithm doesn't know how to proceed.";
612 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
613 break;
614 case ERROR_IN_STEP_COMPUTATION:
615 solutionDescription = "ERROR_IN_STEP_COMPUTATION[IPOPT]: An unrecoverable error occurred while IPOPT tried to compute the search direction.";
616 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
617 break;
618 case INVALID_NUMBER_DETECTED:
619 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.";
620 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
621 break;
622 case INTERNAL_ERROR:
623 solutionDescription = "INTERNAL_ERROR[IPOPT]: An unknown internal error occurred. Please contact the IPOPT authors through the mailing list.";
624 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
625 break;
626 default:
627 solutionDescription = "OTHER[IPOPT]: other unknown solution status from Ipopt solver";
628 osresult->setSolutionStatus(solIdx, "other", solutionDescription);
629 }
631 delete osrlwriter;
632 delete[] mdObjValues;
633 osrlwriter = NULL;
634
635 }
636 catch(const ErrorClass& eclass){
637#ifdef DEBUG
638 cout << "error in OSIpoptSolver, line 636:\n" << eclass.errormsg << endl;
639#endif
642 std::string osrl = osrlwriter->writeOSrL( osresult);
643 delete osrlwriter;
644 osrlwriter = NULL;
645 throw ErrorClass( osrl) ;
646 delete[] mdObjValues;
647 mdObjValues = NULL;
648 }
650}
651
652
654 try{
655 /* set the default options */
656 this->bSetSolverOptions = true;
657 app->Options()->SetNumericValue("tol", 1e-9);
658 app->Options()->SetIntegerValue("print_level", 0);
659 app->Options()->SetIntegerValue("max_iter", 20000);
660 app->Options()->SetStringValue("mu_strategy", "adaptive");
661 app->Options()->SetStringValue("output_file", "ipopt.out");
662 app->Options()->SetStringValue("check_derivatives_for_naninf", "yes");
663 /* end of the default options */
664
665 if(osoption == NULL && osol.length() > 0)
666 {
667 m_osolreader = new OSoLReader();
669 }
670
671 if( osoption != NULL && osoption->getNumberOfSolverOptions() > 0 ){
672 std::cout << "number of solver options " << osoption->getNumberOfSolverOptions() << std::endl;
673 std::vector<SolverOption*> optionsVector;
674 optionsVector = osoption->getSolverOptions( "ipopt");
675 char *pEnd;
676 int i;
677 int num_ipopt_options = optionsVector.size();
678 for(i = 0; i < num_ipopt_options; i++){
679 std::cout << "ipopt solver option " << optionsVector[ i]->name << std::endl;
680 if(optionsVector[ i]->type == "numeric" ){
681 std::cout << "FOUND A NUMERIC OPTION " << os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) << std::endl;
682 app->Options()->SetNumericValue(optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
683 }
684 else if(optionsVector[ i]->type == "integer" ){
685 std::cout << "FOUND AN INTEGER OPTION " << atoi( optionsVector[ i]->value.c_str() ) << std::endl;
686 app->Options()->SetIntegerValue(optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ) );
687 }
688 else if(optionsVector[ i]->type == "string" ){
689 std::cout << "FOUND A STRING OPTION " << optionsVector[ i]->value.c_str() << std::endl;
690 app->Options()->SetStringValue(optionsVector[ i]->name, optionsVector[ i]->value);
691 }
692 }
693 }
694 }
695 catch(const ErrorClass& eclass){
696#ifdef DEBUG
697 cout << "error in OSIpoptSolver, line 695:\n" << eclass.errormsg << endl;
698#endif
699 std::cout << "THERE IS AN ERROR" << std::endl;
703 throw ErrorClass( osrl) ;
704 }
705}//end setSolverOptions()
706
707
708
710 try{
711
712 if(osil.length() == 0 && osinstance == NULL) throw ErrorClass("there is no instance");
713 if(osinstance == NULL){
714 m_osilreader = new OSiLReader();
716 }
717 // Create a new instance of your nlp
719 app = new IpoptApplication();
720 this->bCallbuildSolverInstance = true;
721 }
722 catch(const ErrorClass& eclass){
723#ifdef DEBUG
724 cout << "error in OSIpoptSolver, line 722:\n" << eclass.errormsg << endl;
725#endif
726 std::cout << "THERE IS AN ERROR" << std::endl;
730 throw ErrorClass( osrl) ;
731 }
732}//end buildSolverInstance()
733
734
735
736void IpoptSolver::solve() throw (ErrorClass) {
737 if( this->bCallbuildSolverInstance == false) buildSolverInstance();
738 if( this->bSetSolverOptions == false) setSolverOptions();
739 try{
740 clock_t start, finish;
741 double duration;
742 start = clock();
743 //OSiLWriter osilwriter;
744 //cout << osilwriter.writeOSiL( osinstance) << endl;
745 if(osinstance->getVariableNumber() <= 0)throw ErrorClass("Ipopt requires decision variables");
746 finish = clock();
747 duration = (double) (finish - start) / CLOCKS_PER_SEC;
748 cout << "Parsing took (seconds): " << duration << endl;
749 //dataEchoCheck();
750 /***************now the ipopt invokation*********************/
751 // see if we have a linear program
752 if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Ipopt NEEDS AN OBJECTIVE FUNCTION");
754 app->Options()->SetStringValue("hessian_approximation", "limited-memory");
755 // if it is a max problem call scaling and set to -1
756 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") != 0){
757 app->Options()->SetStringValue("nlp_scaling_method", "user-scaling");
758 }
759 // Intialize the IpoptApplication and process the options
760 std::cout << "Call Ipopt Initialize" << std::endl;
761 app->Initialize();
762 std::cout << "Finished Ipopt Initialize" << std::endl;
763 //nlp->osinstance = this->osinstance;
764 // Ask Ipopt to solve the problem
765 std::cout << "Call Ipopt Optimize" << std::endl;
766 ApplicationReturnStatus status = app->OptimizeTNLP( nlp);
767 std::cout << "Finish Ipopt Optimize" << std::endl;
769 std::cout << "Finish writing the osrl" << std::endl;
770 //if (status != Solve_Succeeded) {
771 if (status < -2) {
772 throw ErrorClass("Ipopt FAILED TO SOLVE THE PROBLEM: " + ipoptErrorMsg);
773 }
774 }
775 catch(const ErrorClass& eclass){
776#ifdef DEBUG
777 cout << "error in OSIpoptSolver, line 775:\n" << eclass.errormsg << endl;
778#endif
782 throw ErrorClass( osrl) ;
783 }
784}//solve
785
786
788
789 int i;
790
791 // print out problem parameters
792 cout << "This is problem: " << osinstance->getInstanceName() << endl;
793 cout << "The problem source is: " << osinstance->getInstanceSource() << endl;
794 cout << "The problem description is: " << osinstance->getInstanceDescription() << endl;
795 cout << "number of variables = " << osinstance->getVariableNumber() << endl;
796 cout << "number of Rows = " << osinstance->getConstraintNumber() << endl;
797
798 // print out the variable information
799 if(osinstance->getVariableNumber() > 0){
800 for(i = 0; i < osinstance->getVariableNumber(); i++){
801 if(osinstance->getVariableNames() != NULL) cout << "variable Names " << osinstance->getVariableNames()[ i] << endl;
802 if(osinstance->getVariableTypes() != NULL) cout << "variable Types " << osinstance->getVariableTypes()[ i] << endl;
803 if(osinstance->getVariableLowerBounds() != NULL) cout << "variable Lower Bounds " << osinstance->getVariableLowerBounds()[ i] << endl;
804 if(osinstance->getVariableUpperBounds() != NULL) cout << "variable Upper Bounds " << osinstance->getVariableUpperBounds()[i] << endl;
805 }
806 }
807
808 // print out objective function information
810 if( osinstance->getObjectiveMaxOrMins()[0] == "min") cout << "problem is a minimization" << endl;
811 else cout << "problem is a maximization" << endl;
812 for(i = 0; i < osinstance->getVariableNumber(); i++){
813 cout << "OBJ COEFFICIENT = " << osinstance->getDenseObjectiveCoefficients()[0][i] << endl;
814 }
815 }
816 // print out constraint information
818 for(i = 0; i < osinstance->getConstraintNumber(); i++){
819 if(osinstance->getConstraintNames() != NULL) cout << "row name = " << osinstance->getConstraintNames()[i] << endl;
820 if(osinstance->getConstraintLowerBounds() != NULL) cout << "row lower bound = " << osinstance->getConstraintLowerBounds()[i] << endl;
821 if(osinstance->getConstraintUpperBounds() != NULL) cout << "row upper bound = " << osinstance->getConstraintUpperBounds()[i] << endl;
822 }
823 }
824
825 // print out linear constraint data
826 cout << endl;
827 cout << "number of nonzeros = " << osinstance->getLinearConstraintCoefficientNumber() << endl;
828 for(i = 0; i <= osinstance->getVariableNumber(); i++){
829 cout << "Start Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->starts[ i] << endl;
830 }
831 cout << endl;
832 for(i = 0; i < osinstance->getLinearConstraintCoefficientNumber(); i++){
833 cout << "Index Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->indexes[i] << endl;
834 cout << "Nonzero Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->values[i] << endl;
835 }
836
837 // print out quadratic data
838 cout << "number of qterms = " << osinstance->getNumberOfQuadraticTerms() << endl;
839 for(int i = 0; i < osinstance->getNumberOfQuadraticTerms(); i++){
840 cout << "Row Index = " << osinstance->getQuadraticTerms()->rowIndexes[i] << endl;
841 cout << "Var Index 1 = " << osinstance->getQuadraticTerms()->varOneIndexes[ i] << endl;
842 cout << "Var Index 2 = " << osinstance->getQuadraticTerms()->varTwoIndexes[ i] << endl;
843 cout << "Coefficient = " << osinstance->getQuadraticTerms()->coefficients[ i] << endl;
844 }
845} // end dataEchoCheck
846
847
848IpoptProblem::IpoptProblem(OSInstance *osinstance_, OSOption *osoption_, OSResult *osresult_) {
849 osinstance = osinstance_;
850 osoption = osoption_;
851 osresult = osresult_;
852}
853
855
856}
857
858
859
std::string os_dtoa_format(double x)
double os_strtod(const char *s00, char **se)
Definition OSdtoa.cpp:2541
std::string osol
osol holds the options for the solver
bool bSetSolverOptions
bSetSolverOptions is set to true if setSolverOptions has been called, false otherwise
std::string osrl
osrl holds the solution or result of the model
OSInstance * osinstance
osinstance holds the problem instance in-memory as an OSInstance object
bool bCallbuildSolverInstance
bCallbuildSolverInstance is set to true if buildSolverService has been called
std::string osil
osil holds the problem instance as a std::string
OSOption * osoption
osoption holds the solver options in-memory as an OSOption object
OSResult * osresult
osresult holds the solution or result of the model in-memory as an OSResult object
used for throwing exceptions.
std::string errormsg
errormsg is the error that is causing the exception to be thrown
the InitVarValue class.
Definition OSOption.h:1160
double value
initial value
Definition OSOption.h:1170
int idx
variable index
Definition OSOption.h:1164
Variables * variables
variables is a pointer to a Variables object
Objectives * objectives
objectives is a pointer to a Objectives object
virtual bool get_starting_point(Ipopt::Index n, bool init_x, Ipopt::Number *x, bool init_z, Ipopt::Number *z_L, Ipopt::Number *z_U, Ipopt::Index m, bool init_lambda, Ipopt::Number *lambda)
Method to return the starting point for the algorithm.
OSOption * osoption
OSInstance * osinstance
virtual bool get_scaling_parameters(Ipopt::Number &obj_scaling, bool &use_x_scaling, Ipopt::Index n, Ipopt::Number *x_scaling, bool &use_g_scaling, Ipopt::Index m, Ipopt::Number *g_scaling)
std::string * ipoptErrorMsg
virtual bool eval_grad_f(Ipopt::Index n, const Ipopt::Number *x, bool new_x, Ipopt::Number *grad_f)
Method to return the gradient of the objective.
IpoptProblem(OSInstance *osinstance_, OSOption *osoption_, OSResult *osresult_, std::string *ipoptErrorMsg_)
the IpoptProblemclass constructor
virtual bool eval_g(Ipopt::Index n, const Ipopt::Number *x, bool new_x, Ipopt::Index m, Ipopt::Number *g)
Method to return the constraint residuals.
virtual ~IpoptProblem()
the IpoptProblem class destructor
virtual bool get_nlp_info(Ipopt::Index &n, Ipopt::Index &m, Ipopt::Index &nnz_jac_g, Ipopt::Index &nnz_h_lag, IndexStyleEnum &index_style)
IPOpt specific methods for defining the nlp problem.
virtual bool eval_jac_g(Ipopt::Index n, const Ipopt::Number *x, bool new_x, Ipopt::Index m, Ipopt::Index nele_jac, Ipopt::Index *iRow, Ipopt::Index *jCol, Ipopt::Number *values)
Method to return: 1) The structure of the jacobian (if "values" is NULL) 2) The values of the jacobia...
OSResult * osresult
virtual bool eval_f(Ipopt::Index n, const Ipopt::Number *x, bool new_x, Ipopt::Number &obj_value)
Method to return the objective value.
virtual bool get_bounds_info(Ipopt::Index n, Ipopt::Number *x_l, Ipopt::Number *x_u, Ipopt::Index m, Ipopt::Number *g_l, Ipopt::Number *g_u)
Method to return the bounds for my problem.
virtual bool eval_h(Ipopt::Index n, const Ipopt::Number *x, bool new_x, Ipopt::Number obj_factor, Ipopt::Index m, const Ipopt::Number *lambda, bool new_lambda, Ipopt::Index nele_hess, Ipopt::Index *iRow, Ipopt::Index *jCol, Ipopt::Number *values)
Method to return: 1) The structure of the hessian of the lagrangian (if "values" is NULL) 2) The valu...
virtual void finalize_solution(Ipopt::SolverReturn status, Ipopt::Index n, const Ipopt::Number *x, const Ipopt::Number *z_L, const Ipopt::Number *z_U, Ipopt::Index m, const Ipopt::Number *g, const Ipopt::Number *lambda, Ipopt::Number obj_value, const Ipopt::IpoptData *ip_data, Ipopt::IpoptCalculatedQuantities *ip_cq)
This method is called when the algorithm is complete so the TNLP can store/write the solution.
std::string * ipoptErrorMsg
OSoLReader * m_osolreader
m_osolreader is an OSoLReader object used to create an osoption from an osol string if needed
OSrLWriter * osrlwriter
OSiLReader * m_osilreader
m_osilreader is an OSiLReader object used to create an osinstance from an osil string if needed
virtual void setSolverOptions()
The implementation of the virtual functions.
Ipopt::SmartPtr< Ipopt::TNLP > nlp
Ipopt::SmartPtr< Ipopt::IpoptApplication > app
virtual void solve()
solve results in an instance being read into the Ipopt data structures and optimize
IpoptSolver()
the IpoptSolver class constructor
~IpoptSolver()
the IpoptSolver class destructor
void dataEchoCheck()
use this for debugging, print out the instance that the solver thinks it has and compare this with th...
virtual void buildSolverInstance()
The implementation of the virtual functions.
The in-memory representation of an OSiL instance..
SparseJacobianMatrix * calculateAllConstraintFunctionGradients(double *x, double *objLambda, double *conLambda, bool new_x, int highestOrder)
Calculate the gradient of all constraint functions.
double * getConstraintLowerBounds()
Get constraint lower bounds.
int getNumberOfQuadraticTerms()
Get the number of specified (usually nonzero) qTerms in the quadratic coefficients.
double * getVariableUpperBounds()
Get variable upper bounds.
SparseJacobianMatrix * getJacobianSparsityPattern()
bool bUseExpTreeForFunEval
bUseExpTreeForFunEval is set to true if you wish to use the OS Expression Tree for function evaluatio...
std::string getInstanceDescription()
Get instance description.
std::string getInstanceSource()
Get instance source.
int getConstraintNumber()
Get number of constraints.
double * calculateAllConstraintFunctionValues(double *x, double *objLambda, double *conLambda, bool new_x, int highestOrder)
Calculate all of the constraint function values.
int getLinearConstraintCoefficientNumber()
Get number of specified (usually nonzero) linear constraint coefficient values.
char * getVariableTypes()
Get variable initial values.
double * calculateAllObjectiveFunctionValues(double *x, double *objLambda, double *conLambda, bool new_x, int highestOrder)
Calculate all of the objective function values.
SparseMatrix * getLinearConstraintCoefficientsInColumnMajor()
Get linear constraint coefficients in column major.
double * calculateObjectiveFunctionGradient(double *x, double *objLambda, double *conLambda, int objIdx, bool new_x, int highestOrder)
Calculate the gradient of the objective function indexed by objIdx.
double ** getDenseObjectiveCoefficients()
getDenseObjectiveCoefficients.
bool initForAlgDiff()
This should be called by nonlinear solvers using callback functions.
InstanceData * instanceData
A pointer to an InstanceData object.
int getNumberOfNonlinearExpressions()
Get number of nonlinear expressions.
SparseHessianMatrix * getLagrangianHessianSparsityPattern()
QuadraticTerms * getQuadraticTerms()
Get all the quadratic terms in the instance.
double * getVariableLowerBounds()
Get variable lower bounds.
int getVariableNumber()
Get number of variables.
std::string * getVariableNames()
Get variable names.
std::string getInstanceName()
Get instance name.
SparseHessianMatrix * calculateLagrangianHessian(double *x, double *objLambda, double *conLambda, bool new_x, int highestOrder)
Calculate the Hessian of the Lagrangian Expression Tree This method will build the CppAD expression t...
std::string * getConstraintNames()
Get constraint names.
std::string * getObjectiveMaxOrMins()
Get objective maxOrMins.
double * getConstraintUpperBounds()
Get constraint upper bounds.
int getObjectiveNumber()
Get number of objectives.
The Option Class.
Definition OSOption.h:3565
InitVarValue ** getInitVarValuesSparse()
Get the initial values associated with the variables in sparse form.
int getNumberOfInitVarValues()
Get the number of initial variable values.
std::vector< SolverOption * > getSolverOptions(std::string solver_name)
Get the options associated with a given solver.
int getNumberOfSolverOptions()
Get the number of solver options.
The Result Class.
Definition OSResult.h:2549
bool setGeneralMessage(std::string message)
Set the general message.
bool setAnOtherVariableResult(int solIdx, int otherIdx, std::string name, std::string description, int *indexes, std::string *s, int n)
Set the [i]th optimization solution's other (non-standard/solver specific)variable-related results,...
bool setPrimalVariableValues(int solIdx, double *x, int n)
Set the [i]th optimization solution's primal variable values, where i equals the given solution index...
bool setDualVariableValues(int solIdx, double *lbValues, double *ubValues, int n)
Set the [i]th optimization solution's dual variable values, where i equals the given solution index.
bool setSolutionNumber(int number)
set the number of solutions.
bool setInstanceName(std::string instanceName)
Set instance name.
bool setNumberOfOtherVariableResults(int solIdx, int numberOfOtherVariableResults)
Set the [i]th optimization solution's other (non-standard/solver specific) variable-related results,...
bool setGeneralStatusType(std::string type)
Set the general status type, which can be: success, error, warning.
bool setObjectiveNumber(int objectiveNumber)
Set the objective number.
bool setObjectiveValues(int solIdx, double *objectiveValues, int n)
Set the [i]th optimization solution's objective values, where i equals the given solution index.
bool setServiceName(std::string serviceName)
Set service name.
bool setVariableNumber(int variableNumber)
Set the variable number.
bool setSolutionStatus(int solIdx, std::string type, std::string description)
Set the [i]th optimization solution status, where i equals the given solution index.
bool setConstraintNumber(int constraintNumber)
Set the constraint number.
Used to read an OSiL string.
Definition OSiLReader.h:38
OSInstance * readOSiL(const std::string &osil)
parse the OSiL model instance.
Used to read an OSoL string.
Definition OSoLReader.h:38
OSOption * readOSoL(const std::string &osol)
parse the OSoL solver options.
Take an OSResult object and write a string that validates against OSrL.
Definition OSrLWriter.h:31
std::string writeOSrL(OSResult *theosresult)
create an osrl string from an OSResult object
std::string maxOrMin
declare the objective function to be a max or a min
Definition OSInstance.h:157
int numberOfObjectives
numberOfObjectives is the number of objective functions in the instance
Definition OSInstance.h:201
Objective ** obj
coef is pointer to an array of ObjCoef object pointers
Definition OSInstance.h:205
int * varTwoIndexes
varTwoIndexes holds an integer array of the second variable indexes of all the quadratic terms.
Definition OSGeneral.h:450
int * rowIndexes
rowIndexes holds an integer array of row indexes of all the quadratic terms.
Definition OSGeneral.h:440
double * coefficients
coefficients holds a double array all the quadratic term coefficients.
Definition OSGeneral.h:455
int * varOneIndexes
varOneIndexes holds an integer array of the first variable indexes of all the quadratic terms.
Definition OSGeneral.h:445
The in-memory representation of a SparseHessianMatrix..
Definition OSGeneral.h:377
int * hessRowIdx
hessRowIdx is an integer array of row indices in the range 0, ..., n - 1.
Definition OSGeneral.h:394
int hessDimension
hessDimension is the number of nonzeros in each array.
Definition OSGeneral.h:389
double * hessValues
hessValues is a double array of the Hessian values.
Definition OSGeneral.h:404
int * hessColIdx
hessColIdx is an integer array of column indices in the range 0, ..., n - 1.
Definition OSGeneral.h:399
a sparse Jacobian matrix data structure
Definition OSGeneral.h:301
int * indexes
indexes holds an integer array of variable indices.
Definition OSGeneral.h:335
int valueSize
valueSize is the dimension of the values array
Definition OSGeneral.h:318
int * starts
starts holds an integer array of start elements, each start element points to the start of partials f...
Definition OSGeneral.h:324
double * values
values holds a double array of nonzero partial derivatives
Definition OSGeneral.h:340
int * indexes
indexes holds an integer array of rowIdx (or colIdx) elements in coefMatrix (AMatrix).
Definition OSGeneral.h:258
int * starts
starts holds an integer array of start elements in coefMatrix (AMatrix), which points to the start of...
Definition OSGeneral.h:252
double * values
values holds a double array of value elements in coefMatrix (AMatrix), which contains nonzero element...
Definition OSGeneral.h:264
double ub
ub corresponds to the optional attribute that holds the variable upper bound.
Definition OSInstance.h:61
double lb
lb corresponds to the optional attribute that holds the variable lower bound.
Definition OSInstance.h:56
Variable ** var
Here we define a pointer to an array of var pointers.
Definition OSInstance.h:97
#define OSDBL_MAX