/home/coin/SVN-release/OS-2.4.1/Couenne/src/main/BonCouenne.cpp

Go to the documentation of this file.
00001 // $Id: BonCouenne.cpp 732 2011-07-03 20:06:50Z pbelotti $
00002 //
00003 // (C) Copyright International Business Machines Corporation and Carnegie Mellon University 2006, 2007
00004 // All Rights Reserved.
00005 // This code is published under the Eclipse Public License (EPL).
00006 //
00007 // Authors :
00008 // Pietro Belotti, Carnegie Mellon University,
00009 // Pierre Bonami, International Business Machines Corporation
00010 //
00011 // Date : 12/19/2006
00012 
00013 #if defined(_MSC_VER)
00014 // Turn off compiler warning about long names
00015 #  pragma warning(disable:4786)
00016 #endif
00017 
00018 #include <iomanip>
00019 #include <fstream>
00020 
00021 #include <stdlib.h>
00022 
00023 #include "CoinTime.hpp"
00024 #include "CoinError.hpp"
00025 #include "BonCouenneInterface.hpp"
00026 
00027 #include "BonCouenneSetup.hpp"
00028 
00029 #include "BonCbc.hpp"
00030 
00031 #include "CbcCutGenerator.hpp"
00032 #include "CouenneProblem.hpp"
00033 #include "CouenneCutGenerator.hpp"
00034 
00035 #include "CouenneRecordBestSol.hpp"
00036 
00037 using namespace Couenne;
00038 
00039 // the maximum difference between a printed optimum and a CouNumber
00040 #define PRINTED_PRECISION 1e-5
00041 
00042 #include "CouenneExprVar.hpp"
00043 #include "CouenneExprConst.hpp"
00044 #include "CouenneExprSum.hpp"
00045 #include "CouenneExprClone.hpp"
00046 #include "CouenneProblemElem.hpp"
00047 #include "CouenneProblem.hpp"
00048 #include "CouenneJournalist.hpp"
00049 
00050 #ifdef COIN_HAS_NTY
00051 int nOrbBr = 0; // FIXME: horrible global variable. Brrr.
00052 #endif
00053 
00054 #include "CoinSignal.hpp"
00055 
00056 #if 0
00057 extern "C" {
00058 
00059   static int nInterrupts = 0;
00060   static void signal_handler (int sig) {
00061 
00062     if (!nInterrupts) {
00063       std::cerr << "[BREAK]" << std::endl;
00064       abort ();
00065     }
00066     return;
00067   }
00068 }
00069 #endif
00070 
00071 //#define FM_FRES
00072 
00073 int main (int argc, char *argv[]) {
00074 
00075   //CoinSighandler_t saveSignal = SIG_DFL;
00076   //saveSignal = signal (SIGINT, signal_handler);
00077 
00078     printf ("Couenne %s --  an Open-Source solver for Mixed Integer Nonlinear Optimization\n\
00079 Mailing list: couenne@list.coin-or.org\n\
00080 Instructions: http://www.coin-or.org/Couenne\n", 
00081             strcmp (COUENNE_VERSION, "trunk") ? COUENNE_VERSION : "");
00082 
00083   WindowsErrorPopupBlocker();
00084   using namespace Ipopt;
00085 
00086   char * pbName = NULL;
00087 
00088   const int infeasible = 1;
00089 
00090   try {
00091 
00092     Bonmin::Bab bb;
00093     bb.setUsingCouenne (true);
00094 
00095     CouenneProblem *p = NULL;
00096     CouenneInterface *ci = NULL;
00097 
00098 #if 0
00099     //ci = new CouenneInterface;
00100     p = new CouenneProblem;
00101 
00102     p -> addVariable (false, p -> domain ());
00103     p -> addVariable (false, p -> domain ());
00104     p -> addVariable (false, p -> domain ());
00105     p -> addVariable (false, p -> domain ());
00106 
00107     p -> addObjective    (new exprSum (new exprClone (p->Var (1)), new exprClone (p->Var (2))), "min");
00108     p -> addLEConstraint (new exprSum (new exprClone (p->Var (0)), new exprClone (p->Var (2))), 
00109                           new exprConst (1));
00110     p -> addEQConstraint (new exprSum (new exprClone (p->Var (1)), new exprClone (p->Var (2))), 
00111                           new exprConst (1));
00112     p -> addEQConstraint (new exprSum (new exprClone (p->Var (1)), new exprClone (p->Var (3))), 
00113                           new exprConst (1));
00114     p -> addEQConstraint (new exprSum (new exprClone (p->Var (2)), new exprClone (p->Var (3))), 
00115                           new exprConst (1));
00116 #endif
00117 
00118     CouenneSetup couenne;
00119 
00120     if (!couenne.InitializeCouenne (argv, p, NULL, ci, &bb))
00121       throw infeasible;
00122 
00124     CouenneCutGenerator *cg = NULL;
00125 
00126     // there is only one CouenneCutGenerator object; scan array until
00127     // dynamic_cast returns non-NULL
00128 
00129     if (couenne. cutGenerators () . size () > 0) {
00130 
00131       for (std::list <Bonmin::BabSetupBase::CuttingMethod>::iterator 
00132                    i  = couenne.cutGenerators () . begin ();
00133            !cg && (i != couenne.cutGenerators () . end ()); 
00134            ++i) 
00135 
00136         cg = dynamic_cast <CouenneCutGenerator *> (i -> cgl);
00137     }
00138 
00139     // This assumes that first cut generator is CouenneCutGenerator
00140     // CouenneCutGenerator *cg = dynamic_cast<CouenneCutGenerator *> 
00141     //   (couenne.cutGenerators().begin()->cgl);
00142 
00143     if (cg)
00144 
00145       cg -> setBabPtr (&bb);
00146 
00147     else {
00148       printf ("main(): ### ERROR: Can not get CouenneCutGenerator\n");
00149       exit (1);
00150     }
00151 
00152     // initial printout
00153 
00154     ConstJnlstPtr jnlst = couenne. couennePtr () -> Jnlst ();
00155 
00156     CouenneProblem *prob = couenne. couennePtr () -> Problem ();
00157 
00158     jnlst -> Printf (J_ERROR, J_COUENNE, "\
00159 Loaded instance \"%s\"\n\
00160 Constraints:     %8d\n\
00161 Variables:       %8d (%d integer)\n\
00162 Auxiliaries:     %8d (%d integer)\n\n",
00163                      prob -> problemName ().c_str (),
00164                      prob -> nOrigCons (),
00165                      prob -> nOrigVars (),
00166                      prob -> nOrigIntVars (),
00167                      prob -> nVars    () - prob -> nOrigVars (),
00168                      CoinMax (0, prob -> nIntVars () - prob -> nOrigIntVars ()));
00169 
00170     double time_start = CoinCpuTime();
00171 
00172 #if 0
00173     CouenneFeasibility feasibility;
00174     bb.model().setProblemFeasibility (feasibility);
00175 #endif
00176 
00178     double timeLimit = 0;
00179     couenne.options () -> GetNumericValue ("time_limit", timeLimit, "couenne.");
00180     couenne.setDoubleParameter (Bonmin::BabSetupBase::MaxTime, 
00181                                 CoinMax (1., timeLimit - time_start));
00182 
00183     //jnlst -> Printf (J_ERROR, J_COUENNE, "Starting branch-and-bound\n");
00184 
00186 
00187     bb (couenne); // do branch and bound
00188 
00189 #ifdef COIN_HAS_NTY
00190     if (nOrbBr)
00191       printf ("%d orbital nontrivial branchings\n", nOrbBr);
00192 #endif
00193 
00195 
00196     std::cout.precision (10);
00197 
00199     int nr=-1, nt=-1;
00200     double st=-1;
00201 
00202     if (cg) cg -> getStats (nr, nt, st);
00203     else printf ("Warning, could not get pointer to CouenneCutGenerator\n");
00204 
00205     CouenneProblem *cp = cg ? cg -> Problem () : NULL;
00206 
00207 #if defined (FM_TRACE_OPTSOL) || defined (FM_FRES)
00208     double cbcLb = bb.model ().getBestPossibleObjValue();
00209     double printObj = 0;
00210     bool foundSol = false;
00211 #endif
00212 
00213 #ifdef FM_TRACE_OPTSOL
00214 
00215     FILE *fSol = fopen ("bidon.sol", "w");
00216 
00217     if(fSol == NULL) {
00218       printf("### ERROR: can not open bidon.sol\n");
00219       exit(1);
00220     }
00221 
00222     if(cp != NULL) {
00223       double cbcObjVal = bb.model().getObjValue();
00224       int modelNvars = bb.model().getNumCols();
00225 
00226       CouenneRecordBestSol *rs = cp->getRecordBestSol(); 
00227       const double *cbcSol = bb.model().getColSolution();
00228       double *modCbcSol = new double[modelNvars];
00229       double modCbcSolVal= 1e100, modCbcSolMaxViol = 0;
00230       bool cbcSolIsFeas = false;
00231 
00232       if(modelNvars != cp->nVars()) {
00233         printf("### ERROR: modelNvars: %d nVars: %d\n", 
00234                modelNvars, cp->nVars());
00235         exit(1);
00236       }
00237 
00238       if(cbcObjVal < 1e49) {
00239 
00240 #ifdef FM_CHECKNLP2
00241         int cMS = rs->getCardModSol();
00242         cbcSolIsFeas = cp->checkNLP2(cbcSol, 0, false, // do not care about obj
00243                                      false, // do not stop at first viol 
00244                                      true, // checkAll 
00245                                      cp->getFeasTol());
00246         CoinCopyN(rs->getModSol(cMS), cMS, modCbcSol);
00247         modCbcSolVal = rs->getModSolVal();
00248         modCbcSolMaxViol = rs->getModSolMaxViol();
00249 #else /* not FM_CHECKNLP2 */
00250         int cMS = cp->nVars();
00251         cbcSolIsFeas = cp->checkNLP(cbcSol, modCbcSolVal, true);
00252         CoinCopyN(cbcSol, cMS, modCbcSol);
00253         modCbcSolMaxViol = cp->getFeasTol();
00254 #endif /* not FM_CHECKNLP2 */
00255         foundSol = true;
00256       }
00257 
00258       const double *couenneSol = rs->getSol();
00259       double *modCouenneSol = new double[modelNvars];
00260       double modCouenneSolVal= 1e100, modCouenneSolMaxViol = 0;
00261       bool couenneSolIsFeas = false;
00262 
00263       if(couenneSol != NULL) {
00264 #ifdef FM_CHECKNLP2
00265         couenneSolIsFeas = cp->checkNLP2(couenneSol, 0, false, 
00266                                          false, true, 
00267                                          cp->getFeasTol());
00268         int cMS = rs->getCardModSol();
00269         CoinCopyN(rs->getModSol(cMS), cMS, modCouenneSol);
00270         modCouenneSolVal = rs->getModSolVal();
00271         modCouenneSolMaxViol = rs->getModSolMaxViol();
00272 #else /* not FM_CHECKNLP2 */
00273         couenneSolIsFeas = cp->checkNLP(couenneSol, modCouenneSolVal, true);
00274         int cMS = cp->nVars();
00275         CoinCopyN(couenneSol, cMS, modCouenneSol);
00276         modCouenneSolMaxViol = cp->getFeasTol();
00277 #endif /* not FM_CHECKNLP2 */
00278         foundSol = true;
00279       }
00280 
00281       int retcomp = rs->compareAndSave(modCbcSol, modCbcSolVal, 
00282                                        modCbcSolMaxViol, cbcSolIsFeas,
00283                                        modCouenneSol, modCouenneSolVal, 
00284                                        modCouenneSolMaxViol, couenneSolIsFeas, 
00285                                        modelNvars, cp->getFeasTol());
00286       switch (retcomp) {
00287       case -1: printf("No solution found\n"); break;
00288       case 0: printf("Best solution found by Cbc  Value: %10.4f  Tolerance: %10g\n", modCbcSolVal, modCbcSolMaxViol); break;
00289       case 1: printf("Best solution found by Couenne  Value: %10.4f  Tolerance: %10g\n", modCouenneSolVal, modCouenneSolMaxViol); break;
00290       default: break; // never happens
00291       }
00292 
00293       if(rs->getHasSol()) {
00294         if(cbcLb > rs->getVal()) { // Best sol found by Couenne and not
00295                                    // transmitted to Cbc
00296           cbcLb = rs->getVal();
00297         }
00298         printObj = rs->getVal();
00299         rs->printSol(fSol);
00300       }
00301       delete[] modCbcSol;
00302       delete[] modCouenneSol;
00303     }
00304     fclose(fSol);
00305 #endif /* FM_TRACE_OPTSOL */
00306 
00307 #ifdef FM_FRES
00308     if(cp != NULL) {
00309       FILE *f_res = NULL;
00310       f_res = fopen("fres.xxx", "r");
00311       if(f_res == NULL) {
00312         f_res = fopen("fres.xxx", "w");
00313         fprintf(f_res, "END_OF_HEADER\n");
00314       }
00315       else {
00316         fclose(f_res);
00317         f_res = fopen("fres.xxx", "a");
00318       }
00319       char *pbName, shortName[256];  
00320       
00321       pbName = strdup(cp -> problemName ().c_str ());
00322       char *f_name_pos = strrchr(pbName, '/');
00323       if(f_name_pos != NULL) {
00324         strcpy(shortName, &(f_name_pos[1]));
00325       }
00326       else {
00327         strcpy(shortName, pbName);
00328       }
00329       
00330       fprintf(f_res, "%20s %10.4f", shortName, cbcLb);
00331       if(foundSol) {
00332         fprintf(f_res, " %10.4f", printObj);
00333       }
00334       else {
00335         fprintf(f_res, "         *");
00336       }
00337       fprintf(f_res, " %10d %10.4f\n", bb.numNodes (),
00338               CoinCpuTime () - time_start);
00339       fclose(f_res);
00340     }
00341 #endif
00342 
00343     // retrieve test value to check
00344     double global_opt;
00345     couenne.options () -> GetNumericValue ("couenne_check", global_opt, "couenne.");
00346 
00347     double 
00348       ub = bb.model (). getObjValue (),
00349       lb = bb.model (). getBestPossibleObjValue ();
00350 
00351     char *gapstr = new char [80];
00352 
00353     sprintf (gapstr, "%.2f%%", 100. * (ub - lb) / (1. + fabs (lb)));
00354 
00355     jnlst -> Printf (J_ERROR, J_COUENNE, "\n\
00356 Linearization cuts added at root node:   %8d\n\
00357 Linearization cuts added in total:       %8d  (separation time: %gs)\n\
00358 Total solving time:                      %8gs (%gs in branch-and-bound)\n\
00359 Lower bound:                           %10g\n\
00360 Upper bound:                           %10g  (gap: %s)\n\
00361 Branch-and-bound nodes:                  %8d\n\n",
00362                      nr, nt, st, 
00363                      CoinCpuTime () - time_start,
00364                      cg ? (CoinCpuTime () - cg -> rootTime ()) : CoinCpuTime (),
00365                      lb, 
00366                      ub,
00367                      (ub > COUENNE_INFINITY/1e4) ? "inf" : gapstr,
00368                      bb.numNodes ());
00369 
00370     delete [] gapstr;
00371 
00372     if (global_opt < COUENNE_INFINITY) { // some value found in couenne.opt
00373 
00374       double opt = bb.model (). getBestPossibleObjValue ();
00375 
00376       printf ("Global Optimum Test on %-40s %s\n", 
00377               cp ? cp -> problemName ().c_str () : "unknown", 
00378               (fabs (opt - global_opt) / 
00379                (1. + CoinMax (fabs (opt), fabs (global_opt))) < PRINTED_PRECISION) ? 
00380               (const char *) "OK" : (const char *) "FAILED");
00381               //opt, global_opt,
00382               //fabs (opt - global_opt));
00383 
00384     } else // good old statistics
00385 
00386     if (couenne.displayStats ()) { // print statistics
00387 
00388       if (cg && !cp) printf ("Warning, could not get pointer to problem\n");
00389       else
00390         printf ("Stats: %-15s %4d [var] %4d [int] %4d [con] %4d [aux] "
00391                 "%6d [root] %8d [tot] %6g [sep] %8g [time] %8g [bb] "
00392                 "%20e [lower] %20e [upper] %7d [nodes]\n",// %s %s\n",
00393                 cp ? cp -> problemName (). c_str () : "unknown",
00394                 (cp) ? cp -> nOrigVars     () : -1, 
00395                 (cp) ? cp -> nOrigIntVars  () : -1, 
00396                 (cp) ? cp -> nOrigCons     () : -1,
00397                 (cp) ? (cp -> nVars     () - 
00398                         cp -> nOrigVars ()): -1,
00399                 nr, nt, st, 
00400                 CoinCpuTime () - time_start,
00401                 cg ? (CoinCpuTime () - cg -> rootTime ()) : CoinCpuTime (),
00402                 bb.model (). getBestPossibleObjValue (),
00403                 bb.model (). getObjValue (),
00404                 //bb.bestBound (),
00405                 //bb.bestObj (),
00406                 bb.numNodes ());
00407                 //bb.iterationCount ());
00408                 //status.c_str (), message.c_str ());
00409     }
00410 
00411 //    nlp_and_solver -> writeAmplSolFile (message, bb.bestSolution (), NULL);
00412   }
00413   catch(Bonmin::TNLPSolver::UnsolvedError *E) {
00414      E->writeDiffFiles();
00415      E->printError(std::cerr);
00416     //There has been a failure to solve a problem with Ipopt.
00417     //And we will output file with information on what has been changed in the problem to make it fail.
00418     //Now depending on what algorithm has been called (B-BB or other) the failed problem may be at different place.
00419     //    const OsiSolverInterface &si1 = (algo > 0) ? nlpSolver : *model.solver();
00420   }
00421   catch (Bonmin::OsiTMINLPInterface::SimpleError &E) {
00422     std::cerr<<E.className()<<"::"<<E.methodName()
00423              <<std::endl
00424              <<E.message()<<std::endl;
00425   }
00426   catch (CoinError &E) {
00427     std::cerr<<E.className()<<"::"<<E.methodName()
00428              <<std::endl
00429              <<E.message()<<std::endl;
00430   }
00431   catch (Ipopt::OPTION_INVALID &E)
00432   {
00433    std::cerr<<"Ipopt exception : "<<E.Message()<<std::endl;
00434   }
00435   catch (int generic_error) {
00436     if (generic_error == infeasible)
00437       printf ("problem infeasible\n");
00438   }
00439 
00440 //  catch(...) {
00441 //    std::cerr<<pbName<<" unrecognized excpetion"<<std::endl;
00442 //    std::cerr<<pbName<<"\t Finished \t exception"<<std::endl;
00443 //    throw;
00444 //  }
00445 
00446   delete [] pbName;
00447   return 0;
00448 }

Generated on Thu Nov 10 03:05:45 2011 by  doxygen 1.4.7