00001
00020
00021
00022 #include <iostream>
00023
00024
00025 #include "OSBonminSolver.h"
00026 #include "OSDataStructures.h"
00027 #include "OSParameters.h"
00028 #include "OSMathUtil.h"
00029 #include "CoinFinite.hpp"
00030 #include "BonOsiTMINLPInterface.hpp"
00031 #include "BonTMINLP.hpp"
00032
00033
00034 #include "CoinTime.hpp"
00035
00036 using std::cout;
00037 using std::endl;
00038 using std::ostringstream;
00039
00040 BonminSolver::BonminSolver() {
00041 osrlwriter = new OSrLWriter();
00042 osresult = new OSResult();
00043 m_osilreader = NULL;
00044 m_osolreader = NULL;
00045 bonminErrorMsg = "";
00046
00047 }
00048
00049 BonminSolver::~BonminSolver() {
00050 #ifdef DEBUG
00051 cout << "inside BonminSolver destructor" << endl;
00052 #endif
00053 if(m_osilreader != NULL) delete m_osilreader;
00054 m_osilreader = NULL;
00055 if(m_osolreader != NULL) delete m_osolreader;
00056 m_osolreader = NULL;
00057
00058 if(osrlwriter != NULL )delete osrlwriter;
00059 osrlwriter = NULL;
00060 if(osresult != NULL ){
00061 delete osresult;
00062 osresult = NULL;
00063
00064 }
00065 #ifdef DEBUG
00066 cout << "leaving BonminSolver destructor" << endl;
00067 #endif
00068 }
00069
00070
00071
00072
00073
00074
00075 bool BonminProblem::get_variables_types(Index n, VariableType* var_types){
00076 int i = 0;
00077 char *varType;
00078 varType = osinstance->getVariableTypes();
00079 n = osinstance->getVariableNumber();
00080 for(i = 0; i < n; i++){
00081 if( varType[i] == 'B') {
00082 var_types[i] = BINARY;
00083 }
00084 else{
00085 if( varType[i] == 'I') {
00086 var_types[i] = INTEGER;
00087 }
00088 else{
00089 var_types[i] = CONTINUOUS;
00090 }
00091 }
00092 }
00093 return true;
00094 }
00095
00096 bool BonminProblem::get_variables_linearity(Index n, Ipopt::TNLP::LinearityType* var_types){
00097
00098
00099
00100 std::cout << "Initialize Nonlinear Structures" << std::endl;
00101 try{
00102 osinstance->initForAlgDiff( );
00103 }
00104 catch(const ErrorClass& eclass){
00105 bonminErrorMsg = eclass.errormsg;
00106 throw;
00107 }
00113 std::map<int, int> varIndexMap;
00114 std::map<int, int>::iterator posVarIndexMap;
00115 varIndexMap = osinstance->getAllNonlinearVariablesIndexMap( );
00116
00117
00118 int i;
00119
00120 for(i = 0; i < n; i++){
00121 var_types[ i] = Ipopt::TNLP::LINEAR;
00122 }
00127 for(posVarIndexMap = varIndexMap.begin(); posVarIndexMap != varIndexMap.end(); ++posVarIndexMap){
00128 std::cout << "Variable Index = " << posVarIndexMap->first << std::endl ;
00129 var_types[ posVarIndexMap->first] = Ipopt::TNLP::NON_LINEAR;
00130 }
00131 std::cout << "Number of nonlinear variables = " << varIndexMap.size() << std::endl;
00132
00133 return true;
00134 }
00135
00136 bool BonminProblem::get_constraints_linearity(Index m, Ipopt::TNLP::LinearityType* const_types){
00137
00138 int i;
00139
00140 for(i = 0; i < m; i++){
00141 const_types[ i] = Ipopt::TNLP::LINEAR;
00142 }
00143
00144
00145 int mm = osinstance->getNumberOfNonlinearExpressionTreeModIndexes();
00146
00147
00148
00149 for(i = 0; i < mm; i++){
00150 if(osinstance->getNonlinearExpressionTreeModIndexes()[ i] >= 0){
00151 std::cout << osinstance->getNonlinearExpressionTreeModIndexes()[ i] << std::endl;
00152 const_types[ osinstance->getNonlinearExpressionTreeModIndexes()[ i] ] = Ipopt::TNLP::NON_LINEAR;
00153
00154 }
00155
00156 }
00157
00158 return true;
00159 }
00160
00161
00162 bool BonminProblem::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00163 Index& nnz_h_lag, TNLP::IndexStyleEnum& index_style)
00164 {
00165 if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Bonmin NEEDS AN OBJECTIVE FUNCTION");
00166
00167 n = osinstance->getVariableNumber();
00168
00169 m = osinstance->getConstraintNumber();
00170 #ifdef DEBUG
00171 cout << "Bonmin number variables !!!!!!!!!!!!!!!!!!!!!!!!!!!" << n << endl;
00172 cout << "Bonmin number constraints !!!!!!!!!!!!!!!!!!!!!!!!!!!" << m << endl;
00173 #endif
00174 try{
00175 osinstance->initForAlgDiff( );
00176 }
00177 catch(const ErrorClass& eclass){
00178 bonminErrorMsg = eclass.errormsg;
00179 throw;
00180 }
00181
00182 osinstance->bUseExpTreeForFunEval = true;
00183
00184 SparseJacobianMatrix *sparseJacobian = NULL;
00185 try{
00186 sparseJacobian = osinstance->getJacobianSparsityPattern();
00187
00188 }
00189 catch(const ErrorClass& eclass){
00190 bonminErrorMsg = eclass.errormsg;
00191 throw;
00192 }
00193
00194 nnz_jac_g = sparseJacobian->valueSize;
00195 #ifdef DEBUG
00196 cout << "nnz_jac_g !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_jac_g << endl;
00197 #endif
00198
00199
00200 if( (osinstance->getNumberOfNonlinearExpressions() == 0) && (osinstance->getNumberOfQuadraticTerms() == 0) ) {
00201 cout << "This is a linear program" << endl;
00202 nnz_h_lag = 0;
00203 }
00204 else{
00205
00206 SparseHessianMatrix *sparseHessian = osinstance->getLagrangianHessianSparsityPattern();
00207
00208 nnz_h_lag = sparseHessian->hessDimension;
00209 }
00210 #ifdef DEBUG
00211 cout << "nnz_h_lag !!!!!!!!!!!!!!!!!!!!!!!!!!!" << nnz_h_lag << endl;
00212 #endif
00213
00214 index_style = TNLP::C_STYLE;
00215
00216 return true;
00217 }
00218
00219
00220 bool BonminProblem::get_bounds_info(Index n, Number* x_l, Number* x_u,
00221 Index m, Number* g_l, Number* g_u){
00222 int i;
00223
00224 double * mdVarLB = osinstance->getVariableLowerBounds();
00225
00226
00227 double * mdVarUB = osinstance->getVariableUpperBounds();
00228
00229 for(i = 0; i < n; i++){
00230 x_l[ i] = mdVarLB[ i];
00231 x_u[ i] = mdVarUB[ i];
00232 #ifdef DEBUG
00233 cout << "x_l !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_l[i] << endl;
00234 cout << "x_u !!!!!!!!!!!!!!!!!!!!!!!!!!!" << x_u[i] << endl;
00235 #endif
00236 }
00237
00238
00239
00240
00241
00242
00243 double * mdConLB = osinstance->getConstraintLowerBounds();
00244
00245 double * mdConUB = osinstance->getConstraintUpperBounds();
00246
00247 for(int i = 0; i < m; i++){
00248 g_l[ i] = mdConLB[ i];
00249 g_u[ i] = mdConUB[ i];
00250 #ifdef DEBUG
00251 cout << "lower !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_l[i] << endl;
00252 cout << "upper !!!!!!!!!!!!!!!!!!!!!!!!!!!" << g_u[i] << endl;
00253 #endif
00254 }
00255 return true;
00256 }
00257
00258
00259
00260 bool BonminProblem::get_starting_point(Index n, bool init_x, Number* x,
00261 bool init_z, Number* z_L, Number* z_U, Index m, bool init_lambda,
00262 Number* lambda) {
00263
00264
00265
00266 assert(init_x == true);
00267 assert(init_z == false);
00268 assert(init_lambda == false);
00269 int i, m1, n1;
00270
00271
00272 #ifdef DEBUG
00273 cout << "get initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00274 #endif
00275
00276
00277 #ifdef DEBUG
00278 cout << "get number of initial values !!!!!!!!!!!!!!!!!!!!!!!!!! " << endl;
00279 #endif
00280 int k;
00281 if (osoption != NULL)
00282 m1 = osoption->getNumberOfInitVarValues();
00283 else
00284 m1 = 0;
00285 #ifdef DEBUG
00286 cout << "number of variables initialed: " << m1 << endl;
00287 #endif
00288
00289 n1 = osinstance->getVariableNumber();
00290 bool* initialed;
00291 initialed = new bool[n1];
00292 #ifdef DEBUG
00293 cout << "number of variables in total: " << n1 << endl;
00294 #endif
00295
00296
00297 for(k = 0; k < n1; k++)
00298 initialed[k] = false;
00299
00300 if (m1 > 0)
00301 {
00302 #ifdef DEBUG
00303 cout << "get initial values " << endl;
00304 #endif
00305
00306 InitVarValue** initVarVector = osoption->getInitVarValuesSparse();
00307 #ifdef DEBUG
00308 cout << "done " << endl;
00309 #endif
00310
00311 double initval;
00312 try
00313 {
00314 for(k = 0; k < m1; k++)
00315 {
00316 i = initVarVector[k]->idx;
00317 if (initVarVector[k]->idx > n1)
00318 throw ErrorClass ("Illegal index value in variable initialization");
00319
00320 initval = initVarVector[k]->value;
00321 if (osinstance->instanceData->variables->var[i]->ub == OSDBL_MAX)
00322 { if (osinstance->instanceData->variables->var[i]->lb > initval)
00323 throw ErrorClass ("Initial value outside of bounds");
00324 }
00325 else
00326 if (osinstance->instanceData->variables->var[i]->lb == -OSDBL_MAX)
00327 { if (osinstance->instanceData->variables->var[i]->ub < initval)
00328 throw ErrorClass ("Initial value outside of bounds");
00329 }
00330 else
00331 { if ((osinstance->instanceData->variables->var[i]->lb > initval) ||
00332 (osinstance->instanceData->variables->var[i]->ub < initval))
00333 throw ErrorClass ("Initial value outside of bounds");
00334 }
00335
00336 x[initVarVector[k]->idx] = initval;
00337 initialed[initVarVector[k]->idx] = true;
00338 }
00339 }
00340 catch(const ErrorClass& eclass){
00341 cout << "Error in BonminProblem::get_starting_point (OSBonminSolver.cpp)";
00342 cout << endl << endl << endl;
00343 }
00344 }
00345
00346 double default_initval;
00347 default_initval = 1.7171;
00348
00349
00350 for(k = 0; k < n1; k++)
00351 {
00352 if (!initialed[k])
00353 if (osinstance->instanceData->variables->var[k]->ub == OSDBL_MAX)
00354 if (osinstance->instanceData->variables->var[k]->lb <= default_initval)
00355 x[k] = default_initval;
00356 else
00357 x[k] = osinstance->instanceData->variables->var[k]->lb;
00358 else
00359 if (osinstance->instanceData->variables->var[k]->lb == -OSDBL_MAX)
00360 if (osinstance->instanceData->variables->var[k]->ub >= default_initval)
00361 x[k] = default_initval;
00362 else
00363 x[k] = osinstance->instanceData->variables->var[k]->ub;
00364 else
00365 if ((osinstance->instanceData->variables->var[k]->lb <= default_initval) &&
00366 (osinstance->instanceData->variables->var[k]->ub >= default_initval))
00367 x[k] = default_initval;
00368 else
00369 if (osinstance->instanceData->variables->var[k]->lb > default_initval)
00370 x[k] = osinstance->instanceData->variables->var[k]->lb;
00371 else
00372 x[k] = osinstance->instanceData->variables->var[k]->ub;
00373 }
00374 for(i = 0; i < n1; i++){
00375 #ifdef DEBUG
00376 std::cout << "INITIAL VALUE !!!!!!!!!!!!!!!!!!!! " << x[ i] << std::endl;
00377 #endif
00378 }
00379
00380 osinstance->calculateAllObjectiveFunctionValues( x, true);
00381 delete[] initialed;
00382 return true;
00383 }
00384
00385
00386
00387 bool BonminProblem::eval_f(Index n, const Number* x, bool new_x, Number& obj_value){
00388
00389 try{
00390 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") == 0){
00391 obj_value = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 )[ 0];
00392 }else{
00393 obj_value = -osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 )[ 0];
00394 }
00395
00396 }
00397 catch(const ErrorClass& eclass){
00398 bonminErrorMsg = eclass.errormsg;
00399 throw;
00400 }
00401 if( CoinIsnan( (double)obj_value) ) return false;
00402 return true;
00403 }
00404
00405 bool BonminProblem::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f){
00406 int i;
00407 double *objGrad;
00408 try{
00409
00410
00411
00412 objGrad = osinstance->calculateObjectiveFunctionGradient( const_cast<double*>(x), NULL, NULL, -1, new_x, 1);
00413 }
00414 catch(const ErrorClass& eclass){
00415 bonminErrorMsg = eclass.errormsg;
00416 throw;
00417 }
00418 for(i = 0; i < n; i++){
00419 if( osinstance->instanceData->objectives->obj[ 0]->maxOrMin.compare("min") == 0){
00420 grad_f[ i] = objGrad[ i];
00421 }else{
00422 grad_f[ i] = -objGrad[ i];
00423 }
00424
00425 }
00426
00427 return true;
00428 }
00429
00430
00431 bool BonminProblem::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g) {
00432 try{
00433 double *conVals = osinstance->calculateAllConstraintFunctionValues( const_cast<double*>(x), NULL, NULL, new_x, 0 );
00434 int i;
00435 for(i = 0; i < m; i++){
00436 if( CoinIsnan( (double)conVals[ i] ) ) return false;
00437 g[i] = conVals[ i] ;
00438 }
00439 return true;
00440 }
00441 catch(const ErrorClass& eclass){
00442 bonminErrorMsg = eclass.errormsg;
00443 throw;
00444 }
00445 }
00446
00447
00448
00449 bool BonminProblem::eval_jac_g(Index n, const Number* x, bool new_x,
00450 Index m, Index nele_jac, Index* iRow, Index *jCol,
00451 Number* values){
00452 SparseJacobianMatrix *sparseJacobian;
00453 if (values == NULL) {
00454
00455
00456
00457
00458
00459 try{
00460 sparseJacobian = osinstance->getJacobianSparsityPattern();
00461 }
00462 catch(const ErrorClass& eclass){
00463 bonminErrorMsg = eclass.errormsg;
00464 throw;
00465 }
00466 int i = 0;
00467 int k, idx;
00468 for(idx = 0; idx < m; idx++){
00469 for(k = *(sparseJacobian->starts + idx); k < *(sparseJacobian->starts + idx + 1); k++){
00470 iRow[i] = idx;
00471 jCol[i] = *(sparseJacobian->indexes + k);
00472
00473
00474 i++;
00475 }
00476 }
00477 }
00478 else {
00479
00480 try{
00481 sparseJacobian = osinstance->calculateAllConstraintFunctionGradients( const_cast<double*>(x), NULL, NULL, new_x, 1);
00482 }
00483 catch(const ErrorClass& eclass){
00484 bonminErrorMsg = eclass.errormsg;
00485 throw;
00486 }
00487
00488 for(int i = 0; i < nele_jac; i++){
00489 values[ i] = sparseJacobian->values[i];
00490
00491
00492
00493 }
00494 }
00495 return true;
00496 }
00497
00498
00499 bool BonminProblem::eval_h(Index n, const Number* x, bool new_x,
00500 Number obj_factor, Index m, const Number* lambda,
00501 bool new_lambda, Index nele_hess, Index* iRow,
00502 Index* jCol, Number* values){
00503
00505 SparseHessianMatrix *sparseHessian;
00506
00507 int i;
00508 if (values == NULL) {
00509
00510
00511 try{
00512 sparseHessian = osinstance->getLagrangianHessianSparsityPattern( );
00513 }
00514 catch(const ErrorClass& eclass){
00515 bonminErrorMsg = eclass.errormsg;
00516 throw;
00517 }
00518
00519 for(i = 0; i < nele_hess; i++){
00520 iRow[i] = *(sparseHessian->hessColIdx + i);
00521 jCol[i] = *(sparseHessian->hessRowIdx + i);
00522
00523
00524 }
00525 }
00526 else {
00527
00528
00529 double* objMultipliers = new double[1];
00530 objMultipliers[0] = obj_factor;
00531 try{
00532 sparseHessian = osinstance->calculateLagrangianHessian( const_cast<double*>(x), objMultipliers, const_cast<double*>(lambda) , new_x, 2);
00533 delete[] objMultipliers;
00534 }
00535 catch(const ErrorClass& eclass){
00536 bonminErrorMsg = eclass.errormsg;
00537 delete[] objMultipliers;
00538 throw;
00539 }
00540 for(i = 0; i < nele_hess; i++){
00541 values[ i] = *(sparseHessian->hessValues + i);
00542 }
00543 }
00545 return true;
00546 }
00547
00548 bool BonminProblem::get_scaling_parameters(Number& obj_scaling,
00549 bool& use_x_scaling, Index n,
00550 Number* x_scaling,
00551 bool& use_g_scaling, Index m,
00552 Number* g_scaling){
00553
00554
00555
00556
00557
00558
00559
00560
00561
00562 return true;
00563 }
00564
00565
00566
00567
00568 void BonminProblem::finalize_solution(TMINLP::SolverReturn status_,
00569 Index n, const Number* x, Number obj_value)
00570 {
00571
00572 status = status_;
00573 #ifdef DEBUG
00574 std::cout << "FINALIZE OBJ SOLUTION VALUE = " << obj_value << std::endl;
00575 #endif
00576
00577 }
00578
00579
00580
00581
00582 void BonminSolver::buildSolverInstance() throw (ErrorClass) {
00583 try{
00584
00585 if(osil.length() == 0 && osinstance == NULL) throw ErrorClass("there is no instance");
00586 if(osinstance == NULL){
00587 m_osilreader = new OSiLReader();
00588 osinstance = m_osilreader->readOSiL( osil);
00589 }
00590
00591 tminlp = new BonminProblem( osinstance, osoption);
00592 this->bCallbuildSolverInstance = true;
00593
00594
00595 }
00596 catch(const ErrorClass& eclass){
00597 std::cout << "THERE IS AN ERROR" << std::endl;
00598 osresult->setGeneralMessage( eclass.errormsg);
00599 osresult->setGeneralStatusType( "error");
00600 osrl = osrlwriter->writeOSrL( osresult);
00601 throw ErrorClass( osrl) ;
00602 }
00603 }
00604
00605
00606
00607 void BonminSolver::setSolverOptions() throw (ErrorClass) {
00608 try{
00609 this->bSetSolverOptions = true;
00610 bonminSetup.initializeOptionsAndJournalist();
00611
00612 bonminSetup.roptions()->AddStringOption2("print_solution","Do we print the solution or not?",
00613 "yes",
00614 "no", "No, we don't.",
00615 "yes", "Yes, we do.",
00616 "A longer comment can be put here");
00617
00618
00619 bonminSetup.options()->SetNumericValue("bonmin.time_limit", 5000);
00620
00621
00622
00623
00624
00625
00626
00627
00628
00629
00630
00631
00632
00633
00634
00635
00636
00637 bonminSetup.readOptionsString("bonmin.algorithm B-BB\n");
00638
00639
00640 bonminSetup.options()->SetIntegerValue("bonmin.bb_log_level", 0 );
00641 bonminSetup.options()->SetIntegerValue("bonmin.nlp_log_level", 0 );
00642
00643
00644 int printSolution;
00645 bonminSetup.options()->GetEnumValue("print_solution", printSolution,"");
00646 if(printSolution == 1){
00647 tminlp->printSolutionAtEndOfAlgorithm();
00648 }
00649
00650 if(osoption == NULL && osol.length() > 0)
00651 {
00652 m_osolreader = new OSoLReader();
00653 osoption = m_osolreader->readOSoL( osol);
00654 }
00655
00656 if(osoption != NULL && osoption->getNumberOfSolverOptions() > 0 ){
00657 char *pEnd;
00658 int i;
00659 std::vector<SolverOption*> optionsVector;
00660 optionsVector = osoption->getSolverOptions( "bonmin");
00661 int num_bonmin_options = optionsVector.size();
00662 for(i = 0; i < num_bonmin_options; i++){
00663 if(optionsVector[ i]->type == "numeric" ){
00664 std::cout << "FOUND A NUMERIC OPTION " << os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) << std::endl;
00665 if(optionsVector[ i]->category == "ipopt"){
00666 bonminSetup.options()->SetNumericValue(optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
00667 }else{
00668 if(optionsVector[ i]->category == "cbc" ){
00669 bonminSetup.options()->SetNumericValue("milp_solver."+optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
00670 }
00671 else{
00672 bonminSetup.options()->SetNumericValue("bonmin."+optionsVector[ i]->name, os_strtod( optionsVector[ i]->value.c_str(), &pEnd ) );
00673 }
00674 }
00675 }
00676 else if(optionsVector[ i]->type == "integer" ){
00677 std::cout << "FOUND AN INTEGER OPTION " <<optionsVector[ i]->name << std::endl;
00678 if(optionsVector[ i]->category == "ipopt"){
00679 bonminSetup.options()->SetIntegerValue(optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ) );
00680 }else{
00681 if(optionsVector[ i]->category == "cbc" ){
00682 std::cout << "SETTING INTEGER CBC OPTION" << std::endl;
00683 bonminSetup.options()->SetIntegerValue("milp_solver."+optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ));
00684 }
00685 else{
00686 bonminSetup.options()->SetIntegerValue("bonmin."+optionsVector[ i]->name, atoi( optionsVector[ i]->value.c_str() ) );
00687 }
00688 }
00689 }
00690 else if(optionsVector[ i]->type == "string" ){
00691 std::cout << "FOUND A STRING OPTION " <<optionsVector[ i]->name << std::endl;
00692 if(optionsVector[ i]->category == "ipopt"){
00693 bonminSetup.options()->SetStringValue(optionsVector[ i]->name, optionsVector[ i]->value );
00694 }else{
00695 if(optionsVector[ i]->category == "cbc" ){
00696 bonminSetup.options()->SetStringValue("milp_solver."+optionsVector[ i]->name, optionsVector[ i]->value);
00697 }
00698 else{
00699 bonminSetup.options()->SetStringValue("bonmin."+optionsVector[ i]->name, optionsVector[ i]->value);
00700 }
00701 }
00702
00703 }
00704 }
00705 }
00706
00707 }
00708
00709 catch(const ErrorClass& eclass){
00710 std::cout << "THERE IS AN ERROR" << std::endl;
00711 osresult->setGeneralMessage( eclass.errormsg);
00712 osresult->setGeneralStatusType( "error");
00713 osrl = osrlwriter->writeOSrL( osresult);
00714 throw ErrorClass( osrl) ;
00715 }
00716 }
00717
00718
00719
00720 void BonminSolver::solve() throw (ErrorClass) {
00721 if( this->bCallbuildSolverInstance == false) buildSolverInstance();
00722 if( this->bSetSolverOptions == false) setSolverOptions();
00723 try{
00724
00725 if(osinstance->getObjectiveNumber() <= 0) throw ErrorClass("Bonmin NEEDS AN OBJECTIVE FUNCTION");
00726
00727
00728
00729 if(osinstance->getVariableNumber() <= 0)throw ErrorClass("Bonmin requires decision variables");
00730
00731
00732 try {
00733 bonminSetup.initialize( GetRawPtr(tminlp) );
00734
00735 bb( bonminSetup);
00736
00737 }
00738 catch(TNLPSolver::UnsolvedError *E) {
00739
00740 std::cerr<<"Ipopt has failed to solve a problem"<<std::endl;
00741 }
00742 catch(OsiTMINLPInterface::SimpleError &E) {
00743 std::cerr<<E.className()<<"::"<<E.methodName()
00744 <<std::endl
00745 <<E.message()<<std::endl;
00746 }
00747 catch(CoinError &E) {
00748 std::cerr<<E.className()<<"::"<<E.methodName()
00749 <<std::endl
00750 <<E.message()<<std::endl;
00751 }
00752
00753 if(( bb.model().isContinuousUnbounded() == true) && (osinstance->getNumberOfIntegerVariables() + osinstance->getNumberOfBinaryVariables() <= 0) ){
00754 std::string solutionDescription = "";
00755 std::string message = "Success";
00756 int solIdx = 0;
00757 if(osresult->setServiceName( "Bonin solver service") != true)
00758 throw ErrorClass("OSResult error: setServiceName");
00759 if(osresult->setInstanceName( osinstance->getInstanceName()) != true)
00760 throw ErrorClass("OSResult error: setInstanceName");
00761 if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00762 throw ErrorClass("OSResult error: setVariableNumer");
00763 if(osresult->setObjectiveNumber( 1) != true)
00764 throw ErrorClass("OSResult error: setObjectiveNumber");
00765 if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00766 throw ErrorClass("OSResult error: setConstraintNumber");
00767 if(osresult->setSolutionNumber( 1) != true)
00768 throw ErrorClass("OSResult error: setSolutionNumer");
00769 if(osresult->setGeneralMessage( message) != true)
00770 throw ErrorClass("OSResult error: setGeneralMessage");
00771 solutionDescription = "The problem is unbounded";
00772 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00773 osresult->setGeneralStatusType("normal");
00774 osrl = osrlwriter->writeOSrL( osresult);
00775 return;
00776 }
00777
00778
00779 if(( bb.model().isProvenInfeasible() == true) ){
00780 std::string solutionDescription = "";
00781 std::string message = "Success";
00782 int solIdx = 0;
00783 if(osresult->setServiceName( "Bonin solver service") != true)
00784 throw ErrorClass("OSResult error: setServiceName");
00785 if(osresult->setInstanceName( osinstance->getInstanceName()) != true)
00786 throw ErrorClass("OSResult error: setInstanceName");
00787 if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00788 throw ErrorClass("OSResult error: setVariableNumer");
00789 if(osresult->setObjectiveNumber( 1) != true)
00790 throw ErrorClass("OSResult error: setObjectiveNumber");
00791 if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00792 throw ErrorClass("OSResult error: setConstraintNumber");
00793 if(osresult->setSolutionNumber( 1) != true)
00794 throw ErrorClass("OSResult error: setSolutionNumer");
00795 if(osresult->setGeneralMessage( message) != true)
00796 throw ErrorClass("OSResult error: setGeneralMessage");
00797 solutionDescription = "The problem is infeasible";
00798 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00799 osresult->setGeneralStatusType("normal");
00800 osrl = osrlwriter->writeOSrL( osresult);
00801 return;
00802 }
00803 status = tminlp->status;
00804 writeResult();
00805
00806 }
00807 catch(const ErrorClass& eclass){
00808 osresult->setGeneralMessage( eclass.errormsg);
00809 osresult->setGeneralStatusType( "error");
00810 osrl = osrlwriter->writeOSrL( osresult);
00811 throw ErrorClass( osrl) ;
00812 }
00813 }
00814
00815
00816
00817
00818 void BonminSolver::writeResult(){
00819 double *x = NULL;
00820 double *z = NULL;
00821 int i = 0;
00822 int solIdx = 0;
00823 std::string solutionDescription = "";
00824 std::string message = "Bonmin solver finishes to the end.";
00825
00826
00827 try{
00828 x = new double[osinstance->getVariableNumber() ];
00829 z = new double[1];
00830
00831 if(osresult->setServiceName( "Bonmin solver service") != true)
00832 throw ErrorClass("OSResult error: setServiceName");
00833 if(osresult->setInstanceName( osinstance->getInstanceName()) != true)
00834 throw ErrorClass("OSResult error: setInstanceName");
00835
00836
00837
00838
00839 if(osresult->setVariableNumber( osinstance->getVariableNumber()) != true)
00840 throw ErrorClass("OSResult error: setVariableNumer");
00841 if(osresult->setObjectiveNumber( 1) != true)
00842 throw ErrorClass("OSResult error: setObjectiveNumber");
00843 if(osresult->setConstraintNumber( osinstance->getConstraintNumber()) != true)
00844 throw ErrorClass("OSResult error: setConstraintNumber");
00845 if(osresult->setSolutionNumber( 1) != true)
00846 throw ErrorClass("OSResult error: setSolutionNumer");
00847 if(osresult->setGeneralMessage( message) != true)
00848 throw ErrorClass("OSResult error: setGeneralMessage");
00849
00850 switch( status){
00851 case TMINLP::SUCCESS:
00852 solutionDescription = "SUCCESS[BONMIN]: Algorithm terminated normally at a locally optimal point, satisfying the convergence tolerances.";
00853
00854 osresult->setSolutionStatus(solIdx, "locallyOptimal", solutionDescription);
00855
00856 *(z + 0) = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(bb.bestSolution()), true)[ 0];
00857 osresult->setObjectiveValuesDense(solIdx, z);
00858 for(i=0; i < osinstance->getVariableNumber(); i++){
00859 *(x + i) = bb.bestSolution()[i];
00860
00861 }
00862 osresult->setPrimalVariableValuesDense(solIdx, x);
00863 break;
00864
00865 case TMINLP::LIMIT_EXCEEDED:
00866 solutionDescription = "LIMIT_EXCEEDED[BONMIN]: A resource limit was exceeded, we provide the current solution.";
00867
00868 osresult->setSolutionStatus(solIdx, "stoppedByLimit", solutionDescription);
00869
00870
00871
00872 *(z + 0) = osinstance->calculateAllObjectiveFunctionValues( const_cast<double*>(bb.model().getColSolution()), true)[ 0];
00873 osresult->setObjectiveValuesDense(solIdx, z);
00874 for(i=0; i < osinstance->getVariableNumber(); i++){
00875 *(x + i) = bb.model().getColSolution()[i];
00876
00877 }
00878 osresult->setPrimalVariableValuesDense(solIdx, x);
00879 break;
00880
00881 case TMINLP::MINLP_ERROR:
00882 solutionDescription = "MINLP_ERROR [BONMIN]: Algorithm stopped with unspecified error.";
00883
00884 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00885
00886 break;
00887
00888 case TMINLP::CONTINUOUS_UNBOUNDED:
00889 solutionDescription = "CONTINUOUS_UNBOUNDED [BONMIN]: The continuous relaxation is unbounded, the MINLP may or may not be unbounded.";
00890
00891 osresult->setSolutionStatus(solIdx, "error", solutionDescription);
00892
00893 break;
00894
00895
00896 case TMINLP::INFEASIBLE:
00897 solutionDescription = "INFEASIBLE [BONMIN]: Problem may be infeasible.";
00898
00899 osresult->setSolutionStatus(solIdx, "infeasible", solutionDescription);
00900 break;
00901
00902
00903 default:
00904 solutionDescription = "OTHER[BONMIN]: other unknown solution status from Bonmin solver";
00905
00906 osresult->setSolutionStatus(solIdx, "other", solutionDescription);
00907 }
00908 osresult->setGeneralStatusType("normal");
00909 osrl = osrlwriter->writeOSrL( osresult);
00910 delete[] x;
00911 x = NULL;
00912 delete[] z;
00913 z = NULL;
00914 }
00915
00916
00917 catch(const ErrorClass& eclass){
00918 delete[] x;
00919 x = NULL;
00920 delete[] z;
00921 z = NULL;
00922 osresult->setGeneralMessage( eclass.errormsg);
00923 osresult->setGeneralStatusType( "error");
00924 osrl = osrlwriter->writeOSrL( osresult);
00925 throw ErrorClass( osrl) ;
00926 }
00927
00928
00929 }
00930
00931
00932 void BonminSolver::dataEchoCheck(){
00933
00934 int i;
00935
00936
00937 cout << "This is problem: " << osinstance->getInstanceName() << endl;
00938 cout << "The problem source is: " << osinstance->getInstanceSource() << endl;
00939 cout << "The problem description is: " << osinstance->getInstanceDescription() << endl;
00940 cout << "number of variables = " << osinstance->getVariableNumber() << endl;
00941 cout << "number of Rows = " << osinstance->getConstraintNumber() << endl;
00942
00943
00944 if(osinstance->getVariableNumber() > 0){
00945 for(i = 0; i < osinstance->getVariableNumber(); i++){
00946 if(osinstance->getVariableNames() != NULL) cout << "variable Names " << osinstance->getVariableNames()[ i] << endl;
00947 if(osinstance->getVariableTypes() != NULL) cout << "variable Types " << osinstance->getVariableTypes()[ i] << endl;
00948 if(osinstance->getVariableLowerBounds() != NULL) cout << "variable Lower Bounds " << osinstance->getVariableLowerBounds()[ i] << endl;
00949 if(osinstance->getVariableUpperBounds() != NULL) cout << "variable Upper Bounds " << osinstance->getVariableUpperBounds()[i] << endl;
00950 }
00951 }
00952
00953
00954 if(osinstance->getVariableNumber() > 0 || osinstance->instanceData->objectives->obj != NULL || osinstance->instanceData->objectives->numberOfObjectives > 0){
00955 if( osinstance->getObjectiveMaxOrMins()[0] == "min") cout << "problem is a minimization" << endl;
00956 else cout << "problem is a maximization" << endl;
00957 for(i = 0; i < osinstance->getVariableNumber(); i++){
00958 cout << "OBJ COEFFICIENT = " << osinstance->getDenseObjectiveCoefficients()[0][i] << endl;
00959 }
00960 }
00961
00962 if(osinstance->getConstraintNumber() > 0){
00963 for(i = 0; i < osinstance->getConstraintNumber(); i++){
00964 if(osinstance->getConstraintNames() != NULL) cout << "row name = " << osinstance->getConstraintNames()[i] << endl;
00965 if(osinstance->getConstraintLowerBounds() != NULL) cout << "row lower bound = " << osinstance->getConstraintLowerBounds()[i] << endl;
00966 if(osinstance->getConstraintUpperBounds() != NULL) cout << "row upper bound = " << osinstance->getConstraintUpperBounds()[i] << endl;
00967 }
00968 }
00969
00970
00971 cout << endl;
00972 cout << "number of nonzeros = " << osinstance->getLinearConstraintCoefficientNumber() << endl;
00973 for(i = 0; i <= osinstance->getVariableNumber(); i++){
00974 cout << "Start Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->starts[ i] << endl;
00975 }
00976 cout << endl;
00977 for(i = 0; i < osinstance->getLinearConstraintCoefficientNumber(); i++){
00978 cout << "Index Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->indexes[i] << endl;
00979 cout << "Nonzero Value = " << osinstance->getLinearConstraintCoefficientsInColumnMajor()->values[i] << endl;
00980 }
00981
00982
00983 cout << "number of qterms = " << osinstance->getNumberOfQuadraticTerms() << endl;
00984 for(int i = 0; i < osinstance->getNumberOfQuadraticTerms(); i++){
00985 cout << "Row Index = " << osinstance->getQuadraticTerms()->rowIndexes[i] << endl;
00986 cout << "Var Index 1 = " << osinstance->getQuadraticTerms()->varOneIndexes[ i] << endl;
00987 cout << "Var Index 2 = " << osinstance->getQuadraticTerms()->varTwoIndexes[ i] << endl;
00988 cout << "Coefficient = " << osinstance->getQuadraticTerms()->coefficients[ i] << endl;
00989 }
00990 }
00991
00992
00993 BonminProblem::BonminProblem(OSInstance *osinstance_, OSOption *osoption_) {
00994 osinstance = osinstance_;
00995 osoption = osoption_;
00996 printSol_ = false;
00997 }
00998
00999 BonminProblem::~BonminProblem() {
01000
01001 }
01002
01003
01004