/home/coin/SVN-release/OS-2.4.2/Bonmin/src/Interfaces/Ampl/BonAmplTMINLP.cpp

Go to the documentation of this file.
00001 // (C) Copyright International Business Machines Corporation and
00002 // Carnegie Mellon University 2004-2011
00003 //
00004 // All Rights Reserved.
00005 // This code is published under the Common Public License.
00006 //
00007 // Authors :
00008 // Carl D. Laird, Carnegie Mellon University,
00009 // Andreas Waechter, International Business Machines Corporation
00010 // Pierre Bonami, Carnegie Mellon University,
00011 //
00012 // Date : 12/01/2004
00013 #include "IpBlas.hpp"
00014 
00015 #include <list>
00016 
00017 #include "AmplTNLP.hpp"
00018 #include "BonAmplTMINLP.hpp"
00019 #include <iostream>
00020 #include <fstream>
00021 
00022 #include "CoinHelperFunctions.hpp"
00023 #include "BonExitCodes.hpp"
00024 
00025 using namespace Ipopt;
00026 
00027 namespace Bonmin{
00028   void 
00029   RegisteredOptions::fillAmplOptionList(ExtraCategoriesInfo which, Ipopt::AmplOptionsList * amplOptList){
00030       std::list<int> test;
00031       std::list< Ipopt::RegisteredOption * > options;
00032       chooseOptions(which, options);
00033       for(std::list< Ipopt::RegisteredOption * >::iterator i = options.begin();
00034            i != options.end() ; i++)
00035        {
00036            std::string name = "bonmin.";
00037            name += (*i)->Name();
00038            Ipopt::RegisteredOptionType T = (*i)->Type();
00039            Ipopt::AmplOptionsList::AmplOptionType  type;
00040            switch(T){
00041              case Ipopt::OT_Number: type = Ipopt::AmplOptionsList::Number_Option;
00042                   break;
00043              case Ipopt::OT_Integer: type = Ipopt::AmplOptionsList::Integer_Option;
00044                   break;
00045              case Ipopt::OT_String: type = Ipopt::AmplOptionsList::String_Option;
00046                   break;
00047              case Ipopt::OT_Unknown:
00048              default:
00049                 throw CoinError("RegisteredOptions","fillAmplOptionList","Unknown option type");
00050            }
00051            amplOptList->AddAmplOption(name, name, type, (*i)->ShortDescription());
00052        }
00053 }
00054 }
00055 #include "asl.h"
00056 #include "asl_pfgh.h"
00057 #include "getstub.h"
00058 
00059 namespace ampl_utils
00060 {
00061   void sos_kludge(int nsos, int *sosbeg, double *sosref);
00062 }
00063 namespace Bonmin
00064 {
00065 
00066   AmplTMINLP::AmplTMINLP()
00067       :
00068       TMINLP(),
00069       appName_(),
00070       upperBoundingObj_(-1),
00071       ampl_tnlp_(NULL),
00072       branch_(),
00073       sos_(),
00074       suffix_handler_(NULL),
00075       constraintsConvexities_(NULL),
00076       numberNonConvex_(0),
00077       nonConvexConstraintsAndRelaxations_(NULL),
00078       numberSimpleConcave_(0),
00079       simpleConcaves_(NULL),
00080       hasLinearObjective_(false)
00081   {}
00082 
00083 
00084   AmplTMINLP::AmplTMINLP(const SmartPtr<const Journalist>& jnlst,
00085       const SmartPtr<Bonmin::RegisteredOptions> roptions,
00086       const SmartPtr<OptionsList> options,
00087       char**& argv,
00088       AmplSuffixHandler* suffix_handler /*=NULL*/,
00089       const std::string& appName,
00090       std::string* nl_file_content /* = NULL */
00091                         )
00092       :
00093       TMINLP(),
00094       appName_(),
00095       upperBoundingObj_(-1),
00096       ampl_tnlp_(NULL),
00097       branch_(),
00098       sos_(),
00099       suffix_handler_(NULL),
00100       constraintsConvexities_(NULL),
00101       numberNonConvex_(0),
00102       nonConvexConstraintsAndRelaxations_(NULL),
00103       numberSimpleConcave_(0),
00104       simpleConcaves_(NULL),
00105       hasLinearObjective_(false)
00106   {
00107     Initialize(jnlst, roptions, options, argv, suffix_handler, appName, nl_file_content);
00108   }
00109 
00110   void
00111   AmplTMINLP::Initialize(const SmartPtr<const Journalist>& jnlst,
00112       const SmartPtr<Bonmin::RegisteredOptions> roptions,
00113       const SmartPtr<OptionsList> options,
00114       char**& argv,
00115       AmplSuffixHandler* suffix_handler /*=NULL*/,
00116       const std::string& appName,
00117       std::string* nl_file_content /* = NULL */
00118                         )
00119   {
00120     appName_ = appName;
00121     options->GetEnumValue("file_solution",writeAmplSolFile_,"bonmin.");
00122     jnlst_ = jnlst;
00123 
00124     if (suffix_handler==NULL)
00125       suffix_handler_ = suffix_handler = new AmplSuffixHandler();
00126 
00127     // Add the suffix handler for scaling
00128     suffix_handler->AddAvailableSuffix("scaling_factor", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00129     suffix_handler->AddAvailableSuffix("scaling_factor", AmplSuffixHandler::Constraint_Source, AmplSuffixHandler::Number_Type);
00130     suffix_handler->AddAvailableSuffix("scaling_factor", AmplSuffixHandler::Objective_Source, AmplSuffixHandler::Number_Type);
00131 
00132     // Modified for warm-start from AMPL
00133     suffix_handler->AddAvailableSuffix("ipopt_zL_out", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00134     suffix_handler->AddAvailableSuffix("ipopt_zU_out", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00135     suffix_handler->AddAvailableSuffix("ipopt_zL_in", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00136     suffix_handler->AddAvailableSuffix("ipopt_zU_in", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00137 
00138     // priority suffix
00139     suffix_handler->AddAvailableSuffix("priority", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Index_Type);
00140     suffix_handler->AddAvailableSuffix("direction", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00141     suffix_handler->AddAvailableSuffix("downPseudocost", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00142     suffix_handler->AddAvailableSuffix("upPseudocost", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00143 
00144 
00145 
00146     // sos suffixes
00147     suffix_handler->AddAvailableSuffix("ref", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00148     suffix_handler->AddAvailableSuffix("sos",AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Index_Type);
00149     suffix_handler->AddAvailableSuffix("sos",AmplSuffixHandler::Constraint_Source, AmplSuffixHandler::Index_Type);
00150     suffix_handler->AddAvailableSuffix("sosno",AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00151     suffix_handler->AddAvailableSuffix("sosref",AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00152     suffix_handler->AddAvailableSuffix("sstatus", AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Index_Type);
00153     suffix_handler->AddAvailableSuffix("sstatus", AmplSuffixHandler::Constraint_Source, AmplSuffixHandler::Index_Type);
00154 
00155 
00156     // For marking convex/nonconvex constraints
00157     suffix_handler->AddAvailableSuffix("id",AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Index_Type);
00158     suffix_handler->AddAvailableSuffix("primary_var",AmplSuffixHandler::Constraint_Source, AmplSuffixHandler::Index_Type);
00159 
00160     // For objectives
00161     suffix_handler->AddAvailableSuffix("UBObj", AmplSuffixHandler::Objective_Source, AmplSuffixHandler::Index_Type);
00162 
00163 
00164     // Perturbation radius
00165     suffix_handler->AddAvailableSuffix("perturb_radius",AmplSuffixHandler::Variable_Source, AmplSuffixHandler::Number_Type);
00166 
00167     SmartPtr<AmplOptionsList> ampl_options_list = new AmplOptionsList();
00168     roptions->fillAmplOptionList(RegisteredOptions::BonminCategory, GetRawPtr(ampl_options_list));
00169     roptions->fillAmplOptionList(RegisteredOptions::FilterCategory, GetRawPtr(ampl_options_list));
00170     roptions->fillAmplOptionList(RegisteredOptions::BqpdCategory, GetRawPtr(ampl_options_list));
00171     fillApplicationOptions(GetRawPtr(ampl_options_list) );
00172     std::string options_id = appName + "_options";
00173     ampl_tnlp_ = new AmplTNLP(jnlst, options, argv, suffix_handler, true,
00174         ampl_options_list, options_id.c_str(),
00175         appName.c_str(), appName.c_str(), nl_file_content);
00176 
00177 
00178     /* Read suffixes */
00179     read_obj_suffixes();
00180     read_priorities();
00181     read_convexities();
00182     read_sos();
00183 
00184 
00185     /* Determine if objective is linear.*/
00186     Index n_non_linear_b= 0;
00187     Index n_non_linear_bi= 0;
00188     Index n_non_linear_c= 0;
00189     Index n_non_linear_ci = 0;
00190     Index n_non_linear_o= 0;
00191     Index n_non_linear_oi = 0;
00192     Index n_binaries = 0;
00193     Index n_integers = 0;
00194     ampl_tnlp_->get_discrete_info(n_non_linear_b, n_non_linear_bi, n_non_linear_c,
00195         n_non_linear_ci, n_non_linear_o, n_non_linear_oi,
00196         n_binaries, n_integers);
00197     if (n_non_linear_b == 0 && n_non_linear_o == 0) {
00198       hasLinearObjective_ = true;
00199     }
00200   }
00201 
00202   AmplTMINLP::~AmplTMINLP()
00203   {
00204     delete [] constraintsConvexities_;
00205     delete [] simpleConcaves_;
00206     delete [] nonConvexConstraintsAndRelaxations_;
00207     delete ampl_tnlp_;
00208   }
00209 
00210   void
00211   AmplTMINLP::read_priorities()
00212   {
00213     int numcols, m, dummy1, dummy2;
00214     TNLP::IndexStyleEnum index_style;
00215     ampl_tnlp_->get_nlp_info(numcols, m, dummy1, dummy2, index_style);
00216 
00217     const AmplSuffixHandler * suffix_handler = GetRawPtr(suffix_handler_);
00218 
00219     const Index* pri = suffix_handler->GetIntegerSuffixValues("priority", AmplSuffixHandler::Variable_Source);
00220     const Index* brac = suffix_handler->GetIntegerSuffixValues("direction", AmplSuffixHandler::Variable_Source);
00221     const Number* upPs = suffix_handler->GetNumberSuffixValues("upPseudocost", AmplSuffixHandler::Variable_Source);
00222     const Number* dwPs = suffix_handler->GetNumberSuffixValues("downPseudocost", AmplSuffixHandler::Variable_Source);
00223 
00224 
00225     branch_.gutsOfDestructor();
00226     branch_.size = numcols;
00227     if (pri) {
00228       branch_.priorities = new int[numcols];
00229       for (int i = 0 ; i < numcols ; i++) {
00230         branch_.priorities [i] = -pri[i] + 9999;
00231       }
00232     }
00233     if (brac) {
00234       branch_.branchingDirections = CoinCopyOfArray(brac,numcols);
00235     }
00236     if (upPs && !dwPs) dwPs = upPs;
00237     else if (dwPs && !upPs) upPs = dwPs;
00238 
00239     if (upPs) {
00240       branch_.upPsCosts = CoinCopyOfArray(upPs,numcols);
00241     }
00242     if (dwPs) {
00243       branch_.downPsCosts = CoinCopyOfArray(dwPs,numcols);
00244     }
00245 
00246     const double* perturb_radius =
00247       suffix_handler->GetNumberSuffixValues("perturb_radius", AmplSuffixHandler::Variable_Source);
00248     perturb_info_.SetPerturbationArray(numcols, perturb_radius);
00249   }
00250 
00251   void
00252   AmplTMINLP::read_sos()
00253   {
00254     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00255     int i = 0;//ASL_suf_sos_explict_free;
00256     int copri[2], **p_sospri;
00257     copri[0] = 0;
00258     copri[1] = 0;
00259     int * starts = NULL;
00260     int * indices = NULL;
00261     char * types = NULL;
00262     double * weights = NULL;
00263     int * priorities = NULL;
00264     p_sospri = &priorities;
00265 
00266     sos_.gutsOfDestructor();
00267     int m = n_con;
00268     sos_.num = suf_sos(i, &sos_.numNz, &types, p_sospri, copri,
00269         &starts, &indices, &weights);
00270     if(m != n_con){
00271       throw CoinError("number of constraints changed by suf_sos. Not supported.",
00272                        "read_sos","Bonmin::AmplTMINLP");
00273    }
00274     if (sos_.num) {
00275       //Copy sos information
00276       sos_.priorities = CoinCopyOfArray(priorities,sos_.num);
00277       sos_.starts = CoinCopyOfArray(starts, sos_.num + 1);
00278       sos_.indices = CoinCopyOfArray(indices, sos_.numNz);
00279       sos_.types = CoinCopyOfArray(types, sos_.num);
00280       sos_.weights = CoinCopyOfArray(weights, sos_.numNz);
00281 
00282       ampl_utils::sos_kludge(sos_.num, sos_.starts, sos_.weights);
00283       for (int ii=0;ii<sos_.num;ii++) {
00284         int ichar = sos_.types[ii] - '0';
00285         if (ichar != 1 && ichar != 2) {
00286           std::cerr<<"Unsuported type of sos constraint: "<<sos_.types[ii]<<std::endl;
00287           throw;
00288         }
00289         sos_.types[ii]= static_cast<char>(ichar);
00290       }
00291     }
00292   }
00293 
00294   void
00295   AmplTMINLP::read_obj_suffixes()
00296   {
00297     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00298     DBG_ASSERT(asl);
00299     //Get the values
00300     if (n_obj < 2) return;
00301 
00302     const AmplSuffixHandler * suffix_handler = GetRawPtr(suffix_handler_);
00303 
00304     const Index* UBObj = suffix_handler->GetIntegerSuffixValues("UBObj", AmplSuffixHandler::Objective_Source);
00305     if (UBObj) {
00306       //get index of lower bounding objective
00307       int lowerBoundingObj = (UBObj[0] == 1)? 1:0;
00308       // Pass information to Ipopt
00309       ampl_tnlp_->set_active_objective(lowerBoundingObj);
00310 
00311       //get the index of upper bounding objective
00312       for (int i = 0; i < n_obj; i++) {
00313         if (UBObj[i]==1) {
00314           if (upperBoundingObj_ != -1) {
00315             jnlst_->Printf(J_ERROR, J_MAIN,
00316                 "Too many objectives for upper-bounding");
00317           }
00318           upperBoundingObj_ = i;
00319         }
00320       }
00321     }
00322     else {
00323       ampl_tnlp_->set_active_objective(0);
00324     }
00325   }
00327   struct NonConvexSuff
00328   {
00330     NonConvexSuff():
00331         cIdx(-1),relIdx(-1),scXIdx(-1),scYIdx(-1)
00332     {}
00334     NonConvexSuff(const NonConvexSuff& other):
00335         cIdx(other.cIdx), relIdx(other.relIdx),
00336         scXIdx(other.scXIdx), scYIdx(other.scYIdx)
00337     {}
00339     int cIdx;
00341     int relIdx;
00343     int scXIdx;
00345     int scYIdx;
00346   };
00347   void AmplTMINLP::read_convexities()
00348   {
00349     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00350     DBG_ASSERT(asl);
00351 
00352     const AmplSuffixHandler * suffix_handler = GetRawPtr(suffix_handler_);
00353     const Index * id = suffix_handler->GetIntegerSuffixValues("id", AmplSuffixHandler::Variable_Source);
00354     const Index * primary_var = suffix_handler->GetIntegerSuffixValues("primary_var", AmplSuffixHandler::Constraint_Source);
00355 
00356 
00357     if (primary_var!= NULL) {
00358       if (constraintsConvexities_ != NULL) {
00359         delete [] constraintsConvexities_;
00360       }
00361       constraintsConvexities_ = new TMINLP::Convexity[n_con];
00362       if (id == NULL) {
00363         std::cerr<<"Incorrect suffixes description in ampl model. id's are not declared "<<std::endl;
00364         exit(ERROR_IN_AMPL_SUFFIXES);
00365       }
00366       int numberSimpleConcave = 0;
00367       std::map<int, int> id_map;
00368 
00369       for (int i = 0 ; i < n_var ; i++) {
00370         id_map[id[i]] = i;
00371       }
00372 
00373 
00374       for (int i = 0 ; i < n_con ; i++) {
00375         if (primary_var[i] != 0) {
00376           constraintsConvexities_[i] = TMINLP::SimpleConcave;
00377           numberSimpleConcave++;
00378         }
00379         else constraintsConvexities_[i] = TMINLP::Convex;
00380       }
00381       simpleConcaves_ = new SimpleConcaveConstraint[numberSimpleConcave];
00382       nonConvexConstraintsAndRelaxations_ = new MarkedNonConvex[numberSimpleConcave];
00383       numberSimpleConcave = 0;
00384       int * jCol = new int[n_var];
00385       for (int i = 0 ; i < n_con ; i++) {
00386         if (primary_var[i] != 0) {
00387           nonConvexConstraintsAndRelaxations_[numberSimpleConcave].cIdx = i;
00388           nonConvexConstraintsAndRelaxations_[numberSimpleConcave].cRelaxIdx = -1;
00389           simpleConcaves_[numberSimpleConcave].cIdx = i;
00390           simpleConcaves_[numberSimpleConcave].yIdx = id_map[primary_var[i]];
00391 
00392           //Now get gradient of i to get xIdx.
00393           int nnz;
00394           int & yIdx = simpleConcaves_[numberSimpleConcave].yIdx;
00395           int & xIdx = simpleConcaves_[numberSimpleConcave].xIdx;
00396           eval_grad_gi(n_var, NULL, false, i, nnz, jCol, NULL);
00397           if (nnz != 2) {//Error in ampl model
00398             std::cout<<"Incorrect suffixes description in ampl model. Constraint with id "
00399             <<id<<" is simple concave and should have only two nonzero elements"<<std::endl;
00400             exit(ERROR_IN_AMPL_SUFFIXES);
00401           }
00402           if (jCol[0] - 1== yIdx) {
00403             xIdx = jCol[1] - 1;
00404           }
00405           else {
00406             if (jCol[1] - 1!= yIdx) {//Error in ampl model
00407               std::cout<<"Incorrect suffixes description in ampl model. Constraint with id "
00408               <<id<<" : variable marked as y does not appear in the constraint."<<std::endl;
00409               exit(ERROR_IN_AMPL_SUFFIXES);
00410             }
00411             xIdx = jCol[0] - 1;
00412           }
00413           numberSimpleConcave++;
00414         }
00415       }
00416       delete [] jCol;
00417       numberSimpleConcave_ = numberSimpleConcave;
00418       numberNonConvex_ = numberSimpleConcave;
00419     }
00420 
00421   }
00422 
00423   bool AmplTMINLP::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g, Index& nnz_h_lag, TNLP::IndexStyleEnum& index_style)
00424   {
00425     return ampl_tnlp_->get_nlp_info(n, m, nnz_jac_g, nnz_h_lag, index_style);
00426   }
00427 
00428   bool AmplTMINLP::get_variables_types(Index n, VariableType* var_types)
00429   {
00430     // The variables are sorted by type in AMPL, so all we need to
00431     // know are the counts of each type.
00432 
00433 
00434     Index n_non_linear_b= 0;
00435     Index n_non_linear_bi= 0;
00436     Index n_non_linear_c= 0;
00437     Index n_non_linear_ci = 0;
00438     Index n_non_linear_o= 0;
00439     Index n_non_linear_oi = 0;
00440     Index n_binaries = 0;
00441     Index n_integers = 0;
00442     ampl_tnlp_->get_discrete_info(n_non_linear_b, n_non_linear_bi, n_non_linear_c,
00443         n_non_linear_ci, n_non_linear_o, n_non_linear_oi,
00444         n_binaries, n_integers);
00445     int totalNumberOfNonContinuous = 0;
00446 
00447 
00448     //begin with variables non-linear both in the objective function and in some constraints
00449     // The first ones are continuous
00450     Index start = 0;
00451     Index end = n_non_linear_b - n_non_linear_bi;
00452     for (int i=start; i<end; i++) {
00453       var_types[i] = CONTINUOUS;
00454     }
00455 
00456     // The second ones are integers
00457     start = end;
00458     end = start + n_non_linear_bi;
00459     for (int i=start; i < end; i++) {
00460       var_types[i] = INTEGER;
00461       totalNumberOfNonContinuous++;
00462     }
00463     //next variables non-linear in some constraints
00464     // The first ones are continuous
00465     start = end;
00466     end = std::max(start,end + n_non_linear_c - n_non_linear_ci - n_non_linear_b);
00467     for (int i=start; i<end; i++) {
00468       var_types[i] = CONTINUOUS;
00469     }
00470 
00471     // The second ones are integers
00472     start = end;
00473     end = start + n_non_linear_ci;
00474     for (int i=start; i < end; i++) {
00475       var_types[i] = INTEGER;
00476       totalNumberOfNonContinuous++;
00477     }
00478 
00479     //next variables non-linear in the objective function
00480     // The first ones are continuous
00481     start = end;
00482     end = std::max(start,end + n_non_linear_o - std::max(n_non_linear_b, n_non_linear_c) - n_non_linear_oi);
00483     for (int i=start; i<end; i++) {
00484       var_types[i] = CONTINUOUS;
00485     }
00486 
00487     // The second ones are integers
00488     start = end;
00489     end = start + n_non_linear_oi;
00490     for (int i=start; i < end; i++) {
00491       var_types[i] = INTEGER;
00492       totalNumberOfNonContinuous++;
00493     }
00494 
00495     //At last the linear variables
00496     // The first ones are continuous
00497     start = end;
00498     end = n - n_binaries - n_integers;
00499     for (int i=start; i<end; i++) {
00500       var_types[i] = CONTINUOUS;
00501     }
00502 
00503     // The second ones are binaries
00504     start = end;
00505     end = start + n_binaries;
00506     for (int i=start; i < end; i++) {
00507       var_types[i] = BINARY;
00508       totalNumberOfNonContinuous++;
00509     }
00510 
00511     // The third ones are integers
00512     start = end;
00513     end = start + n_integers;
00514     for (int i=start; i < end; i++) {
00515       var_types[i] = INTEGER;
00516       totalNumberOfNonContinuous++;
00517     }
00518     return true;
00519   }
00520 
00521   bool AmplTMINLP::get_variables_linearity(Index n, Ipopt::TNLP::LinearityType* var_types)
00522   {
00523     // The variables are sorted by type in AMPL, so all we need to
00524     // know are the counts of each type.
00525 
00526 
00527     Index n_non_linear_b= 0;
00528     Index n_non_linear_bi= 0;
00529     Index n_non_linear_c= 0;
00530     Index n_non_linear_ci = 0;
00531     Index n_non_linear_o= 0;
00532     Index n_non_linear_oi = 0;
00533     Index n_binaries = 0;
00534     Index n_integers = 0;
00535     ampl_tnlp_->get_discrete_info(n_non_linear_b, n_non_linear_bi, n_non_linear_c,
00536         n_non_linear_ci, n_non_linear_o, n_non_linear_oi,
00537         n_binaries, n_integers);
00538 
00539     //Compute the number of non linear variables:
00540     int n_non_linear = std::max(n_non_linear_c, n_non_linear_o);//n_non_linear_c + n_non_linear_o - n_non_linear_b;
00541 
00542     //printf("n_non_linear_c %i n_non_linear_o %i n_non_linear_b %i\n", n_non_linear_c, n_non_linear_o, n_non_linear_b);
00543 
00544     int start = 0;
00545     int end = n_non_linear;
00546     for (int i=start; i<end; i++) {
00547       var_types[i] = Ipopt::TNLP::NON_LINEAR;
00548     }
00549 
00550     //At last the linear variables
00551     // The first ones are continuous
00552     start = end;
00553     end = n;
00554     for (int i=start; i<end; i++) {
00555       var_types[i] = Ipopt::TNLP::LINEAR;
00556     }
00557     return true;
00558   }
00559 
00562   bool
00563   AmplTMINLP::get_constraints_linearity(Index m,
00564       Ipopt::TNLP::LinearityType* const_types)
00565   {
00566     return ampl_tnlp_->get_constraints_linearity(m, const_types);
00567   }
00568 
00569   bool AmplTMINLP::get_bounds_info(Index n, Number* x_l, Number* x_u,
00570       Index m, Number* g_l, Number* g_u)
00571   {
00572     return ampl_tnlp_->get_bounds_info(n, x_l, x_u, m, g_l, g_u);
00573   }
00574 
00575   bool AmplTMINLP::get_starting_point(Index n, bool init_x, Number* x,
00576       bool init_z, Number* z_L, Number* z_U,
00577       Index m, bool init_lambda, Number* lambda)
00578   {
00579     return ampl_tnlp_->get_starting_point(n, init_x, x,
00580         init_z, z_L, z_U,
00581         m, init_lambda, lambda);
00582   }
00583 
00584   bool AmplTMINLP::eval_f(Index n, const Number* x, bool new_x, Number& obj_value)
00585   {
00586     return ampl_tnlp_->eval_f(n, x, new_x, obj_value);
00587   }
00588 
00589   bool AmplTMINLP::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f)
00590   {
00591     return ampl_tnlp_->eval_grad_f(n, x, new_x, grad_f);
00592   }
00593 
00594   bool AmplTMINLP::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g)
00595   {
00596     return ampl_tnlp_->eval_g(n, x, new_x, m, g);
00597   }
00598 
00599   bool AmplTMINLP::eval_jac_g(Index n, const Number* x, bool new_x,
00600       Index m, Index nele_jac, Index* iRow,
00601       Index *jCol, Number* values)
00602   {
00603     return ampl_tnlp_->eval_jac_g(n, x, new_x,
00604         m, nele_jac, iRow, jCol,
00605         values);
00606   }
00607 
00608   bool AmplTMINLP::eval_h(Index n, const Number* x, bool new_x,
00609       Number obj_factor, Index m, const Number* lambda,
00610       bool new_lambda, Index nele_hess, Index* iRow,
00611       Index* jCol, Number* values)
00612   {
00613     return ampl_tnlp_->eval_h(n, x, new_x, obj_factor,
00614         m, lambda, new_lambda, nele_hess, iRow,
00615         jCol, values);
00616   }
00617 
00618   bool AmplTMINLP::eval_gi(Index n, const Number* x, bool new_x,
00619       Index i, Number& gi)
00620   {
00621     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00622 
00623     // ignore new_x for now
00624     xunknown();
00625 
00626     fint nerror = 0;
00627     gi = conival(i, const_cast<real*>(x), &nerror);
00628     if (nerror!=0) {
00629       return false;
00630     }
00631     else {
00632       return true;
00633     }
00634   }
00635 
00636   bool AmplTMINLP::eval_grad_gi(Index n, const Number* x, bool new_x,
00637       Index i, Index& nele_grad_gi, Index* jCol,
00638       Number* values)
00639   {
00640     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00641 
00642     if (jCol) {
00643       // Only compute the number of nonzeros and the indices
00644       DBG_ASSERT(!values);
00645       nele_grad_gi = 0;
00646       for (cgrad* cg=Cgrad[i]; cg; cg = cg->next) {
00647         jCol[nele_grad_gi++] = cg->varno + 1;
00648       }
00649       return true;
00650     }
00651     DBG_ASSERT(values);
00652 
00653     // ignore new_x for now
00654     xunknown();
00655 
00656     asl->i.congrd_mode = 1;
00657     fint nerror = 0;
00658     congrd(i, const_cast<real*>(x), values, &nerror);
00659     if (nerror!=0) {
00660       return false;
00661     }
00662     else {
00663       return true;
00664     }
00665   }
00666 
00667   void AmplTMINLP::finalize_solution(TMINLP::SolverReturn status,
00668       Index n, const Number* x, Number obj_value)
00669   {
00670     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00671     std::string message;
00672     std::string status_str;
00673     if (status == TMINLP::SUCCESS) {
00674       status_str = "\t\"Finished\"";
00675       message = "\n" + appName_ +": Optimal";
00676       solve_result_num = 3;
00677     }
00678     else if (status == TMINLP::INFEASIBLE) {
00679       status_str = "\t\"Finished\"";
00680       message = "\n" + appName_ + ": Infeasible problem";
00681       solve_result_num = 220;
00682     }
00683     else if (status == TMINLP::CONTINUOUS_UNBOUNDED) {
00684       status_str = "\t\"Finished\"";
00685       message = "\n" + appName_ +" Continuous relaxation is unbounded.";
00686       solve_result_num = 300;
00687     }
00688     else if (status == TMINLP::LIMIT_EXCEEDED) {
00689       status_str = "\t\"Not finished\"";
00690       message = "\n" + appName_ + ": Optimization interupted on limit.";
00691       if(x)
00692         solve_result_num = 421; /* Limit reached with integer feasible solution.*/
00693       else
00694         solve_result_num = 410; /* Limit reached without solution.*/
00695     }
00696     else if (status == TMINLP::MINLP_ERROR) {
00697       status_str = "\t\"Aborted\"";
00698       message = "\n" + appName_ + ": Error encountered in optimization.";
00699       solve_result_num = 500;
00700     }
00701     if (writeAmplSolFile_) {
00702       write_solution(message, x);
00703       std::cout<<"\n "<<status_str<<std::endl;
00704     }
00705     else {
00706       std::cout<<status_str<<message<<std::endl;
00707       std::string fName = appName_ + ".sol";
00708       std::ofstream of(fName.c_str());
00709       for (int i = 0 ; i < n ; i++) {
00710         of<<i<<"\t"<<x[i]<<std::endl;
00711       }
00712       of<<"-1\n";
00713     }
00714   }
00715 
00716   void AmplTMINLP::write_solution(const std::string & message, const Number *x_sol)
00717   {
00718     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00719 
00720     DBG_ASSERT(asl);
00721     //    DBG_ASSERT(x_sol);
00722 
00723     // We need to copy the message into a non-const char array to make
00724     // it work with the AMPL C function.
00725     char* cmessage = new char[message.length()+1];
00726     strcpy(cmessage, message.c_str());
00727 
00728     // In order to avoid casting into non-const, we copy the data here...
00729     double* x_sol_copy = NULL;
00730     if (x_sol) {
00731       x_sol_copy = new double[n_var];
00732       for (int i=0; i<n_var; i++) {
00733         x_sol_copy[i] = x_sol[i];
00734       }
00735     }
00736     write_sol(cmessage, x_sol_copy, NULL, NULL);
00737 
00738     delete [] x_sol_copy;
00739     delete [] cmessage;
00740 
00741   }
00742 
00743 
00745   void AmplTMINLP::getLinearPartOfObjective(double * obj)
00746   {
00747     Index n, m, nnz_jac_g, nnz_h_lag;
00748     TNLP::IndexStyleEnum index_style = TNLP::FORTRAN_STYLE;
00749     get_nlp_info( n, m, nnz_jac_g, nnz_h_lag, index_style);
00750     eval_grad_f(n, NULL, 0,obj);
00751     Index n_non_linear_b= 0;
00752     Index n_non_linear_bi= 0;
00753     Index n_non_linear_c= 0;
00754     Index n_non_linear_ci = 0;
00755     Index n_non_linear_o= 0;
00756     Index n_non_linear_oi = 0;
00757     Index n_binaries = 0;
00758     Index n_integers = 0;
00759     ampl_tnlp_->get_discrete_info(n_non_linear_b, n_non_linear_bi, n_non_linear_c,
00760         n_non_linear_ci, n_non_linear_o, n_non_linear_oi,
00761         n_binaries, n_integers);
00762 
00763     // Now get the coefficients of variables wich are linear in the objective
00764     // The first ones are not
00765     Index start = 0;
00766     Index end = n_non_linear_b;
00767     for (int i=start; i<end; i++) {
00768       obj[i] = 0.;
00769     }
00770 
00771     //next variables should be linear in the objective just skip them
00772     // (from current end to (end + n_non_linear_c - n_non_linear_ci - n_non_linear_b;)
00773 
00774 
00775     // Those are nonlinear in the objective
00776     start = end + n_non_linear_c;
00777     end = start + n_non_linear_o;
00778     for (int i=start; i < end; i++) {
00779       obj[i]=0.;
00780     }
00781     //The rest is linear keep the values of the gradient
00782   }
00783 
00784 
00785 
00788   bool
00789   AmplTMINLP::eval_upper_bound_f(Index n, const Number* x,
00790       Number& obj_value)
00791   {
00792     ASL_pfgh* asl = ampl_tnlp_->AmplSolverObject();
00793     //xknown(x);    // This tells ampl to use a new x
00794     fint nerror = -1;
00795     obj_value = objval(upperBoundingObj_, const_cast<double *>(x), &nerror);
00796     if (nerror > 0) {
00797       jnlst_->Printf(J_ERROR, J_MAIN,
00798           "Error in evaluating upper bounding objecting");
00799       throw -1;
00800     }
00801     return nerror;
00802 
00803   }
00805   const ASL_pfgh*
00806   AmplTMINLP::AmplSolverObject() const
00807   {
00808     return ampl_tnlp_->AmplSolverObject();
00809   }
00810 
00811 } // namespace Bonmin

Generated on Wed Nov 30 03:03:55 2011 by  doxygen 1.4.7