24#include "OSCommonUtil.h"
30using std::ostringstream;
44 cout <<
"inside IpoptSolver destructor" << endl;
57 cout <<
"leaving IpoptSolver destructor" << endl;
63 Index& nnz_h_lag, IndexStyleEnum& index_style)
65 if(
osinstance->getObjectiveNumber() <= 0)
throw ErrorClass(
"Ipopt NEEDS AN OBJECTIVE FUNCTION");
71 cout <<
"number variables !!!!!!!!!!!!!!!!!!!!!!!!!!!" << n << endl;
72 cout <<
"number constraints !!!!!!!!!!!!!!!!!!!!!!!!!!!" << m << endl;
77 catch(
const ErrorClass& eclass){
79 cout <<
"error in OSIpoptSolver, line 78:\n" << eclass.
errormsg << endl;
87 SparseJacobianMatrix *sparseJacobian = NULL;
89 sparseJacobian =
osinstance->getJacobianSparsityPattern();
91 catch(
const ErrorClass& eclass){
93 cout <<
"error in OSIpoptSolver, line 91:\n" << eclass.
errormsg << endl;
101 cout <<
"nnz_jac_g !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_jac_g << endl;
105 if( (
osinstance->getNumberOfNonlinearExpressions() == 0) && (
osinstance->getNumberOfQuadraticTerms() == 0) ) {
107 cout <<
"This is a linear program" << endl;
113 SparseHessianMatrix *sparseHessian =
osinstance->getLagrangianHessianSparsityPattern();
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;
123 index_style = TNLP::C_STYLE;
125 cout <<
"return from get_nlp_info (OSIpoptSolver.cpp)" << nnz_h_lag << endl;
135 Index m, Number* g_l, Number* g_u){
137 double * mdVarLB =
osinstance->getVariableLowerBounds();
140 double * mdVarUB =
osinstance->getVariableUpperBounds();
142 for(i = 0; i < n; i++){
143 x_l[ i] = mdVarLB[ i];
144 x_u[ i] = mdVarUB[ i];
154 double * mdConLB =
osinstance->getConstraintLowerBounds();
156 double * mdConUB =
osinstance->getConstraintUpperBounds();
158 for(
int i = 0; i < m; i++){
159 g_l[ i] = mdConLB[ i];
160 g_u[ i] = mdConUB[ i];
170 bool init_z, Number* z_L, Number* z_U, Index m,
bool init_lambda,
175 assert(init_x ==
true);
176 assert(init_z ==
false);
177 assert(init_lambda ==
false);
181 cout <<
"get initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
186 cout <<
"get number of initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
187 cout <<
"Is osoption = NULL? " << (
osoption == NULL) << endl;
191 m1 =
osoption->getNumberOfInitVarValues();
195 cout <<
"number of variables initialed: " << m1 << endl;
200 initialed =
new bool[n1];
202 cout <<
"number of variables in total: " << n1 << endl;
205 for(k = 0; k < n1; k++)
206 initialed[k] =
false;
211 cout <<
"get initial values " << endl;
214 InitVarValue** initVarVector =
osoption->getInitVarValuesSparse();
216 cout <<
"done " << endl;
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");
227 initval = initVarVector[k]->
value;
229 {
if (
osinstance->instanceData->variables->var[k]->lb > initval)
230 throw ErrorClass (
"Initial value outside of bounds");
234 {
if (
osinstance->instanceData->variables->var[k]->ub < initval)
235 throw ErrorClass (
"Initial value outside of bounds");
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");
243 x[initVarVector[k]->
idx] = initval;
244 initialed[initVarVector[k]->idx] =
true;
247 catch(
const ErrorClass& eclass)
248 { cout <<
"Error in IpoptProblem::get_starting_point (OSIpoptSolver.cpp, line 247)";
249 cout << endl << endl << endl;
253 double default_initval;
254 default_initval = 1.7171;
256 for(k = 0; k < n1; k++)
257 { cout <<
"verify component " << k << endl;
260 if (
osinstance->instanceData->variables->var[k]->lb <= default_initval)
261 x[k] = default_initval;
263 x[k] =
osinstance->instanceData->variables->var[k]->lb;
266 if (
osinstance->instanceData->variables->var[k]->ub >= default_initval)
267 x[k] = default_initval;
269 x[k] =
osinstance->instanceData->variables->var[k]->ub;
271 if ((
osinstance->instanceData->variables->var[k]->lb <= default_initval) &&
272 (
osinstance->instanceData->variables->var[k]->ub >= default_initval))
273 x[k] = default_initval;
275 if (
osinstance->instanceData->variables->var[k]->lb > default_initval)
276 x[k] =
osinstance->instanceData->variables->var[k]->lb;
278 x[k] =
osinstance->instanceData->variables->var[k]->ub;
281 for(i = 0; i < n1; i++){
282 std::cout <<
"INITIAL VALUE !!!!!!!!!!!!!!!!!!!! " << x[ i] << std::endl;
294 obj_value =
osinstance->calculateAllObjectiveFunctionValues(
const_cast<double*
>(x), NULL, NULL, new_x, 0 )[ 0];
296 catch(
const ErrorClass& eclass){
298 cout <<
"error in OSIpoptSolver, line 296:\n" << eclass.
errormsg << endl;
303 if( CommonUtil::ISOSNAN( (
double)obj_value) )
return false;
312 objGrad =
osinstance->calculateObjectiveFunctionGradient(
const_cast<double*
>(x), NULL, NULL, -1, new_x, 1);
314 catch(
const ErrorClass& eclass){
316 cout <<
"error in OSIpoptSolver, line 314:\n" << eclass.
errormsg << endl;
321 for(i = 0; i < n; i++){
322 grad_f[ i] = objGrad[ i];
330 double *conVals =
osinstance->calculateAllConstraintFunctionValues(
const_cast<double*
>(x), NULL, NULL, new_x, 0 );
332 for(i = 0; i < m; i++){
333 if( CommonUtil::ISOSNAN( (
double)conVals[ i] ) )
return false;
338 catch(
const ErrorClass& eclass){
340 cout <<
"error in OSIpoptSolver, line 338:\n" << eclass.
errormsg << endl;
350 Index m, Index nele_jac, Index* iRow, Index *jCol,
352 SparseJacobianMatrix *sparseJacobian;
353 if (values == NULL) {
360 sparseJacobian =
osinstance->getJacobianSparsityPattern();
362 catch(
const ErrorClass& eclass){
364 cout <<
"error in OSIpoptSolver, line 362:\n" << eclass.
errormsg << endl;
371 for(idx = 0; idx < m; idx++){
372 for(k = *(sparseJacobian->
starts + idx); k < *(sparseJacobian->
starts + idx + 1); k++){
374 jCol[i] = *(sparseJacobian->
indexes + k);
384 sparseJacobian =
osinstance->calculateAllConstraintFunctionGradients(
const_cast<double*
>(x), NULL, NULL, new_x, 1);
386 catch(
const ErrorClass& eclass){
388 cout <<
"error in OSIpoptSolver, line 386:\n" << eclass.
errormsg << endl;
394 for(
int i = 0; i < nele_jac; i++){
395 values[ i] = sparseJacobian->
values[i];
406 Number obj_factor, Index m,
const Number* lambda,
407 bool new_lambda, Index nele_hess, Index* iRow,
408 Index* jCol, Number* values){
411 SparseHessianMatrix *sparseHessian;
414 if (values == NULL) {
418 sparseHessian =
osinstance->getLagrangianHessianSparsityPattern( );
420 catch(
const ErrorClass& eclass){
422 cout <<
"error in OSIpoptSolver, line 420:\n" << eclass.
errormsg << endl;
428 for(i = 0; i < nele_hess; i++){
438 double* objMultipliers =
new double[1];
439 objMultipliers[0] = obj_factor;
441 sparseHessian =
osinstance->calculateLagrangianHessian(
const_cast<double*
>(x), objMultipliers, (
double*)lambda , new_x, 2);
442 delete[] objMultipliers;
444 catch(
const ErrorClass& eclass){
446 cout <<
"error in OSIpoptSolver, line 444:\n" << eclass.
errormsg << endl;
449 delete[] objMultipliers;
452 for(i = 0; i < nele_hess; i++){
453 values[ i] = *(sparseHessian->
hessValues + i);
461 bool& use_x_scaling, Index n,
463 bool& use_g_scaling, Index m,
465 if(
osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare(
"min") != 0){
468 else obj_scaling = 1;
469 use_x_scaling =
false;
470 use_g_scaling =
false;
475 Index n,
const Number* x,
const Number* z_L,
const Number* z_U,
476 Index m,
const Number* g,
const Number* lambda,
478 const IpoptData* ip_data,
479 IpoptCalculatedQuantities* ip_cq)
484 OSrLWriter *osrlwriter ;
485 osrlwriter =
new OSrLWriter();
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]);
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]);
496 for (Index i=0; i<n; i++) {
497 printf(
"z_U[%d] = %e\n", i, z_U[i]);
500 printf(
"\nObjective value f(x*) = %e\n", obj_value);
503 int numberOfOtherVariableResults;
505 ostringstream outStr;
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 =
"";
514 if(
osresult->setServiceName(
"Ipopt solver service") !=
true)
515 throw ErrorClass(
"OSResult error: setServiceName");
517 throw ErrorClass(
"OSResult error: setInstanceName");
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");
533 if(
osresult->setGeneralMessage( message) !=
true)
534 throw ErrorClass(
"OSResult error: setGeneralMessage");
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);
548 numberOfOtherVariableResults = 2;
549 osresult->setNumberOfOtherVariableResults(solIdx, numberOfOtherVariableResults);
552 rcost =
new std::string[
osinstance->getVariableNumber()];
555 for (Index i = 0; i < n; i++) {
559 osresult->setAnOtherVariableResult(solIdx, otherIdx,
"varL",
"Lagrange Multiplier on the Variable Lower Bound", rcost,
osinstance->getVariableNumber());
562 for (Index i = 0; i < n; i++) {
566 osresult->setAnOtherVariableResult(solIdx, otherIdx,
"varU",
"Lagrange Multiplier on the Variable Upper Bound", rcost,
osinstance->getVariableNumber());
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);
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);
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);
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);
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);
606 case DIVERGING_ITERATES:
607 solutionDescription =
"DIVERGING_ITERATES[IPOPT]: It seems that the iterates diverge.";
608 osresult->setSolutionStatus(solIdx,
"unbounded", solutionDescription);
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);
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);
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);
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);
627 solutionDescription =
"OTHER[IPOPT]: other unknown solution status from Ipopt solver";
628 osresult->setSolutionStatus(solIdx,
"other", solutionDescription);
630 osresult->setGeneralStatusType(
"normal");
632 delete[] mdObjValues;
636 catch(
const ErrorClass& eclass){
638 cout <<
"error in OSIpoptSolver, line 636:\n" << eclass.
errormsg << endl;
641 osresult->setGeneralStatusType(
"error");
645 throw ErrorClass( osrl) ;
646 delete[] mdObjValues;
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");
672 std::cout <<
"number of solver options " <<
osoption->getNumberOfSolverOptions() << std::endl;
673 std::vector<SolverOption*> optionsVector;
674 optionsVector =
osoption->getSolverOptions(
"ipopt");
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 ) );
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() ) );
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);
695 catch(
const ErrorClass& eclass){
697 cout <<
"error in OSIpoptSolver, line 695:\n" << eclass.
errormsg << endl;
699 std::cout <<
"THERE IS AN ERROR" << std::endl;
701 osresult->setGeneralStatusType(
"error");
703 throw ErrorClass(
osrl) ;
712 if(
osil.length() == 0 &&
osinstance == NULL)
throw ErrorClass(
"there is no instance");
719 app =
new IpoptApplication();
722 catch(
const ErrorClass& eclass){
724 cout <<
"error in OSIpoptSolver, line 722:\n" << eclass.
errormsg << endl;
726 std::cout <<
"THERE IS AN ERROR" << std::endl;
728 osresult->setGeneralStatusType(
"error");
730 throw ErrorClass(
osrl) ;
740 clock_t start, finish;
745 if(
osinstance->getVariableNumber() <= 0)
throw ErrorClass(
"Ipopt requires decision variables");
747 duration = (double) (finish - start) / CLOCKS_PER_SEC;
748 cout <<
"Parsing took (seconds): " << duration << endl;
752 if(
osinstance->getObjectiveNumber() <= 0)
throw ErrorClass(
"Ipopt NEEDS AN OBJECTIVE FUNCTION");
753 if( (
osinstance->getNumberOfNonlinearExpressions() == 0) && (
osinstance->getNumberOfQuadraticTerms() == 0) )
754 app->Options()->SetStringValue(
"hessian_approximation",
"limited-memory");
756 if(
osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare(
"min") != 0){
757 app->Options()->SetStringValue(
"nlp_scaling_method",
"user-scaling");
760 std::cout <<
"Call Ipopt Initialize" << std::endl;
762 std::cout <<
"Finished Ipopt Initialize" << std::endl;
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;
772 throw ErrorClass(
"Ipopt FAILED TO SOLVE THE PROBLEM: " +
ipoptErrorMsg);
775 catch(
const ErrorClass& eclass){
777 cout <<
"error in OSIpoptSolver, line 775:\n" << eclass.
errormsg << endl;
780 osresult->setGeneralStatusType(
"error");
782 throw ErrorClass(
osrl) ;
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;
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;
809 if(
osinstance->getVariableNumber() > 0 ||
osinstance->instanceData->objectives->obj != NULL ||
osinstance->instanceData->objectives->numberOfObjectives > 0){
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;
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;
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;
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;
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;
std::string os_dtoa_format(double x)
double os_strtod(const char *s00, char **se)
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
double value
initial value
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.
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...
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
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..
Take an OSResult object and write a string that validates against OSrL.
std::string writeOSrL(OSResult *theosresult)
create an osrl string from an OSResult object
int * hessRowIdx
hessRowIdx is an integer array of row indices in the range 0, ..., n - 1.
int hessDimension
hessDimension is the number of nonzeros in each array.
double * hessValues
hessValues is a double array of the Hessian values.
int * hessColIdx
hessColIdx is an integer array of column indices in the range 0, ..., n - 1.
int * indexes
indexes holds an integer array of variable indices.
int valueSize
valueSize is the dimension of the values array
int * starts
starts holds an integer array of start elements, each start element points to the start of partials f...
double * values
values holds a double array of nonzero partial derivatives