Changeset 42


Ignore:
Timestamp:
Sep 15, 2006 11:32:56 PM (13 years ago)
Author:
ladanyi
Message:

B&B running fine; preparing to distribute pseudocosts

Location:
trunk/Bonmin/experimental/Bcp
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • trunk/Bonmin/experimental/Bcp/BM.cpp

    r30 r42  
    6868}
    6969
     70/****************************************************************************/
     71
    7072template <>
    7173void BCP_parameter_set<BM_par>::set_default_entries() {
     
    113115//#############################################################################
    114116
    115 void
    116 BM_tm::readIpopt()
    117 {
    118     if (par.entry(BM_par::IpoptParamfile).length() != 0) {
    119         FILE* ipopt = fopen(par.entry(BM_par::IpoptParamfile).c_str(), "r");
    120         if (!ipopt) {
    121             throw BCP_fatal_error("Non-existent IpoptParamfile\n");
    122         }
    123         fseek(ipopt, 0, SEEK_END);
    124         int len = ftell(ipopt) + 1;
    125         fseek(ipopt, 0, SEEK_SET);
    126         char* ipopt_content = new char[len+1];
    127         len = fread(ipopt_content, 1, len, ipopt);
    128         ipopt_file_content.assign(ipopt_content, len);
    129         delete[] ipopt_content;
    130     }
    131 
    132     FILE* nl = fopen(par.entry(BM_par::NL_filename).c_str(), "r");
    133     if (!nl) {
    134         throw BCP_fatal_error("Non-existent NL_filename\n");
    135     }
    136     fseek(nl, 0, SEEK_END);
    137     int len = ftell(nl) + 1;
    138     fseek(nl, 0, SEEK_SET);
    139     char* nl_content = new char[len+1];
    140     len = fread(nl_content, 1, len, nl);
    141     nl_file_content.assign(nl_content, len);
    142     delete[] nl_content;
    143 }
    144 
    145 /****************************************************************************/
    146 
    147 void
    148 BM_tm::initialize_core(BCP_vec<BCP_var_core*>& vars,
    149                        BCP_vec<BCP_cut_core*>& cuts,
    150                        BCP_lp_relax*& matrix)
    151 {
    152     BonminAmplInterface nlpSolver;
    153     char* argv_[3];
    154     char** argv = argv_;
    155     argv[0] = NULL;
    156     argv[1] = strdup(par.entry(BM_par::NL_filename).c_str());
    157     argv[2] = NULL;
    158     nlpSolver.readAmplNlFile(argv, 0, 0);
    159     free(argv[1]);
    160  
    161     nlpSolver.extractInterfaceParams();
    162  
    163     OsiClpSolverInterface clp;
    164     int addObjVar = par.entry(BM_par::PureBranchAndBound) ? 0 : 1;
    165     nlpSolver.extractLinearRelaxation(clp, addObjVar);
    166  
    167     const int numCols = clp.getNumCols();
    168     const int numRows = clp.getNumRows();
    169 
    170     double* obj = new double[numCols];
    171     if (par.entry(BM_par::PureBranchAndBound)) {
    172         CoinFillN(obj, numCols, 0.0);
    173     }
    174     else {
    175         CoinDisjointCopyN(clp.getObjCoefficients(), numCols, obj);
    176     }
    177 
    178     vars.reserve(numCols);
    179     const double* clb = clp.getColLower();
    180     const double* cub = clp.getColUpper();
    181     for (int i = 0; i < numCols; ++i)
    182         {
    183             BCP_var_t type = BCP_ContinuousVar;
    184             if (clp.isBinary(i)) type = BCP_BinaryVar;
    185             if (clp.isInteger(i)) type = BCP_IntegerVar;
    186             vars.push_back(new BCP_var_core(type, obj[i], clb[i], cub[i]));
    187         }
    188 
    189     if (par.entry(BM_par::PureBranchAndBound)) {
    190         // Just fake something into the core matrix. In this case: 0 <= 1
    191         BCP_vec<double> OBJ(obj, numCols);
    192         BCP_vec<double> CLB(clb, numCols);
    193         BCP_vec<double> CUB(cub, numCols);
    194         BCP_vec<double> RLB(1, 0.0);
    195         BCP_vec<double> RUB(1, 1.0);
    196         BCP_vec<int> VB;
    197         BCP_vec<int> EI;
    198         BCP_vec<double> EV;
    199         VB.push_back(0);
    200         VB.push_back(2);
    201         EI.push_back(0);         EV.push_back(0.0);
    202         EI.push_back(numCols-1); EV.push_back(0.0);
    203         matrix = new BCP_lp_relax(false, VB, EI, EV, OBJ, CLB, CUB, RLB, RUB);
    204         cuts.push_back(new BCP_cut_core(0.0, 1.0));
    205     } else {
    206         cuts.reserve(numRows);
    207         const double* rlb = clp.getRowLower();
    208         const double* rub = clp.getRowUpper();
    209         for (int i = 0; i < numRows; ++i)
    210             {
    211                 cuts.push_back(new BCP_cut_core(rlb[i], rub[i]));
    212             }
    213 
    214         matrix = new BCP_lp_relax(true /*column major ordered*/);
    215         matrix->copyOf(*clp.getMatrixByCol(), obj, clb, cub, rlb, rub);
    216     }
    217 }
    218 
    219 /****************************************************************************/
    220 void
    221 BM_tm::create_root(BCP_vec<BCP_var*>& added_vars,
    222                    BCP_vec<BCP_cut*>& added_cuts,
    223                    BCP_user_data*& user_data,
    224                    BCP_pricing_status& pricing_status)
    225 {
    226     BM_node* data = new BM_node;
    227     data->numNlpFailed_ = 0;
    228     user_data = data;
    229 }
    230 
    231 /****************************************************************************/
    232 
    233 void
    234 BM_tm::display_feasible_solution(const BCP_solution* sol)
    235 {
    236     const BM_solution* bs = dynamic_cast<const BM_solution*>(sol);
    237     if (!bs) {
    238         throw BCP_fatal_error("Trying to pack non-BM_solution.\n");
    239     }
    240 
    241     /* Parse again the input file so that we have a nice and clean ampl
    242        setup */
    243     BonminAmplInterface nlpSolver; 
    244     char* argv_[3];
    245     char** argv = argv_;
    246     argv[0] = NULL;
    247     argv[1] = strdup(par.entry(BM_par::NL_filename).c_str());
    248     argv[2] = NULL;
    249     nlpSolver.readAmplNlFile(argv, 0, 0);
    250     free(argv[1]);
    251     OsiClpSolverInterface clp;
    252     int addObjVar = par.entry(BM_par::PureBranchAndBound) ? 0 : 1;
    253     nlpSolver.extractLinearRelaxation(clp, addObjVar);
    254     const int numCols = clp.getNumCols();
    255 
    256     int i;
    257  
    258     /* This will give the vector of core variables we have created in
    259        BM_tm::initialize_core */
    260     const BCP_vec<BCP_var_core*>& vars = getTmProblemPointer()->core->vars;
    261 
    262     /* Create a dense vector with the value of each variable. Round the
    263        integer vars (they were tested to be near enough to integrality) so
    264        in the printouts we won't have 9.99999991234, etc.) */
    265     double* dsol = new double[numCols];
    266     for (i = 0; i < numCols; ++i) {
    267         dsol[i] = 0.0;
    268     }
    269     const int size = bs->_ind.size();
    270     for (i = 0; i < size; ++i) {
    271         const int ind = bs->_ind[i];
    272         const double v = bs->_values[i];
    273         const BCP_var_t type = vars[ind]->var_type();
    274         dsol[ind] = (type == BCP_ContinuousVar) ? v : floor(v+0.5);
    275     }
    276 
    277     /* Display the solution on the screen */
    278     printf("bonmin: feasible solution found.  Objective value: %f\n",
    279            bs->_objective);
    280     for (i = 0; i < size; ++i) {
    281         printf("    index: %5i   value: %f\n", bs->_ind[i], dsol[bs->_ind[i]]);
    282     }
    283     printf("\n");
    284 
    285     /* create the AMPL solfile */
    286     nlpSolver.writeAmplSolFile("\nbon-min: Optimal solution", dsol, NULL);
    287     delete[] dsol;
    288 }
    289 
    290 //#############################################################################
    291 
    292 BM_lp::BM_lp() :
    293     BCP_lp_user(),
    294     sos(),
    295     babSolver_(3),
    296     nlp(),
    297     feasChecker_(NULL),
    298     in_strong(0)
    299 {
    300     babSolver_.setSolver(&nlp);
    301     setOsiBabSolver(&babSolver_);
    302 }
    303 
    304 /****************************************************************************/
    305 
    306 BM_lp::~BM_lp()
    307 {
    308     delete feasChecker_;
    309     delete[] primal_solution_;
    310     delete[] branching_priority_;
    311     for (int i = sos.size() - 1; i >= 0; --i) {
    312         delete sos[i];
    313     }
    314 }
    315 
    316 /****************************************************************************/
    317 
    318 OsiSolverInterface *
    319 BM_lp::initialize_solver_interface()
    320 {
    321     OsiClpSolverInterface * clp = new OsiClpSolverInterface;
    322     OsiBabSolver babSolver(3);
    323     babSolver.setSolver(clp);
    324     clp->setAuxiliaryInfo(&babSolver);
    325     clp->messageHandler()->setLogLevel(0);
    326     setOsiBabSolver(dynamic_cast<OsiBabSolver *>(clp->getAuxiliaryInfo()));
    327     return clp;
    328 }
    329 
    330 /****************************************************************************/
    331 
    332 void
    333 BM_lp::initialize_new_search_tree_node(const BCP_vec<BCP_var*>& vars,
    334                                        const BCP_vec<BCP_cut*>& cuts,
    335                                        const BCP_vec<BCP_obj_status>& vs,
    336                                        const BCP_vec<BCP_obj_status>& cs,
    337                                        BCP_vec<int>& var_changed_pos,
    338                                        BCP_vec<double>& var_new_bd,
    339                                        BCP_vec<int>& cut_changed_pos,
    340                                        BCP_vec<double>& cut_new_bd)
    341 {
    342     int i;
    343     // First copy the bounds into nlp. That way all the branching decisions
    344     // will be transferred over.
    345     for (i = sos.size() - 1; i >= 0; --i) {
    346         sos[i]->first = 0;
    347         sos[i]->last = sos[i]->length;
    348     }
    349     OsiSolverInterface * osi = getLpProblemPointer()->lp_solver;
    350     const int numCols = osi->getNumCols();
    351     const double* clb = osi->getColLower();
    352     const double* cub = osi->getColUpper();
    353     for (i = 0; i < numCols; ++i) {
    354         const BCP_var_core* v =
    355             dynamic_cast<const BCP_var_core*>(vars[i]);
    356         if (v) {
    357             nlp.setColLower(i, clb[i]);
    358             nlp.setColUpper(i, cub[i]);
    359             continue;
    360         }
    361         const BM_branching_var* bv =
    362             dynamic_cast<const BM_branching_var*>(vars[i]);
    363         if (bv) {
    364             const int ind = bv->sos_index;
    365             const int split = bv->split;
    366             const char type = sos[ind]->type;
    367             const int *indices = sos[ind]->indices;
    368             if (bv->ub() == 0.0) {
    369                 const int last = sos[ind]->last;
    370                 for (int j = split + type; j < last; ++j) {
    371                     nlp.setColLower(indices[j], 0.0);
    372                     nlp.setColUpper(indices[j], 0.0);
    373                 }
    374                 sos[ind]->last = split + type;
    375             } else {
    376                 const int first = sos[ind]->first;
    377                 for (int j = first; j < split; ++j) {
    378                     nlp.setColLower(indices[j], 0.0);
    379                     nlp.setColUpper(indices[j], 0.0);
    380                 }
    381                 sos[ind]->first = split;
    382             }
    383             continue;
    384         }
    385         throw BCP_fatal_error("BM_lp: unrecognized variable type\n");
    386     }
    387 
    388     in_strong = 0;
    389 }
    390 
    391 /************************************************************************/
    392 void
    393 BM_lp::modify_lp_parameters(OsiSolverInterface* lp, bool in_strong_branching)
    394 
    395     // Called each time the node LP is solved
    396 
    397 {
    398     if (in_strong_branching) {
    399         in_strong = 1;
    400         //     lp->setIntParam(OsiMaxNumIterationHotStart, 50);
    401     }
    402 }
    403 
    404 /****************************************************************************/
    405 
    406 BCP_solution*
    407 BM_lp::test_feasibility(const BCP_lp_result& lp_result,
    408                         const BCP_vec<BCP_var*>& vars,
    409                         const BCP_vec<BCP_cut*>& cuts)
    410 {
    411     BM_solution* sol = NULL;
    412 
    413     if (par.entry(BM_par::PureBranchAndBound)) {
    414         /* PURE_BB TODO:
    415            Do whatever must be done to:
    416            1) compute a lower bound using NLP and save it in "lower_bound_"
    417            2) save the solution of the NLP in "primal_solution_"
    418            3) if the optimum (well, local opt anyway...) happens to be
    419            feasible (or along the NLP optimization we have blundered into
    420            a feas sol) then create "sol"
    421         */
    422         nlp.initialSolve();
    423         if (nlp.isProvenOptimal()) {
    424             const int numCols = nlp.getNumCols();
    425             lower_bound_ = nlp.getObjValue();
    426             CoinDisjointCopyN(nlp.getColSolution(), numCols, primal_solution_);
    427             numNlpFailed_ = 0;
    428             Ipopt::SmartPtr<Ipopt::OptionsList> options = nlp.retrieve_options();
    429             double intTol;
    430             options->GetNumericValue("integer_tolerance",intTol,"bonmin.");
    431 
    432             int i;
    433             for (i = 0; i < numCols; ++i) {
    434                 if (vars[i]->var_type() == BCP_ContinuousVar)
    435                     continue;
    436                 const double psol = CoinMin(CoinMax(vars[i]->lb(),
    437                                                     primal_solution_[i]),
    438                                             vars[i]->ub());
    439                 const double frac = psol - floor(psol);
    440                 const double dist = std::min(frac, 1.0-frac);
    441                 if (dist >= intTol)
    442                     break;
    443             }
    444             if (i == numCols) {
    445                 /* yipee! a feasible solution! */
    446                 sol = new BM_solution;
    447                 //Just copy the solution stored in solver to sol
    448                 for (i = 0 ; i < numCols ; i++) {
    449                     if (primal_solution_[i] > lp_result.primalTolerance())
    450                         sol->add_entry(i, primal_solution_[i]);
    451                 }
    452                 sol->setObjective(lower_bound_);
    453             }
    454         }
    455         else if (nlp.isProvenPrimalInfeasible()) {
    456             // prune it!
    457             lower_bound_ = 1e200;
    458             numNlpFailed_ = 0;
    459             printf("\
    460 BM_lp: At node %i : will fathom because of infeasibility\n",
    461                    current_index());
    462         }
    463         else if (nlp.isAbandoned()) {
    464             // nlp failed
    465             nlp.forceBranchable();
    466             lower_bound_ = nlp.getObjValue();
    467             CoinDisjointCopyN(nlp.getColSolution(), vars.size(),
    468                               primal_solution_);
    469             numNlpFailed_ += 1;
    470             // will branch, but save in the user data how many times we have
    471             // failed, and if we fail too many times in a row then just abandon
    472             // the node.
    473         }
    474         else {
    475             // complain loudly about a bug
    476             throw BCP_fatal_error("Impossible outcome by nlp.initialSolve()\n");
    477         }
    478     }
    479     else {
    480         /* TODO:
    481            don't test feasibility in every node
    482         */
    483         double integerTolerance;
    484         double cutOffIncrement;
    485         nlp.retrieve_options()->GetNumericValue("integer_tolerance",
    486                                                 integerTolerance,"bonmin.");
    487         /* FIXME: cutoff_decr was called cutoff_incr??? */
    488         nlp.retrieve_options()->GetNumericValue("cutoff_decr",
    489                                                 cutOffIncrement,"bonmin.");
    490         OsiSolverInterface * osi = getLpProblemPointer()->lp_solver;
    491    
    492         // Don't test feasibility if the node LP was infeasible
    493         if (osi->isProvenPrimalInfeasible() ) {   
    494             return sol;
    495         }
    496 
    497         if (!feasChecker_) {
    498             feasChecker_ = new IpCbcOACutGenerator2(&nlp, NULL, NULL,
    499                                                     cutOffIncrement,
    500                                                     integerTolerance, 0, 1);
    501             feasChecker_->setLocalSearchNodeLimit(0);
    502             feasChecker_->setMaxLocalSearch(0);
    503             feasChecker_->setMaxLocalSearchPerNode(0);
    504         }
    505 
    506         // The babSolver info used is the one containted in osi
    507         OsiBabSolver * babSolver =
    508             dynamic_cast<OsiBabSolver *> (osi->getAuxiliaryInfo());
    509         //assert(babSolver == &babSolver_);
    510         babSolver->setSolver(nlp);
    511         feasChecker_->generateCuts(*osi, cuts_);
    512         const int numvar = vars.size();
    513         double* solverSol = new double[numvar];
    514         double objValue = 1e200;
    515    
    516         if (babSolver->solution(objValue, solverSol,numvar)) {
    517             sol = new BM_solution;
    518             //Just copy the solution stored in solver to sol
    519             for (int i = 0 ; i < numvar ; i++) {
    520                 if (solverSol[i] > lp_result.primalTolerance())
    521                     sol->add_entry(i, solverSol[i]);
    522             }
    523             sol->setObjective(objValue);
    524         }
    525         delete[] solverSol;
    526     }
    527     return sol;
    528 }
    529 
    530 /****************************************************************************/
    531 
    532 void
    533 BM_lp::generate_cuts_in_lp(const BCP_lp_result& lpres,
    534                            const BCP_vec<BCP_var*>& vars,
    535                            const BCP_vec<BCP_cut*>& cuts,
    536                            BCP_vec<BCP_cut*>& new_cuts,
    537                            BCP_vec<BCP_row*>& new_rows)
    538 {
    539     if (! par.entry(BM_par::PureBranchAndBound)) {
    540         /* TODO:
    541            generate cuts with various other Cgl methods (Gomory, etc)
    542         */
    543         // fill cuts
    544         // eventually fill rows if not in strong branching
    545         int numCuts = cuts_.sizeRowCuts();
    546         for(int i = 0 ; i < numCuts ; i++)
    547             {
    548                 const OsiRowCut& cut = cuts_.rowCut(i);
    549                 new_cuts.push_back(new BB_cut(cut));
    550                 const CoinPackedVector& row = cut.row();
    551                 new_rows.push_back(new BCP_row(row, cut.lb(), cut.ub()));
    552             }
    553     }
    554     cuts_.dumpCuts();
    555 }
    556 
    557 /****************************************************************************/
    558 
    559 void
    560 BM_lp::cuts_to_rows(const BCP_vec<BCP_var*>& vars, // on what to expand
    561                     BCP_vec<BCP_cut*>& cuts,       // what to expand
    562                     BCP_vec<BCP_row*>& rows,       // the expanded rows
    563                     // things that the user can use for lifting cuts if allowed
    564                     const BCP_lp_result& lpres,
    565                     BCP_object_origin origin, bool allow_multiple)
    566 {
    567     const int numCuts = cuts.size();
    568     for (int i = 0; i < numCuts; ++i)
    569         {
    570             BB_cut* cut = dynamic_cast<BB_cut*>(cuts[i]);
    571             if (!cut) {
    572                 throw BCP_fatal_error("Non-\"BB_cut\" in cuts_to_rows!\n");
    573             }
    574             const CoinPackedVector& row = cut->row();
    575             rows.push_back(new BCP_row(row, cuts[i]->lb(), cuts[i]->ub()));
    576         }
    577 }
    578 
    579 /****************************************************************************/
    580 
    581 void
    582 BM_lp::vars_to_cols(const BCP_vec<BCP_cut*>& cuts,
    583                     BCP_vec<BCP_var*>& vars,
    584                     BCP_vec<BCP_col*>& cols,
    585                     const BCP_lp_result& lpres,
    586                     BCP_object_origin origin, bool allow_multiple)
    587 {
    588     /* All vars better be BM_branching_var, and their columns are empty */
    589     const int numVars = vars.size();
    590     for (int i = 0; i < numVars; ++i) {
    591         BM_branching_var* bv = dynamic_cast<BM_branching_var*>(vars[i]);
    592         if (!bv) {
    593             throw BCP_fatal_error("Non-\"BM_branching_var\" in vars_to_cols!\n");
    594         }
    595         cols.push_back(new BCP_col());
    596     }
    597 }
    598 
    599 /****************************************************************************/
    600 
    601 double
    602 BM_lp::compute_lower_bound(const double old_lower_bound,
    603                            const BCP_lp_result& lpres,
    604                            const BCP_vec<BCP_var*>& vars,
    605                            const BCP_vec<BCP_cut*>& cuts)
    606 {
    607     double bd;
    608     if (par.entry(BM_par::PureBranchAndBound)) {
    609         bd = lower_bound_;
    610     } else {
    611         bd = std::max(babSolver_.mipBound(), lpres.objval());
    612     }
    613     return std::max(bd, old_lower_bound);
    614 }
    615 
    616 /****************************************************************************/
    617 
    618 BCP_branching_decision
    619 BM_lp::select_branching_candidates(const BCP_lp_result& lpres,
    620                                    const BCP_vec<BCP_var*>& vars,
    621                                    const BCP_vec<BCP_cut*>& cuts,
    622                                    const BCP_lp_var_pool& local_var_pool,
    623                                    const BCP_lp_cut_pool& local_cut_pool,
    624                                    BCP_vec<BCP_lp_branching_object*>& cands)
    625 {
    626     Ipopt::SmartPtr<Ipopt::OptionsList> options = nlp.retrieve_options();
    627     double intTol;
    628     double cutoffIncr;
    629     options->GetNumericValue("integer_tolerance", intTol, "bonmin.");
    630     options->GetNumericValue("cutoff_decr", cutoffIncr, "bonmin.");
    631 
    632     if (lower_bound_ >= upper_bound() + cutoffIncr) {
    633         return BCP_DoNotBranch_Fathomed;
    634     }
    635 
    636     if (numNlpFailed_ >= par.entry(BM_par::NumNlpFailureMax)) {
    637         printf("WARNING! Too many (%i) NLP failures in a row. Abandoning node.",
    638                numNlpFailed_);
    639         return BCP_DoNotBranch_Fathomed;
    640     }
    641 
    642     if (par.entry(BM_par::BranchOnSos) && sos.size() > 0) {
    643         double bestval = 2.0;
    644         int bestsos = -1;
    645         int bestsplit = -1;
    646         int bestprio = 0;
    647         for (size_t ind = 0; ind < sos.size(); ++ind) {
    648             const BmSosInfo& sosi = *sos[ind];
    649             const int prio = sosi.priority;
    650             const int type = sosi.type;
    651             if ((type == 0 && par.entry(BM_par::BranchOnSos) == 2) ||
    652                 (type == 1 && par.entry(BM_par::BranchOnSos) == 1))
    653                 continue;
    654             if (bestsos >= 0 && prio != bestprio)
    655                 break;
    656             const int first = sosi.first;
    657             const int last = sosi.last;
    658             const int *indices = sosi.indices;
    659             // First count how many nonzeros are there
    660             int cnt = 0;
    661             for (int j = first; j < last; ++j) {
    662                 const double prim_i = primal_solution_[indices[j]];
    663                 if (prim_i > intTol && prim_i < 1-intTol) {
    664                     ++cnt;
    665                 }
    666             }
    667             // For type1 there can be at most 1, for type2 at most 2 nonzeros
    668             // Find the next sos if this is satisfied.
    669             if (cnt <= 1+type)
    670                 continue;
    671 
    672             // Find where does the sos goes over 0.5
    673             int half = -10;
    674             double val = 0.0;
    675             double primal = 0.0;
    676             for (half = first; half < last; ++half) {
    677                 val = primal_solution_[indices[half]];
    678                 primal += val;
    679                 if (primal > 0.5) {
    680                     break;
    681                 }
    682             }
    683 
    684             /* pick this over the prev best choice only if the value used to
    685                jump over the 0.5 boundary is smaller then the previous
    686                candidate */
    687             if (val >= bestval)
    688                 continue;
    689 
    690             /* The branches will be as follows:
    691                type1 : [first,half) and [half,last)
    692                type2 : [first,half) and (half,last)
    693                So if half==first then it must be incremented,
    694                and also for type 2 if half==last-1 then it must be decremented.
    695             */
    696             if (half == first)
    697                 ++half;
    698             if (type == 1 && half == last-1)
    699                 --half;
    700 
    701             bestval = val;
    702             bestsos = ind;
    703             bestsplit = half;
    704             bestprio = prio;
    705         }
    706 
    707         if (bestsos >= 0) {
    708             // OK, we can branch on an SOS
    709             // This vector will contain one entry only, but we must pass in a
    710             // vector
    711             BCP_vec<BCP_var*> new_vars;
    712             // This vector contains which vars have their ounds forcibly
    713             // changed. It's going to be the newly added branching var, hence
    714             // it's position is -1
    715             BCP_vec<int> fvp(1, -1);
    716             // This vector contains the new lb/ub pairs for all children. So
    717             // it's going to be {0.0, 0.0, 1.0, 1.0}
    718             BCP_vec<double> fvb(4, 0.0);
    719             fvb[2] = fvb[3] = 1.0;
    720             new_vars.push_back(new BM_branching_var(bestsos, bestsplit));
    721             cands.push_back(new BCP_lp_branching_object(2, // num of children
    722                                                         &new_vars,
    723                                                         NULL, // no new cuts
    724                                                         &fvp,NULL,&fvb,NULL,
    725                                                         NULL,NULL,NULL,NULL));
    726             if (par.entry(BM_par::PrintBranchingInfo)) {
    727                 printf("\
    728 BM_lp: At node %i branching on SOS%i set %i  split at: %i with val: %f\n",
    729                        current_index(), sos[bestsos]->type+1,
    730                        bestsos, bestsplit, bestval);
    731             }
    732             return BCP_DoBranch;;
    733         }
    734     }
    735 
    736     /* So we couldn't do an sos branching. If it's not pure B&B then we might
    737        as well do the built-in branching */
    738     if (! par.entry(BM_par::PureBranchAndBound)) {
    739         return BCP_lp_user::select_branching_candidates(lpres, vars, cuts,
    740                                                         local_var_pool,
    741                                                         local_cut_pool, cands);
    742     }
    743 
    744     /* If PURE_BB then we have the NLP optimum in "primal_solution_".
    745        Try to find a good branching candidate there. */
    746     const int numVars = vars.size();
    747     int besti = -1;
    748     double bestval = -1e200;
    749     for (int i = 0; i < numVars; ++i) {
    750         if (vars[i]->var_type() == BCP_ContinuousVar)
    751             continue;
    752         const double psol = CoinMin(CoinMax(vars[i]->lb(),primal_solution_[i]),
    753                                     vars[i]->ub());
    754         const double frac = psol - floor(psol);
    755         const double dist = std::min(frac, 1.0-frac);
    756        
    757         if (dist < intTol)
    758             continue;
    759         if (par.entry(BM_par::CombinedDistanceAndPriority)) {
    760             if (dist * branching_priority_[i] > bestval) {
    761                 bestval = dist * branching_priority_[i];
    762                 besti = i;
    763             }
    764         }
    765         else {
    766             if (branching_priority_[i] > bestval) {
    767                 bestval = branching_priority_[i];
    768                 besti = i;
    769             }
    770         }
    771     }
    772     if (besti == -1) {
    773         return BCP_DoNotBranch_Fathomed;
    774     }
    775     double psol = CoinMin(CoinMax(vars[besti]->lb(), primal_solution_[besti]),
    776                           vars[besti]->ub());
    777     if (par.entry(BM_par::PrintBranchingInfo)) {
    778         printf("BM_lp: branching on variable %i   value: %f\n", besti, psol);
    779     }
    780     BCP_vec<int> select_pos(1, besti);
    781     append_branching_vars(primal_solution_, vars, select_pos, cands);
    782     return BCP_DoBranch;
    783 }
    784 
    785 /****************************************************************************/
    786 
    787 void
    788 BM_lp::set_user_data_for_children(BCP_presolved_lp_brobj* best,
    789                                   const int selected)
    790 {
    791     BM_node* data = NULL;
    792     data = new BM_node;
    793     data->numNlpFailed_ = numNlpFailed_;
    794     best->user_data()[0] = data;
    795     data = new BM_node;
    796     data->numNlpFailed_ = numNlpFailed_;
    797     best->user_data()[1] = data;
    798 }
    799 
    800 //#############################################################################
    801 
    802 BM_lp::BmSosInfo::BmSosInfo(const TMINLP::SosInfo * sosinfo, int i)
    803 {
    804     priority = sosinfo->priorities ? sosinfo->priorities[i] : 0;
    805     length   = sosinfo->starts[i+1] - sosinfo->starts[i];
    806     type     = sosinfo->types[i] == 1 ? 0 : 1;
    807     indices  = new int[length];
    808     weights  = new double[length];
    809     first    = 0;
    810     last     = length;
    811     CoinDisjointCopyN(sosinfo->indices+sosinfo->starts[i], length, indices);
    812     CoinDisjointCopyN(sosinfo->weights+sosinfo->starts[i], length, weights);
    813     CoinSort_2(weights, weights+length, indices);
    814 }
    815    
    816 /*---------------------------------------------------------------------------*/
    817 
    818 void BM_lp::setSosFrom(const TMINLP::SosInfo * sosinfo)
    819 {
    820     if (!sosinfo || sosinfo->num == 0)
    821         return;
    822 
    823     //we have some sos constraints
    824     sos.reserve(sosinfo->num);
    825     for (int i = 0; i < sosinfo->num; ++i) {
    826         sos.push_back(new BmSosInfo(sosinfo, i));
    827     }
    828 
    829     if (sosinfo->priorities) {
    830         int * priorities = new int[sosinfo->num];
    831         CoinDisjointCopyN(sosinfo->priorities, sosinfo->num, priorities);
    832         /* FIXME: we may need to go down in order! */
    833         if (par.entry(BM_par::SosWithLowPriorityMoreImportant)) {
    834             CoinSort_2(priorities, priorities+sosinfo->num, &sos[0],
    835                        std::less<int>());
    836         } else {
    837             CoinSort_2(priorities, priorities+sosinfo->num, &sos[0],
    838                        std::greater<int>());
    839         }
    840         delete[] priorities;
    841     }
    842 }
    843 
    844 /*---------------------------------------------------------------------------*/
    845 
    846 #if 0
    847 void BM_lp::BmSosInfo::shuffle()
    848 {
    849     if (num == 0)
    850         return;
    851     int pos = 0;
    852     while (pos < num) {
    853         int cnt = 0;
    854         while (priorities[priority_order[pos]] ==
    855                priorities[priority_order[pos+cnt]]) cnt++;
    856         if (cnt > 1) {
    857             std::random_shuffle(priority_order+pos, priority_order+(pos+cnt));
    858         }
    859         pos += cnt;
    860     }
    861 }
    862 #endif
    863 
    864 //#############################################################################
    865 
    866117BCP_MemPool BM_branching_var::memPool(sizeof(BM_branching_var));
  • trunk/Bonmin/experimental/Bcp/BM.hpp

    r30 r42  
    230230    unpack_user_data(BCP_buffer& buf);
    231231
     232    virtual void
     233    process_message(BCP_buffer& buf);
     234
    232235    virtual OsiSolverInterface *
    233236    initialize_solver_interface();
  • trunk/Bonmin/experimental/Bcp/BM_pack.cpp

    r30 r42  
    7575    /* synchronize bonmin & BCP parameters */
    7676    Ipopt::SmartPtr<Ipopt::OptionsList> options = nlp.retrieve_options();
     77
     78    int nlpLogLevel;
     79    options->GetIntegerValue("nlp_log_level", nlpLogLevel, "bonmin.");
     80    nlp.messageHandler()->setLogLevel(nlpLogLevel);
     81
    7782    double bm_intTol;
    7883    double bm_cutoffIncr; // could be negative
  • trunk/Bonmin/experimental/Bcp/Makefile.am

    r21 r42  
    3737########################################################################
    3838
    39 bonminbcp_SOURCES = BB_cut.hpp BB_cut.cpp BM.hpp BM.cpp BM_pack.cpp bm_var.hpp
     39bonminbcp_SOURCES = \
     40        BB_cut.hpp \
     41        BB_cut.cpp \
     42        BM.hpp \
     43        BM.cpp \
     44        BM_tm.cpp \
     45        BM_lp.cpp \
     46        BM_lp_branch.cpp \
     47        BM_pack.cpp \
     48        BM_var.hpp
    4049
    4150bonminbcp_DEPENDENCIES = amplsolver.a
  • trunk/Bonmin/experimental/Bcp/Makefile.in

    r22 r42  
    5151
    5252subdir = experimental/Bcp
    53 DIST_COMMON = $(srcdir)/Makefile.am $(srcdir)/Makefile.in
     53DIST_COMMON = $(srcdir)/Makefile.am $(srcdir)/Makefile.in INSTALL
    5454ACLOCAL_M4 = $(top_srcdir)/aclocal.m4
    5555am__aclocal_m4_deps = $(top_srcdir)/acinclude.m4 \
     
    6363binPROGRAMS_INSTALL = $(INSTALL_PROGRAM)
    6464PROGRAMS = $(bin_PROGRAMS)
    65 am_bonminbcp_OBJECTS = BB_cut.$(OBJEXT) BM.$(OBJEXT) BM_pack.$(OBJEXT)
     65am_bonminbcp_OBJECTS = BB_cut.$(OBJEXT) BM.$(OBJEXT) BM_tm.$(OBJEXT) \
     66        BM_lp.$(OBJEXT) BM_lp_branch.$(OBJEXT) BM_pack.$(OBJEXT)
    6667bonminbcp_OBJECTS = $(am_bonminbcp_OBJECTS)
    6768am__DEPENDENCIES_1 = $(IPOPTOBJDIR)/src/Interfaces/libipopt.la
     
    321322#                                bonmin                                #
    322323########################################################################
    323 bonminbcp_SOURCES = BB_cut.hpp BB_cut.cpp BM.hpp BM.cpp BM_pack.cpp bm_var.hpp
     324bonminbcp_SOURCES = \
     325        BB_cut.hpp \
     326        BB_cut.cpp \
     327        BM.hpp \
     328        BM.cpp \
     329        BM_tm.cpp \
     330        BM_lp.cpp \
     331        BM_lp_branch.cpp \
     332        BM_pack.cpp \
     333        BM_var.hpp
     334
    324335bonminbcp_DEPENDENCIES = amplsolver.a
    325336AMPL_PATHED_FILES = \
     
    445456@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/BB_cut.Po@am__quote@
    446457@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/BM.Po@am__quote@
     458@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/BM_lp.Po@am__quote@
     459@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/BM_lp_branch.Po@am__quote@
    447460@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/BM_pack.Po@am__quote@
     461@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/BM_tm.Po@am__quote@
    448462
    449463.cpp.o:
Note: See TracChangeset for help on using the changeset viewer.