source: trunk/Couenne/src/main/BonCouenne.cpp @ 577

Last change on this file since 577 was 577, checked in by pbelotti, 9 years ago

restoring FP for tests. Changing branching point on integer variables with integer value. Fixing large branching points for variables with large (upper/lower) bounds and large LP values. Added debug code to twoImplBounds for check against optimal solutions. Added fictitious objective function of 0 when none is defined. Redeclared size_t in invmap.cpp for compatibility with MSVS.

  • Property svn:keywords set to Author Date Id Revision
File size: 11.8 KB
Line 
1// $Id: BonCouenne.cpp 577 2011-05-21 20:38:48Z pbelotti $
2//
3// (C) Copyright International Business Machines Corporation and Carnegie Mellon University 2006, 2007
4// All Rights Reserved.
5// This code is published under the Eclipse Public License (EPL).
6//
7// Authors :
8// Pietro Belotti, Carnegie Mellon University,
9// Pierre Bonami, International Business Machines Corporation
10//
11// Date : 12/19/2006
12
13
14#if defined(_MSC_VER)
15// Turn off compiler warning about long names
16#  pragma warning(disable:4786)
17#endif
18#include <iomanip>
19#include <fstream>
20
21#include <stdlib.h>
22
23#include "CoinTime.hpp"
24#include "CoinError.hpp"
25#include "BonCouenneInterface.hpp"
26
27#include "BonCouenneSetup.hpp"
28
29#include "BonCbc.hpp"
30
31#include "CbcCutGenerator.hpp"
32#include "CouenneProblem.hpp"
33#include "CouenneCutGenerator.hpp"
34
35#include "CouenneRecordBestSol.hpp"
36
37using namespace Couenne;
38
39// the maximum difference between a printed optimum and a CouNumber
40#define PRINTED_PRECISION 1e-5
41
42#include "CouenneExprVar.hpp"
43#include "CouenneExprConst.hpp"
44#include "CouenneExprSum.hpp"
45#include "CouenneExprClone.hpp"
46#include "CouenneProblemElem.hpp"
47#include "CouenneProblem.hpp"
48#include "CouenneJournalist.hpp"
49
50#ifdef COIN_HAS_NTY
51int nOrbBr = 0;
52#endif
53
54#include "CoinSignal.hpp"
55
56extern "C" {
57
58  static int nInterrupts = 0;
59  static void signal_handler (int sig) {
60
61    if (!nInterrupts) {
62      std::cerr << "[BREAK]" << std::endl;
63      abort ();
64    }
65    return;
66  }
67}
68
69int main (int argc, char *argv[]) {
70
71  //CoinSighandler_t saveSignal = SIG_DFL;
72  //saveSignal = signal (SIGINT, signal_handler);
73
74    printf ("Couenne %s --  an Open-Source solver for Mixed Integer Nonlinear Optimization\n\
75Mailing list: couenne@list.coin-or.org\n\
76Instructions: http://www.coin-or.org/Couenne\n", 
77            strcmp (COUENNE_VERSION, "trunk") ? COUENNE_VERSION : "");
78
79  WindowsErrorPopupBlocker();
80  using namespace Ipopt;
81
82  char * pbName = NULL;
83
84  const int infeasible = 1;
85
86  try {
87
88    Bonmin::Bab bb;
89    bb.setUsingCouenne (true);
90
91    CouenneProblem *p = NULL;
92    CouenneInterface *ci = NULL;
93
94#if 0
95    //ci = new CouenneInterface;
96    p = new CouenneProblem;
97
98    p -> addVariable (false, p -> domain ());
99    p -> addVariable (false, p -> domain ());
100    p -> addVariable (false, p -> domain ());
101    p -> addVariable (false, p -> domain ());
102
103    p -> addObjective    (new exprSum (new exprClone (p->Var (1)), new exprClone (p->Var (2))), "min");
104    p -> addLEConstraint (new exprSum (new exprClone (p->Var (0)), new exprClone (p->Var (2))),
105                          new exprConst (1));
106    p -> addEQConstraint (new exprSum (new exprClone (p->Var (1)), new exprClone (p->Var (2))),
107                          new exprConst (1));
108    p -> addEQConstraint (new exprSum (new exprClone (p->Var (1)), new exprClone (p->Var (3))),
109                          new exprConst (1));
110    p -> addEQConstraint (new exprSum (new exprClone (p->Var (2)), new exprClone (p->Var (3))),
111                          new exprConst (1));
112#endif
113
114    CouenneSetup couenne;
115    if (!couenne.InitializeCouenne (argv, p, NULL, ci, &bb))
116      throw infeasible;
117
118#ifdef FM_PRINT_INFO
119    // This assumes that first cut generator is CouenneCutGenerator
120    CouenneCutGenerator *ccg = dynamic_cast<CouenneCutGenerator *> 
121      (couenne.cutGenerators().begin()->cgl);
122    if(ccg) {
123      ccg->setBabPtr(&bb);
124    }
125    else {
126      printf("main(): ### ERROR: Can not get CouenneCutGenerator\n");
127      exit(1);
128    }
129#endif
130
131    // initial printout
132
133    ConstJnlstPtr jnlst = couenne. couennePtr () -> Jnlst ();
134
135    CouenneProblem *prob = couenne. couennePtr () -> Problem ();
136
137    jnlst -> Printf (J_ERROR, J_COUENNE, "\
138Loaded instance \"%s\"\n\
139Constraints:     %8d\n\
140Variables:       %8d (%d integer)\n\
141Auxiliaries:     %8d (%d integer)\n\n",
142                     prob -> problemName ().c_str (),
143                     prob -> nOrigCons (),
144                     prob -> nOrigVars (),
145                     prob -> nOrigIntVars (),
146                     prob -> nVars    () - prob -> nOrigVars (),
147                     prob -> nIntVars () - prob -> nOrigIntVars ());
148
149    double time_start = CoinCpuTime();
150
151#if 0
152    CouenneFeasibility feasibility;
153    bb.model().setProblemFeasibility (feasibility);
154#endif
155
156    /// update time limit (read/preprocessing might have taken some)
157    double timeLimit = 0;
158    couenne.options () -> GetNumericValue ("time_limit", timeLimit, "couenne.");
159    couenne.setDoubleParameter (Bonmin::BabSetupBase::MaxTime, 
160                                CoinMax (1., timeLimit - time_start));
161
162    jnlst -> Printf (J_ERROR, J_COUENNE, "Starting branch-and-bound\n");
163
164    //////////////////////////////////
165
166    bb (couenne); // do branch and bound
167
168#ifdef COIN_HAS_NTY
169    if (nOrbBr)
170      printf ("%d orbital nontrivial branchings\n", nOrbBr);
171#endif
172
173    //////////////////////////////////
174
175    std::cout.precision (10);
176
177    //////////////////////////////
178    CouenneCutGenerator *cg = NULL;
179
180    // there is only one cut generator, so scan array until
181    // dynamic_cast returns non-NULL
182
183    if (bb.model (). cutGenerators ()) {
184
185      int nGen = bb.model (). numberCutGenerators ();
186     
187      for (int i=0; !cg && i < nGen; i++)
188        cg = dynamic_cast <CouenneCutGenerator *> 
189          (bb.model (). cutGenerators () [i] -> generator ());
190    }
191
192    ////////////////////////////////
193    int nr=-1, nt=-1;
194    double st=-1;
195
196    if (cg) cg -> getStats (nr, nt, st);
197    else printf ("Warning, could not get pointer to CouenneCutGenerator\n");
198
199    CouenneProblem *cp = cg ? cg -> Problem () : NULL;
200
201    double cbcLb = bb.model ().getBestPossibleObjValue();
202    double printObj = 0;
203    bool foundSol = false;
204
205#ifdef FM_TRACE_OPTSOL
206
207    FILE *fSol = fopen("bidon.sol", "w");
208    if(fSol == NULL) {
209      printf("### ERROR: can not open bidon.sol\n");
210      exit(1);
211    }
212
213    if(cp != NULL) {
214      double cbcObjVal = bb.model().getObjValue();
215      int modelNvars = bb.model().getNumCols();
216
217      CouenneRecordBestSol *rs = cp->getRecordBestSol(); 
218      const double *cbcSol = bb.model().getColSolution();
219      double *modCbcSol = new double[modelNvars];
220      double modCbcSolVal= 1e100, modCbcSolMaxViol = 0;
221      bool cbcSolIsFeas = false;
222
223      if(modelNvars != cp->nVars()) {
224        printf("### ERROR: modelNvars: %d nVars: %d\n", 
225               modelNvars, cp->nVars());
226        exit(1);
227      }
228
229      if(cbcObjVal < 1e49) {
230
231#ifdef FM_CHECKNLP2
232        int cMS = rs->getCardModSol();
233        cbcSolIsFeas = cp->checkNLP2(cbcSol, 0, false, // do not care about obj
234                                     false, // do not stop at first viol
235                                     true, // checkAll
236                                     cp->getFeasTol());
237        CoinCopyN(rs->getModSol(cMS), cMS, modCbcSol);
238        modCbcSolVal = rs->getModSolVal();
239        modCbcSolMaxViol = rs->getModSolMaxViol();
240#else /* not FM_CHECKNLP2 */
241        int cMS = cp->nVars();
242        cbcSolIsFeas = cp->checkNLP(cbcSol, modCbcSolVal, true);
243        CoinCopyN(cbcSol, cMS, modCbcSol);
244        modCbcSolMaxViol = cp->getFeasTol();
245#endif /* not FM_CHECKNLP2 */
246        foundSol = true;
247      }
248
249      const double *couenneSol = rs->getSol();
250      double *modCouenneSol = new double[modelNvars];
251      double modCouenneSolVal= 1e100, modCouenneSolMaxViol = 0;
252      bool couenneSolIsFeas = false;
253
254      if(couenneSol != NULL) {
255#ifdef FM_CHECKNLP2
256        couenneSolIsFeas = cp->checkNLP2(couenneSol, 0, false, 
257                                         false, true, 
258                                         cp->getFeasTol());
259        int cMS = rs->getCardModSol();
260        CoinCopyN(rs->getModSol(cMS), cMS, modCouenneSol);
261        modCouenneSolVal = rs->getModSolVal();
262        modCouenneSolMaxViol = rs->getModSolMaxViol();
263#else /* not FM_CHECKNLP2 */
264        couenneSolIsFeas = cp->checkNLP(couenneSol, modCouenneSolVal, true);
265        int cMS = cp->nVars();
266        CoinCopyN(couenneSol, cMS, modCouenneSol);
267        modCouenneSolMaxViol = cp->getFeasTol();
268#endif /* not FM_CHECKNLP2 */
269        foundSol = true;
270      }
271
272      int retcomp = rs->compareAndSave(modCbcSol, modCbcSolVal, 
273                                       modCbcSolMaxViol, cbcSolIsFeas,
274                                       modCouenneSol, modCouenneSolVal, 
275                                       modCouenneSolMaxViol, couenneSolIsFeas, 
276                                       modelNvars, cp->getFeasTol());
277      switch (retcomp) {
278      case -1: printf("No solution found\n"); break;
279      case 0: printf("Best solution found by Cbc\n"); break;
280      case 1: printf("Best solution found by Couenne\n"); break;
281      default: break; // never happens
282      }
283
284      if(rs->getHasSol()) {
285        if(cbcLb > rs->getVal()) { // Best sol found by Couenne and not
286                                   // transmitted to Cbc
287          cbcLb = rs->getVal();
288        }
289        printObj = rs->getVal();
290        rs->printSol(fSol);
291      }
292      delete[] modCbcSol;
293      delete[] modCouenneSol;
294    }
295    fclose(fSol);
296#endif /* FM_TRACE_OPTSOL */
297
298
299    // retrieve test value to check
300    double global_opt;
301    couenne.options () -> GetNumericValue ("couenne_check", global_opt, "couenne.");
302
303    double 
304      ub = bb.model (). getObjValue (),
305      lb = bb.model (). getBestPossibleObjValue ();
306
307    char *gapstr = new char [80];
308
309    sprintf (gapstr, "%.2f%%", 100. * (ub - lb) / (1. + fabs (lb)));
310
311    jnlst -> Printf (J_ERROR, J_COUENNE, "\n\
312Linearization cuts added at root node:   %8d\n\
313Linearization cuts added in total:       %8d  (separation time: %gs)\n\
314Total solving time:                      %8gs (%gs in branch-and-bound)\n\
315Lower bound:                           %10g\n\
316Upper bound:                           %10g  (gap: %s)\n\
317Branch-and-bound nodes:                  %8d\n\n",
318                     nr, nt, st, 
319                     CoinCpuTime () - time_start,
320                     cg ? (CoinCpuTime () - cg -> rootTime ()) : CoinCpuTime (),
321                     lb, 
322                     ub,
323                     (ub > COUENNE_INFINITY/1e4) ? "inf" : gapstr,
324                     bb.numNodes ());
325
326    delete [] gapstr;
327
328    if (global_opt < COUENNE_INFINITY) { // some value found in couenne.opt
329
330      double opt = bb.model (). getBestPossibleObjValue ();
331
332      printf ("Global Optimum Test on %-40s %s\n", 
333              cp ? cp -> problemName ().c_str () : "unknown", 
334              (fabs (opt - global_opt) / 
335               (1. + CoinMax (fabs (opt), fabs (global_opt))) < PRINTED_PRECISION) ? 
336              (const char *) "OK" : (const char *) "FAILED");
337              //opt, global_opt,
338              //fabs (opt - global_opt));
339
340    } else // good old statistics
341
342    if (couenne.displayStats ()) { // print statistics
343
344      // CAUTION: assuming first cut generator is our CouenneCutGenerator
345
346      if (cg && !cp) printf ("Warning, could not get pointer to problem\n");
347      else
348        printf ("Stats: %-15s %4d [var] %4d [int] %4d [con] %4d [aux] "
349                "%6d [root] %8d [tot] %6g [sep] %8g [time] %8g [bb] "
350                "%20e [lower] %20e [upper] %7d [nodes]\n",// %s %s\n",
351                cp ? cp -> problemName (). c_str () : "unknown",
352                (cp) ? cp -> nOrigVars     () : -1, 
353                (cp) ? cp -> nOrigIntVars  () : -1, 
354                (cp) ? cp -> nOrigCons     () : -1,
355                (cp) ? (cp -> nVars     () - 
356                        cp -> nOrigVars ()): -1,
357                nr, nt, st, 
358                CoinCpuTime () - time_start,
359                cg ? (CoinCpuTime () - cg -> rootTime ()) : CoinCpuTime (),
360                bb.model (). getBestPossibleObjValue (),
361                bb.model (). getObjValue (),
362                //bb.bestBound (),
363                //bb.bestObj (),
364                bb.numNodes ());
365                //bb.iterationCount ());
366                //status.c_str (), message.c_str ());
367    }
368
369//    nlp_and_solver -> writeAmplSolFile (message, bb.bestSolution (), NULL);
370  }
371  catch(Bonmin::TNLPSolver::UnsolvedError *E) {
372     E->writeDiffFiles();
373     E->printError(std::cerr);
374    //There has been a failure to solve a problem with Ipopt.
375    //And we will output file with information on what has been changed in the problem to make it fail.
376    //Now depending on what algorithm has been called (B-BB or other) the failed problem may be at different place.
377    //    const OsiSolverInterface &si1 = (algo > 0) ? nlpSolver : *model.solver();
378  }
379  catch (Bonmin::OsiTMINLPInterface::SimpleError &E) {
380    std::cerr<<E.className()<<"::"<<E.methodName()
381             <<std::endl
382             <<E.message()<<std::endl;
383  }
384  catch (CoinError &E) {
385    std::cerr<<E.className()<<"::"<<E.methodName()
386             <<std::endl
387             <<E.message()<<std::endl;
388  }
389  catch (Ipopt::OPTION_INVALID &E)
390  {
391   std::cerr<<"Ipopt exception : "<<E.Message()<<std::endl;
392  }
393  catch (int generic_error) {
394    if (generic_error == infeasible)
395      printf ("problem infeasible\n");
396  }
397
398//  catch(...) {
399//    std::cerr<<pbName<<" unrecognized excpetion"<<std::endl;
400//    std::cerr<<pbName<<"\t Finished \t exception"<<std::endl;
401//    throw;
402//  }
403
404  delete [] pbName;
405  return 0;
406}
Note: See TracBrowser for help on using the repository browser.