source: trunk/Cbc/src/CbcCutGenerator.cpp

Last change on this file was 2705, checked in by stefan, 9 days ago

disable use of pthread_getcpuclockid on Mac OS X

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 48.7 KB
Line 
1/* $Id: CbcCutGenerator.cpp 2705 2019-10-08 08:59:48Z stefan $ */
2// Copyright (C) 2003, International Business Machines
3// Corporation and others.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6#if defined(_MSC_VER)
7// Turn off compiler warning about long names
8#pragma warning(disable : 4786)
9#endif
10#include "CbcConfig.h"
11#include <cassert>
12#include <cstdlib>
13#include <cmath>
14#include <cfloat>
15
16#ifdef COIN_HAS_CLP
17#include "OsiClpSolverInterface.hpp"
18#else
19#include "OsiSolverInterface.hpp"
20#endif
21//#define CGL_DEBUG 1
22#ifdef CGL_DEBUG
23#include "OsiRowCutDebugger.hpp"
24#endif
25#include "CbcModel.hpp"
26#include "CbcMessage.hpp"
27#include "CbcCutGenerator.hpp"
28#include "CbcBranchDynamic.hpp"
29#include "CglProbing.hpp"
30#include "CoinTime.hpp"
31#ifdef CBC_THREAD
32// need time on a thread by thread basis
33#include <pthread.h>
34#include <time.h>
35#endif
36
37// Default Constructor
38CbcCutGenerator::CbcCutGenerator()
39  : timeInCutGenerator_(0.0)
40  , model_(NULL)
41  , generator_(NULL)
42  , generatorName_(NULL)
43  , whenCutGenerator_(-1)
44  , whenCutGeneratorInSub_(-100)
45  , switchOffIfLessThan_(0)
46  , depthCutGenerator_(-1)
47  , depthCutGeneratorInSub_(-1)
48  , inaccuracy_(0)
49  , numberTimes_(0)
50  , numberCuts_(0)
51  , numberElements_(0)
52  , numberColumnCuts_(0)
53  , numberCutsActive_(0)
54  , numberCutsAtRoot_(0)
55  , numberActiveCutsAtRoot_(0)
56  , numberShortCutsAtRoot_(0)
57  , switches_(1)
58  , maximumTries_(-1)
59{
60}
61// Normal constructor
62CbcCutGenerator::CbcCutGenerator(CbcModel *model, CglCutGenerator *generator,
63  int howOften, const char *name,
64  bool normal, bool atSolution,
65  bool infeasible, int howOftenInSub,
66  int whatDepth, int whatDepthInSub,
67  int switchOffIfLessThan)
68  : timeInCutGenerator_(0.0)
69  , depthCutGenerator_(whatDepth)
70  , depthCutGeneratorInSub_(whatDepthInSub)
71  , inaccuracy_(0)
72  , numberTimes_(0)
73  , numberCuts_(0)
74  , numberElements_(0)
75  , numberColumnCuts_(0)
76  , numberCutsActive_(0)
77  , numberCutsAtRoot_(0)
78  , numberActiveCutsAtRoot_(0)
79  , numberShortCutsAtRoot_(0)
80  , switches_(1)
81  , maximumTries_(-1)
82{
83  if (howOften < -1900) {
84    setGlobalCuts(true);
85    howOften += 2000;
86  } else if (howOften < -900) {
87    setGlobalCutsAtRoot(true);
88    howOften += 1000;
89  }
90  model_ = model;
91  generator_ = generator->clone();
92  generator_->refreshSolver(model_->solver());
93  setNeedsOptimalBasis(generator_->needsOptimalBasis());
94  whenCutGenerator_ = howOften;
95  whenCutGeneratorInSub_ = howOftenInSub;
96  switchOffIfLessThan_ = switchOffIfLessThan;
97  if (name)
98    generatorName_ = CoinStrdup(name);
99  else
100    generatorName_ = CoinStrdup("Unknown");
101  setNormal(normal);
102  setAtSolution(atSolution);
103  setWhenInfeasible(infeasible);
104}
105
106// Copy constructor
107CbcCutGenerator::CbcCutGenerator(const CbcCutGenerator &rhs)
108{
109  model_ = rhs.model_;
110  generator_ = rhs.generator_->clone();
111  //generator_->refreshSolver(model_->solver());
112  whenCutGenerator_ = rhs.whenCutGenerator_;
113  whenCutGeneratorInSub_ = rhs.whenCutGeneratorInSub_;
114  switchOffIfLessThan_ = rhs.switchOffIfLessThan_;
115  depthCutGenerator_ = rhs.depthCutGenerator_;
116  depthCutGeneratorInSub_ = rhs.depthCutGeneratorInSub_;
117  generatorName_ = CoinStrdup(rhs.generatorName_);
118  switches_ = rhs.switches_;
119  maximumTries_ = rhs.maximumTries_;
120  timeInCutGenerator_ = rhs.timeInCutGenerator_;
121  savedCuts_ = rhs.savedCuts_;
122  inaccuracy_ = rhs.inaccuracy_;
123  numberTimes_ = rhs.numberTimes_;
124  numberCuts_ = rhs.numberCuts_;
125  numberElements_ = rhs.numberElements_;
126  numberColumnCuts_ = rhs.numberColumnCuts_;
127  numberCutsActive_ = rhs.numberCutsActive_;
128  numberCutsAtRoot_ = rhs.numberCutsAtRoot_;
129  numberActiveCutsAtRoot_ = rhs.numberActiveCutsAtRoot_;
130  numberShortCutsAtRoot_ = rhs.numberShortCutsAtRoot_;
131}
132
133// Assignment operator
134CbcCutGenerator &
135CbcCutGenerator::operator=(const CbcCutGenerator &rhs)
136{
137  if (this != &rhs) {
138    delete generator_;
139    free(generatorName_);
140    model_ = rhs.model_;
141    generator_ = rhs.generator_->clone();
142    generator_->refreshSolver(model_->solver());
143    whenCutGenerator_ = rhs.whenCutGenerator_;
144    whenCutGeneratorInSub_ = rhs.whenCutGeneratorInSub_;
145    switchOffIfLessThan_ = rhs.switchOffIfLessThan_;
146    depthCutGenerator_ = rhs.depthCutGenerator_;
147    depthCutGeneratorInSub_ = rhs.depthCutGeneratorInSub_;
148    generatorName_ = CoinStrdup(rhs.generatorName_);
149    switches_ = rhs.switches_;
150    maximumTries_ = rhs.maximumTries_;
151    timeInCutGenerator_ = rhs.timeInCutGenerator_;
152    savedCuts_ = rhs.savedCuts_;
153    inaccuracy_ = rhs.inaccuracy_;
154    numberTimes_ = rhs.numberTimes_;
155    numberCuts_ = rhs.numberCuts_;
156    numberElements_ = rhs.numberElements_;
157    numberColumnCuts_ = rhs.numberColumnCuts_;
158    numberCutsActive_ = rhs.numberCutsActive_;
159    numberCutsAtRoot_ = rhs.numberCutsAtRoot_;
160    numberActiveCutsAtRoot_ = rhs.numberActiveCutsAtRoot_;
161    numberShortCutsAtRoot_ = rhs.numberShortCutsAtRoot_;
162  }
163  return *this;
164}
165
166// Destructor
167CbcCutGenerator::~CbcCutGenerator()
168{
169  free(generatorName_);
170  delete generator_;
171}
172
173/* This is used to refresh any inforamtion.
174   It also refreshes the solver in the cut generator
175   in case generator wants to do some work
176*/
177void CbcCutGenerator::refreshModel(CbcModel *model)
178{
179  model_ = model;
180  // added test - helps if generator not thread safe
181  if (whenCutGenerator_ != -100)
182    generator_->refreshSolver(model_->solver());
183}
184/* Generate cuts for the model data contained in si.
185   The generated cuts are inserted into and returned in the
186   collection of cuts cs.
187*/
188bool CbcCutGenerator::generateCuts(OsiCuts &cs, int fullScan, OsiSolverInterface *solver, CbcNode *node)
189{
190  /*
191          Make some decisions about whether we'll generate cuts. First convert
192          whenCutGenerator_ to a set of canonical values for comparison to the node
193          count.
194
195                 0 <    mod 1000000, with a result of 0 forced to 1
196           -99 <= <= 0  convert to 1
197          -100 =        Off, period
198        */
199  int depth;
200  if (node)
201    depth = node->depth();
202  else
203    depth = 0;
204  int howOften = whenCutGenerator_;
205  if (dynamic_cast< CglProbing * >(generator_)) {
206    if (howOften == -100 && model_->doCutsNow(3)) {
207      howOften = 1; // do anyway
208    }
209  }
210  if (howOften == -100)
211    return false;
212  int pass = model_->getCurrentPassNumber() - 1;
213  if (maximumTries_ > 0) {
214    // howOften means what it says
215    if ((pass % howOften) != 0 || depth)
216      return false;
217    else
218      howOften = 1;
219  }
220  if (howOften > 0)
221    howOften = howOften % 1000000;
222  else
223    howOften = 1;
224  if (!howOften)
225    howOften = 1;
226  bool returnCode = false;
227  //OsiSolverInterface * solver = model_->solver();
228  // Reset cuts on first pass
229  if (!pass)
230    savedCuts_ = OsiCuts();
231  /*
232          Determine if we should generate cuts based on node count.
233        */
234  bool doThis = (model_->getNodeCount() % howOften) == 0;
235  /*
236          If the user has provided a depth specification, it will override the node
237          count specification.
238        */
239  if (depthCutGenerator_ > 0) {
240    doThis = (depth % depthCutGenerator_) == 0;
241    if (depth < depthCutGenerator_)
242      doThis = true; // and also at top of tree
243  }
244  /*
245          A few magic numbers ...
246
247          The distinction between -100 and 100 for howOften is that we can override 100
248          with fullScan. -100 means no cuts, period. As does the magic number -200 for
249          whenCutGeneratorInSub_.
250        */
251
252  // But turn off if 100
253  if (howOften == 100)
254    doThis = false;
255  // Switch off if special setting
256  if (whenCutGeneratorInSub_ == -200 && model_->parentModel()) {
257    fullScan = 0;
258    doThis = false;
259  }
260  if (fullScan || doThis) {
261    CoinThreadRandom *randomNumberGenerator = NULL;
262#ifdef COIN_HAS_CLP
263    {
264      OsiClpSolverInterface *clpSolver
265        = dynamic_cast< OsiClpSolverInterface * >(solver);
266      if (clpSolver)
267        randomNumberGenerator = clpSolver->getModelPtr()->randomNumberGenerator();
268    }
269#endif
270    double time1 = 0.0;
271    //#undef CBC_THREAD
272    /* TODO there should be check in configure or the Posix version C preprocessor variable
273     * to decide whether pthread_getcpuclockid is available
274     */
275#if defined(_MSC_VER) || defined(__MSVCRT__) || defined(__APPLE__) || !defined(CBC_THREAD)
276    if (timing()) 
277      time1 = CoinCpuTime();
278#else
279    struct timespec currTime;
280    clockid_t threadClockId;
281    if (timing()) {
282      if (!model_->getNumberThreads()) {
283        time1 = CoinCpuTime();
284      } else {
285        // Get thread clock Id
286        pthread_getcpuclockid(pthread_self(), &threadClockId);
287        // Using thread clock Id get the clock time
288        clock_gettime(threadClockId, &currTime);
289        time1 = static_cast<double>(currTime.tv_sec)
290          +1.0e-9*static_cast<double>(currTime.tv_nsec);
291      }
292    }
293#endif
294    //#define CBC_DEBUG
295    int numberRowCutsBefore = cs.sizeRowCuts();
296    int numberColumnCutsBefore = cs.sizeColCuts();
297#ifdef JJF_ZERO
298    int cutsBefore = cs.sizeCuts();
299#endif
300    CglTreeInfo info;
301    info.level = depth;
302    info.pass = pass;
303    info.formulation_rows = model_->numberRowsAtContinuous();
304    info.inTree = node != NULL;
305    if (model_->parentModel()) {
306      info.parentSolver = model_->parentModel()->continuousSolver();
307      // indicate if doing full search
308      info.hasParent = ((model_->specialOptions() & 67108864) == 0) ? 1 : 2;
309    } else {
310      info.hasParent = 0;
311      info.parentSolver = NULL;
312    }
313    info.originalColumns = model_->originalColumns();
314    info.randomNumberGenerator = randomNumberGenerator;
315    info.options = (globalCutsAtRoot()) ? 8 : 0;
316    if (ineffectualCuts())
317      info.options |= 32;
318    if (globalCuts())
319      info.options |= 16;
320    if (fullScan < 0)
321      info.options |= 128;
322    if (whetherInMustCallAgainMode())
323      info.options |= 1024;
324    // See if we want alternate set of cuts
325    if ((model_->moreSpecialOptions() & 16384) != 0)
326      info.options |= 256;
327    if (model_->parentModel())
328      info.options |= 512;
329    // above had &&!model_->parentModel()&&depth<2)
330    incrementNumberTimesEntered();
331    CglProbing *generator = dynamic_cast< CglProbing * >(generator_);
332    //if (!depth&&!pass)
333    //printf("Cut generator %s when %d\n",generatorName_,whenCutGenerator_);
334    if (!generator) {
335      // Pass across model information in case it could be useful
336      //void * saveData = solver->getApplicationData();
337      //solver->setApplicationData(model_);
338      generator_->generateCuts(*solver, cs, info);
339      //solver->setApplicationData(saveData);
340    } else {
341      // Probing - return tight column bounds
342      CglTreeProbingInfo *info2 = model_->probingInfo();
343      bool doCuts = false;
344      if (info2 && !depth) {
345        info2->options = (globalCutsAtRoot()) ? 8 : 0;
346        info2->level = depth;
347        info2->pass = pass;
348        info2->formulation_rows = model_->numberRowsAtContinuous();
349        info2->inTree = node != NULL;
350        if (model_->parentModel()) {
351          info2->parentSolver = model_->parentModel()->continuousSolver();
352          // indicate if doing full search
353          info2->hasParent = ((model_->specialOptions() & 67108864) == 0) ? 1 : 2;
354        } else {
355          info2->hasParent = 0;
356          info2->parentSolver = NULL;
357        }
358        info2->originalColumns = model_->originalColumns();
359        info2->randomNumberGenerator = randomNumberGenerator;
360        generator->generateCutsAndModify(*solver, cs, info2);
361        doCuts = true;
362      } else if (depth) {
363        /* The idea behind this is that probing may work in a different
364                   way deep in tree.  So every now and then try various
365                   combinations to see what works.
366                */
367#define TRY_NOW_AND_THEN
368#ifdef TRY_NOW_AND_THEN
369        if ((numberTimes_ == 200 || (numberTimes_ > 200 && (numberTimes_ % 2000) == 0))
370          && !model_->parentModel() && info.formulation_rows > 200) {
371          /* In tree, every now and then try various combinations
372                       maxStack, maxProbe (last 5 digits)
373                       123 is special and means CglProbing will try and
374                       be intelligent.
375                    */
376          int test[] = {
377            100123,
378            199999,
379            200123,
380            299999,
381            500123,
382            599999,
383            1000123,
384            1099999,
385            2000123,
386            2099999
387          };
388          int n = static_cast< int >(sizeof(test) / sizeof(int));
389          int saveStack = generator->getMaxLook();
390          int saveNumber = generator->getMaxProbe();
391          int kr1 = 0;
392          int kc1 = 0;
393          int bestStackTree = -1;
394          int bestNumberTree = -1;
395          for (int i = 0; i < n; i++) {
396            //OsiCuts cs2 = cs;
397            int stack = test[i] / 100000;
398            int number = test[i] - 100000 * stack;
399            generator->setMaxLook(stack);
400            generator->setMaxProbe(number);
401            int numberRowCutsBefore = cs.sizeRowCuts();
402            int numberColumnCutsBefore = cs.sizeColCuts();
403            generator_->generateCuts(*solver, cs, info);
404            int numberRowCuts = cs.sizeRowCuts() - numberRowCutsBefore;
405            int numberColumnCuts = cs.sizeColCuts() - numberColumnCutsBefore;
406#ifdef CLP_INVESTIGATE
407            if (numberRowCuts < kr1 || numberColumnCuts < kc1)
408              printf("Odd ");
409#endif
410            if (numberRowCuts > kr1 || numberColumnCuts > kc1) {
411#ifdef CLP_INVESTIGATE
412              printf("*** ");
413#endif
414              kr1 = numberRowCuts;
415              kc1 = numberColumnCuts;
416              bestStackTree = stack;
417              bestNumberTree = number;
418              doCuts = true;
419            }
420#ifdef CLP_INVESTIGATE
421            printf("maxStack %d number %d gives %d row cuts and %d column cuts\n",
422              stack, number, numberRowCuts, numberColumnCuts);
423#endif
424          }
425          generator->setMaxLook(saveStack);
426          generator->setMaxProbe(saveNumber);
427          if (bestStackTree > 0) {
428            generator->setMaxLook(bestStackTree);
429            generator->setMaxProbe(bestNumberTree);
430#ifdef CLP_INVESTIGATE
431            printf("RRNumber %d -> %d, stack %d -> %d\n",
432              saveNumber, bestNumberTree, saveStack, bestStackTree);
433#endif
434          } else {
435            // no good
436            generator->setMaxLook(0);
437#ifdef CLP_INVESTIGATE
438            printf("RRSwitching off number %d -> %d, stack %d -> %d\n",
439              saveNumber, saveNumber, saveStack, 1);
440#endif
441          }
442        }
443#endif
444        if (generator->getMaxLook() > 0 && !doCuts) {
445          generator->generateCutsAndModify(*solver, cs, &info);
446          doCuts = true;
447        }
448      } else {
449        // at root - don't always do
450        if (pass < 15 || (pass & 1) == 0) {
451          generator->generateCutsAndModify(*solver, cs, &info);
452          doCuts = true;
453        }
454      }
455      if (doCuts && generator->tightLower()) {
456        // probing may have tightened bounds - check
457        const double *tightLower = generator->tightLower();
458        const double *lower = solver->getColLower();
459        const double *tightUpper = generator->tightUpper();
460        const double *upper = solver->getColUpper();
461        const double *solution = solver->getColSolution();
462        int j;
463        int numberColumns = solver->getNumCols();
464        double primalTolerance = 1.0e-8;
465        const char *tightenBounds = generator->tightenBounds();
466#ifdef CGL_DEBUG
467        const OsiRowCutDebugger *debugger = solver->getRowCutDebugger();
468        if (debugger && debugger->onOptimalPath(*solver)) {
469          printf("On optimal path CbcCut\n");
470          int nCols = solver->getNumCols();
471          int i;
472          const double *optimal = debugger->optimalSolution();
473          const double *objective = solver->getObjCoefficients();
474          double objval1 = 0.0, objval2 = 0.0;
475          for (i = 0; i < nCols; i++) {
476#if CGL_DEBUG > 1
477            printf("%d %g %g %g %g\n", i, lower[i], solution[i], upper[i], optimal[i]);
478#endif
479            objval1 += solution[i] * objective[i];
480            objval2 += optimal[i] * objective[i];
481            assert(optimal[i] >= lower[i] - 1.0e-5 && optimal[i] <= upper[i] + 1.0e-5);
482            assert(optimal[i] >= tightLower[i] - 1.0e-5 && optimal[i] <= tightUpper[i] + 1.0e-5);
483          }
484          printf("current obj %g, integer %g\n", objval1, objval2);
485        }
486#endif
487        bool feasible = true;
488        if ((model_->getThreadMode() & 2) == 0) {
489          for (j = 0; j < numberColumns; j++) {
490            if (solver->isInteger(j)) {
491              if (tightUpper[j] < upper[j]) {
492                double nearest = floor(tightUpper[j] + 0.5);
493                //assert (fabs(tightUpper[j]-nearest)<1.0e-5); may be infeasible
494                solver->setColUpper(j, nearest);
495                if (nearest < solution[j] - primalTolerance)
496                  returnCode = true;
497              }
498              if (tightLower[j] > lower[j]) {
499                double nearest = floor(tightLower[j] + 0.5);
500                //assert (fabs(tightLower[j]-nearest)<1.0e-5); may be infeasible
501                solver->setColLower(j, nearest);
502                if (nearest > solution[j] + primalTolerance)
503                  returnCode = true;
504              }
505            } else {
506              if (upper[j] > lower[j]) {
507                if (tightUpper[j] == tightLower[j]) {
508                  // fix
509                  //if (tightLower[j]!=lower[j])
510                  solver->setColLower(j, tightLower[j]);
511                  //if (tightUpper[j]!=upper[j])
512                  solver->setColUpper(j, tightUpper[j]);
513                  if (tightLower[j] > solution[j] + primalTolerance || tightUpper[j] < solution[j] - primalTolerance)
514                    returnCode = true;
515                } else if (tightenBounds && tightenBounds[j]) {
516                  solver->setColLower(j, CoinMax(tightLower[j], lower[j]));
517                  solver->setColUpper(j, CoinMin(tightUpper[j], upper[j]));
518                  if (tightLower[j] > solution[j] + primalTolerance || tightUpper[j] < solution[j] - primalTolerance)
519                    returnCode = true;
520                }
521              }
522            }
523            if (upper[j] < lower[j] - 1.0e-3) {
524              feasible = false;
525              break;
526            }
527          }
528        } else {
529          CoinPackedVector lbs;
530          CoinPackedVector ubs;
531          int numberChanged = 0;
532          bool ifCut = false;
533          for (j = 0; j < numberColumns; j++) {
534            if (solver->isInteger(j)) {
535              if (tightUpper[j] < upper[j]) {
536                double nearest = floor(tightUpper[j] + 0.5);
537                //assert (fabs(tightUpper[j]-nearest)<1.0e-5); may be infeasible
538                ubs.insert(j, nearest);
539                numberChanged++;
540                if (nearest < solution[j] - primalTolerance)
541                  ifCut = true;
542              }
543              if (tightLower[j] > lower[j]) {
544                double nearest = floor(tightLower[j] + 0.5);
545                //assert (fabs(tightLower[j]-nearest)<1.0e-5); may be infeasible
546                lbs.insert(j, nearest);
547                numberChanged++;
548                if (nearest > solution[j] + primalTolerance)
549                  ifCut = true;
550              }
551            } else {
552              if (upper[j] > lower[j]) {
553                if (tightUpper[j] == tightLower[j]) {
554                  // fix
555                  lbs.insert(j, tightLower[j]);
556                  ubs.insert(j, tightUpper[j]);
557                  if (tightLower[j] > solution[j] + primalTolerance || tightUpper[j] < solution[j] - primalTolerance)
558                    ifCut = true;
559                } else if (tightenBounds && tightenBounds[j]) {
560                  lbs.insert(j, CoinMax(tightLower[j], lower[j]));
561                  ubs.insert(j, CoinMin(tightUpper[j], upper[j]));
562                  if (tightLower[j] > solution[j] + primalTolerance || tightUpper[j] < solution[j] - primalTolerance)
563                    ifCut = true;
564                }
565              }
566            }
567            if (upper[j] < lower[j] - 1.0e-3) {
568              feasible = false;
569              break;
570            }
571          }
572          if (numberChanged) {
573            OsiColCut cc;
574            cc.setUbs(ubs);
575            cc.setLbs(lbs);
576            if (ifCut) {
577              cc.setEffectiveness(100.0);
578            } else {
579              cc.setEffectiveness(1.0e-5);
580            }
581            cs.insert(cc);
582          }
583        }
584        if (!feasible) {
585          // not feasible -add infeasible cut
586          OsiRowCut rc;
587          rc.setLb(COIN_DBL_MAX);
588          rc.setUb(0.0);
589          cs.insert(rc);
590        }
591      }
592      //if (!solver->basisIsAvailable())
593      //returnCode=true;
594      if (!returnCode) {
595        // bounds changed but still optimal
596#ifdef COIN_HAS_CLP
597        OsiClpSolverInterface *clpSolver
598          = dynamic_cast< OsiClpSolverInterface * >(solver);
599        if (clpSolver) {
600          clpSolver->setLastAlgorithm(2);
601        }
602#endif
603      }
604#ifdef JJF_ZERO
605      // Pass across info to pseudocosts
606      char *mark = new char[numberColumns];
607      memset(mark, 0, numberColumns);
608      int nLook = generator->numberThisTime();
609      const int *lookedAt = generator->lookedAt();
610      const int *fixedDown = generator->fixedDown();
611      const int *fixedUp = generator->fixedUp();
612      for (j = 0; j < nLook; j++)
613        mark[lookedAt[j]] = 1;
614      int numberObjects = model_->numberObjects();
615      for (int i = 0; i < numberObjects; i++) {
616        CbcSimpleIntegerDynamicPseudoCost *obj1 = dynamic_cast< CbcSimpleIntegerDynamicPseudoCost * >(model_->modifiableObject(i));
617        if (obj1) {
618          int iColumn = obj1->columnNumber();
619          if (mark[iColumn])
620            obj1->setProbingInformation(fixedDown[iColumn], fixedUp[iColumn]);
621        }
622      }
623      delete[] mark;
624#endif
625    }
626    CbcCutModifier *modifier = model_->cutModifier();
627    if (modifier) {
628      int numberRowCutsAfter = cs.sizeRowCuts();
629      int k;
630      int nOdd = 0;
631      //const OsiSolverInterface * solver = model_->solver();
632      for (k = numberRowCutsAfter - 1; k >= numberRowCutsBefore; k--) {
633        OsiRowCut &thisCut = cs.rowCut(k);
634        int returnCode = modifier->modify(solver, thisCut);
635        if (returnCode) {
636          nOdd++;
637          if (returnCode == 3)
638            cs.eraseRowCut(k);
639        }
640      }
641      if (nOdd)
642        COIN_DETAIL_PRINT(printf("Cut generator %s produced %d cuts of which %d were modified\n",
643          generatorName_, numberRowCutsAfter - numberRowCutsBefore, nOdd));
644    }
645    {
646      // make all row cuts without test for duplicate
647      int numberRowCutsAfter = cs.sizeRowCuts();
648      int k;
649#ifdef CGL_DEBUG
650      const OsiRowCutDebugger *debugger = solver->getRowCutDebugger();
651#endif
652      //#define WEAKEN_CUTS 1
653#ifdef WEAKEN_CUTS
654      const double *lower = solver->getColLower();
655      const double *upper = solver->getColUpper();
656      const double *solution = solver->getColSolution();
657#endif
658      for (k = numberRowCutsBefore; k < numberRowCutsAfter; k++) {
659        OsiRowCut *thisCut = cs.rowCutPtr(k);
660#ifdef WEAKEN_CUTS
661        // weaken cut if coefficients not integer
662
663        double lb = thisCut->lb();
664        double ub = thisCut->ub();
665        if (lb < -1.0e100 || ub > 1.0e100) {
666          // normal cut
667          CoinPackedVector rpv = thisCut->row();
668          const int n = rpv.getNumElements();
669          const int *indices = rpv.getIndices();
670          const double *elements = rpv.getElements();
671          double bound = 0.0;
672          double sum = 0.0;
673          bool integral = true;
674          int nInteger = 0;
675          for (int k = 0; k < n; k++) {
676            double value = fabs(elements[k]);
677            int column = indices[k];
678            sum += value;
679            if (value != floor(value + 0.5))
680              integral = false;
681            if (solver->isInteger(column)) {
682              nInteger++;
683              double largerBound = CoinMax(fabs(lower[column]),
684                fabs(upper[column]));
685              double solutionBound = fabs(solution[column]) + 10.0;
686              bound += CoinMin(largerBound, solutionBound);
687            }
688          }
689#if WEAKEN_CUTS == 1
690          // leave if all 0-1
691          if (nInteger == bound)
692            integral = true;
693#endif
694          if (!integral) {
695            double weakenBy = 1.0e-7 * (bound + sum);
696#if WEAKEN_CUTS > 2
697            weakenBy *= 10.0;
698#endif
699            if (lb < -1.0e100)
700              thisCut->setUb(ub + weakenBy);
701            else
702              thisCut->setLb(lb - weakenBy);
703          }
704        }
705#endif
706#ifdef CGL_DEBUG
707        if (debugger && debugger->onOptimalPath(*solver)) {
708#if CGL_DEBUG > 1
709          const double *optimal = debugger->optimalSolution();
710          CoinPackedVector rpv = thisCut->row();
711          const int n = rpv.getNumElements();
712          const int *indices = rpv.getIndices();
713          const double *elements = rpv.getElements();
714
715          double lb = thisCut->lb();
716          double ub = thisCut->ub();
717          double sum = 0.0;
718
719          for (int k = 0; k < n; k++) {
720            int column = indices[k];
721            sum += optimal[column] * elements[k];
722          }
723          // is it nearly violated
724          if (sum > ub - 1.0e-8 || sum < lb + 1.0e-8) {
725            double violation = CoinMax(sum - ub, lb - sum);
726            std::cout << generatorName_ << " cut with " << n
727                      << " coefficients, nearly cuts off known solutions by " << violation
728                      << ", lo=" << lb << ", ub=" << ub << std::endl;
729            for (int k = 0; k < n; k++) {
730              int column = indices[k];
731              std::cout << "( " << column << " , " << elements[k] << " ) ";
732              if ((k % 4) == 3)
733                std::cout << std::endl;
734            }
735            std::cout << std::endl;
736            std::cout << "Non zero solution values are" << std::endl;
737            int j = 0;
738            for (int k = 0; k < n; k++) {
739              int column = indices[k];
740              if (fabs(optimal[column]) > 1.0e-9) {
741                std::cout << "( " << column << " , " << optimal[column] << " ) ";
742                if ((j % 4) == 3)
743                  std::cout << std::endl;
744                j++;
745              }
746            }
747            std::cout << std::endl;
748          }
749#endif
750          assert(!debugger->invalidCut(*thisCut));
751          if (debugger->invalidCut(*thisCut))
752            abort();
753        }
754#endif
755        thisCut->mutableRow().setTestForDuplicateIndex(false);
756      }
757    }
758    // Add in saved cuts if violated
759    if (false && !depth) {
760      const double *solution = solver->getColSolution();
761      double primalTolerance = 1.0e-7;
762      int numberCuts = savedCuts_.sizeRowCuts();
763      for (int k = numberCuts - 1; k >= 0; k--) {
764        const OsiRowCut *thisCut = savedCuts_.rowCutPtr(k);
765        double sum = 0.0;
766        int n = thisCut->row().getNumElements();
767        const int *column = thisCut->row().getIndices();
768        const double *element = thisCut->row().getElements();
769        assert(n);
770        for (int i = 0; i < n; i++) {
771          double value = element[i];
772          sum += value * solution[column[i]];
773        }
774        if (sum > thisCut->ub() + primalTolerance) {
775          sum = sum - thisCut->ub();
776        } else if (sum < thisCut->lb() - primalTolerance) {
777          sum = thisCut->lb() - sum;
778        } else {
779          sum = 0.0;
780        }
781        if (sum) {
782          // add to candidates and take out here
783          cs.insert(*thisCut);
784          savedCuts_.eraseRowCut(k);
785        }
786      }
787    }
788    if (!atSolution()) {
789      int numberRowCutsAfter = cs.sizeRowCuts();
790      int k;
791      int nEls = 0;
792      int nCuts = numberRowCutsAfter - numberRowCutsBefore;
793      // Remove NULL cuts!
794      int nNull = 0;
795      const double *solution = solver->getColSolution();
796      bool feasible = true;
797      double primalTolerance = 1.0e-7;
798      int shortCut = (depth) ? -1 : generator_->maximumLengthOfCutInTree();
799      for (k = numberRowCutsAfter - 1; k >= numberRowCutsBefore; k--) {
800        const OsiRowCut *thisCut = cs.rowCutPtr(k);
801        double sum = 0.0;
802        if (thisCut->lb() <= thisCut->ub()) {
803          int n = thisCut->row().getNumElements();
804          if (n <= shortCut)
805            numberShortCutsAtRoot_++;
806          const int *column = thisCut->row().getIndices();
807          const double *element = thisCut->row().getElements();
808          if (n <= 0) {
809            // infeasible cut - give up
810            feasible = false;
811            break;
812          }
813          nEls += n;
814          for (int i = 0; i < n; i++) {
815            double value = element[i];
816            sum += value * solution[column[i]];
817          }
818          if (sum > thisCut->ub() + primalTolerance) {
819            sum = sum - thisCut->ub();
820          } else if (sum < thisCut->lb() - primalTolerance) {
821            sum = thisCut->lb() - sum;
822          } else {
823            sum = 0.0;
824            cs.eraseRowCut(k);
825            nNull++;
826          }
827        }
828      }
829      //if (nNull)
830      //printf("%s has %d cuts and %d elements - %d null!\n",generatorName_,
831      //       nCuts,nEls,nNull);
832      numberRowCutsAfter = cs.sizeRowCuts();
833      nCuts = numberRowCutsAfter - numberRowCutsBefore;
834      nEls = 0;
835      for (k = numberRowCutsBefore; k < numberRowCutsAfter; k++) {
836        const OsiRowCut *thisCut = cs.rowCutPtr(k);
837        int n = thisCut->row().getNumElements();
838        nEls += n;
839      }
840      //printf("%s has %d cuts and %d elements\n",generatorName_,
841      //     nCuts,nEls);
842      CoinBigIndex nElsNow = solver->getMatrixByCol()->getNumElements();
843      int numberColumns = solver->getNumCols();
844      int numberRows = solver->getNumRows();
845      //double averagePerRow = static_cast<double>(nElsNow)/
846      //static_cast<double>(numberRows);
847      CoinBigIndex nAdd;
848      CoinBigIndex nAdd2;
849      CoinBigIndex nReasonable;
850      if (!model_->parentModel() && depth < 2) {
851        if (inaccuracy_ < 3) {
852          nAdd = 10000;
853          if (pass > 0 && numberColumns > -500)
854            nAdd = CoinMin(nAdd, nElsNow + 2 * numberRows);
855        } else {
856          nAdd = 10000;
857          if (pass > 0)
858            nAdd = CoinMin(nAdd, nElsNow + 2 * numberRows);
859        }
860        nAdd2 = 5 * numberColumns;
861        nReasonable = CoinMax(nAdd2, nElsNow / 8 + nAdd);
862        if (!depth && !pass) {
863          // allow more
864          nAdd += nElsNow / 2;
865          nAdd2 += nElsNow / 2;
866          nReasonable += nElsNow / 2;
867        }
868        //if (!depth&&ineffectualCuts())
869        //nReasonable *= 2;
870      } else {
871        nAdd = 200;
872        nAdd2 = 2 * numberColumns;
873        nReasonable = CoinMax(nAdd2, nElsNow / 8 + nAdd);
874      }
875      //#define UNS_WEIGHT 0.1
876#ifdef UNS_WEIGHT
877      const double *colLower = solver->getColLower();
878      const double *colUpper = solver->getColUpper();
879#endif
880      if (/*nEls>CoinMax(nAdd2,nElsNow/8+nAdd)*/ nCuts && feasible) {
881        //printf("need to remove cuts\n");
882        // just add most effective
883#ifndef JJF_ONE
884        CoinBigIndex nDelete = nEls - nReasonable;
885
886        nElsNow = nEls;
887        double *sort = new double[nCuts];
888        int *which = new int[nCuts];
889        // For parallel cuts
890        double *element2 = new double[numberColumns];
891        //#define USE_OBJECTIVE 2
892#ifdef USE_OBJECTIVE
893        const double *objective = solver->getObjCoefficients();
894#if USE_OBJECTIVE > 1
895        double objNorm = 0.0;
896        for (int i = 0; i < numberColumns; i++)
897          objNorm += objective[i] * objective[i];
898        if (objNorm)
899          objNorm = 1.0 / sqrt(objNorm);
900        else
901          objNorm = 1.0;
902        objNorm *= 0.01; // downgrade
903#endif
904#endif
905        CoinZeroN(element2, numberColumns);
906        for (k = numberRowCutsBefore; k < numberRowCutsAfter; k++) {
907          const OsiRowCut *thisCut = cs.rowCutPtr(k);
908          double sum = 0.0;
909          if (thisCut->lb() <= thisCut->ub()) {
910            int n = thisCut->row().getNumElements();
911            const int *column = thisCut->row().getIndices();
912            const double *element = thisCut->row().getElements();
913            assert(n);
914#ifdef UNS_WEIGHT
915            double normU = 0.0;
916            double norm = 1.0e-3;
917            int nU = 0;
918            for (int i = 0; i < n; i++) {
919              double value = element[i];
920              int iColumn = column[i];
921              double solValue = solution[iColumn];
922              sum += value * solValue;
923              value *= value;
924              norm += value;
925              if (solValue > colLower[iColumn] + 1.0e-6 && solValue < colUpper[iColumn] - 1.0e-6) {
926                normU += value;
927                nU++;
928              }
929            }
930#ifdef JJF_ZERO
931            int nS = n - nU;
932            if (numberColumns > 20000) {
933              if (nS > 50) {
934                double ratio = 50.0 / nS;
935                normU /= ratio;
936              }
937            }
938#endif
939            norm += UNS_WEIGHT * (normU - norm);
940#else
941            double norm = 1.0e-3;
942#ifdef USE_OBJECTIVE
943            double obj = 0.0;
944#endif
945            for (int i = 0; i < n; i++) {
946              int iColumn = column[i];
947              double value = element[i];
948              sum += value * solution[iColumn];
949              norm += value * value;
950#ifdef USE_OBJECTIVE
951              obj += value * objective[iColumn];
952#endif
953            }
954#endif
955            if (sum > thisCut->ub()) {
956              sum = sum - thisCut->ub();
957            } else if (sum < thisCut->lb()) {
958              sum = thisCut->lb() - sum;
959            } else {
960              sum = 0.0;
961            }
962#ifdef USE_OBJECTIVE
963            if (sum) {
964#if USE_OBJECTIVE == 1
965              obj = CoinMax(1.0e-6, fabs(obj));
966              norm = sqrt(obj * norm);
967              //sum += fabs(obj)*invObjNorm;
968              //printf("sum %g norm %g normobj %g invNorm %g mod %g\n",
969              //     sum,norm,obj,invObjNorm,obj*invObjNorm);
970              // normalize
971              sum /= sqrt(norm);
972#else
973              // normalize
974              norm = 1.0 / sqrt(norm);
975              sum = (sum + objNorm * obj) * norm;
976#endif
977            }
978#else
979            // normalize
980            sum /= sqrt(norm);
981#endif
982            //sum /= pow(norm,0.3);
983            // adjust for length
984            //sum /= pow(reinterpret_cast<double>(n),0.2);
985            //sum /= sqrt((double) n);
986            // randomize
987            //double randomNumber =
988            //model_->randomNumberGenerator()->randomDouble();
989            //sum *= (0.5+randomNumber);
990          } else {
991            // keep
992            sum = COIN_DBL_MAX;
993          }
994          sort[k - numberRowCutsBefore] = sum;
995          which[k - numberRowCutsBefore] = k;
996        }
997        CoinSort_2(sort, sort + nCuts, which);
998        // Now see which ones are too similar
999        int nParallel = 0;
1000        double testValue = (depth > 1) ? 0.99 : 0.999999;
1001        for (k = 0; k < nCuts; k++) {
1002          int j = which[k];
1003          const OsiRowCut *thisCut = cs.rowCutPtr(j);
1004          if (thisCut->lb() > thisCut->ub())
1005            break; // cut is infeasible
1006          int n = thisCut->row().getNumElements();
1007          const int *column = thisCut->row().getIndices();
1008          const double *element = thisCut->row().getElements();
1009          assert(n);
1010          double norm = 0.0;
1011          double lb = thisCut->lb();
1012          double ub = thisCut->ub();
1013          for (int i = 0; i < n; i++) {
1014            double value = element[i];
1015            element2[column[i]] = value;
1016            norm += value * value;
1017          }
1018          int kkk = CoinMin(nCuts, k + 5);
1019          for (int kk = k + 1; kk < kkk; kk++) {
1020            int jj = which[kk];
1021            const OsiRowCut *thisCut2 = cs.rowCutPtr(jj);
1022            if (thisCut2->lb() > thisCut2->ub())
1023              break; // cut is infeasible
1024            int nB = thisCut2->row().getNumElements();
1025            const int *columnB = thisCut2->row().getIndices();
1026            const double *elementB = thisCut2->row().getElements();
1027            assert(nB);
1028            double normB = 0.0;
1029            double product = 0.0;
1030            for (int i = 0; i < nB; i++) {
1031              double value = elementB[i];
1032              normB += value * value;
1033              product += value * element2[columnB[i]];
1034            }
1035            if (product > 0.0 && product * product > testValue * norm * normB) {
1036              bool parallel = true;
1037              double lbB = thisCut2->lb();
1038              double ubB = thisCut2->ub();
1039              if ((lb < -1.0e20 && lbB > -1.0e20) || (lbB < -1.0e20 && lb > -1.0e20))
1040                parallel = false;
1041              double tolerance;
1042              tolerance = CoinMax(fabs(lb), fabs(lbB)) + 1.0e-6;
1043              if (fabs(lb - lbB) > tolerance)
1044                parallel = false;
1045              if ((ub > 1.0e20 && ubB < 1.0e20) || (ubB > 1.0e20 && ub < 1.0e20))
1046                parallel = false;
1047              tolerance = CoinMax(fabs(ub), fabs(ubB)) + 1.0e-6;
1048              if (fabs(ub - ubB) > tolerance)
1049                parallel = false;
1050              if (parallel) {
1051                nParallel++;
1052                sort[k] = 0.0;
1053                break;
1054              }
1055            }
1056          }
1057          for (int i = 0; i < n; i++) {
1058            element2[column[i]] = 0.0;
1059          }
1060        }
1061        delete[] element2;
1062        CoinSort_2(sort, sort + nCuts, which);
1063        k = 0;
1064        while (nDelete > 0 || !sort[k]) {
1065          int iCut = which[k];
1066          const OsiRowCut *thisCut = cs.rowCutPtr(iCut);
1067          int n = thisCut->row().getNumElements();
1068          // may be best, just to save if short
1069          if (false && n && sort[k]) {
1070            // add to saved cuts
1071            savedCuts_.insert(*thisCut);
1072          }
1073          nDelete -= n;
1074          k++;
1075          if (k >= nCuts)
1076            break;
1077        }
1078        std::sort(which, which + k);
1079        k--;
1080        for (; k >= 0; k--) {
1081          cs.eraseRowCut(which[k]);
1082        }
1083        delete[] sort;
1084        delete[] which;
1085        numberRowCutsAfter = cs.sizeRowCuts();
1086#else
1087        double *norm = new double[nCuts];
1088        int *which = new int[2 * nCuts];
1089        double *score = new double[nCuts];
1090        double *ortho = new double[nCuts];
1091        int nIn = 0;
1092        int nOut = nCuts;
1093        // For parallel cuts
1094        double *element2 = new double[numberColumns];
1095        const double *objective = solver->getObjCoefficients();
1096        double objNorm = 0.0;
1097        for (int i = 0; i < numberColumns; i++)
1098          objNorm += objective[i] * objective[i];
1099        if (objNorm)
1100          objNorm = 1.0 / sqrt(objNorm);
1101        else
1102          objNorm = 1.0;
1103        objNorm *= 0.1; // weight of 0.1
1104        CoinZeroN(element2, numberColumns);
1105        int numberRowCuts = numberRowCutsAfter - numberRowCutsBefore;
1106        int iBest = -1;
1107        double best = 0.0;
1108        int nPossible = 0;
1109        double testValue = (depth > 1) ? 0.7 : 0.5;
1110        for (k = 0; k < numberRowCuts; k++) {
1111          const OsiRowCut *thisCut = cs.rowCutPtr(k + numberRowCutsBefore);
1112          double sum = 0.0;
1113          if (thisCut->lb() <= thisCut->ub()) {
1114            int n = thisCut->row().getNumElements();
1115            const int *column = thisCut->row().getIndices();
1116            const double *element = thisCut->row().getElements();
1117            assert(n);
1118            double normThis = 1.0e-6;
1119            double obj = 0.0;
1120            for (int i = 0; i < n; i++) {
1121              int iColumn = column[i];
1122              double value = element[i];
1123              sum += value * solution[iColumn];
1124              normThis += value * value;
1125              obj += value * objective[iColumn];
1126            }
1127            if (sum > thisCut->ub()) {
1128              sum = sum - thisCut->ub();
1129            } else if (sum < thisCut->lb()) {
1130              sum = thisCut->lb() - sum;
1131            } else {
1132              sum = 0.0;
1133            }
1134            if (sum) {
1135              normThis = 1.0 / sqrt(normThis);
1136              norm[k] = normThis;
1137              sum *= normThis;
1138              obj *= normThis;
1139              score[k] = sum + obj * objNorm;
1140              ortho[k] = 1.0;
1141            }
1142          } else {
1143            // keep and discard others
1144            nIn = 1;
1145            which[0] = k;
1146            for (int j = 0; j < numberRowCuts; j++) {
1147              if (j != k)
1148                which[nOut++] = j;
1149            }
1150            iBest = -1;
1151            break;
1152          }
1153          if (sum) {
1154            if (score[k] > best) {
1155              best = score[k];
1156              iBest = nPossible;
1157            }
1158            which[nPossible++] = k;
1159          } else {
1160            which[nOut++] = k;
1161          }
1162        }
1163        while (iBest >= 0) {
1164          int kBest = which[iBest];
1165          int j = which[nIn];
1166          which[iBest] = j;
1167          which[nIn++] = kBest;
1168          const OsiRowCut *thisCut = cs.rowCutPtr(kBest + numberRowCutsBefore);
1169          int n = thisCut->row().getNumElements();
1170          nReasonable -= n;
1171          if (nReasonable <= 0) {
1172            for (k = nIn; k < nPossible; k++)
1173              which[nOut++] = which[k];
1174            break;
1175          }
1176          // Now see which ones are too similar and choose next
1177          iBest = -1;
1178          best = 0.0;
1179          int nOld = nPossible;
1180          nPossible = nIn;
1181          const int *column = thisCut->row().getIndices();
1182          const double *element = thisCut->row().getElements();
1183          assert(n);
1184          double normNew = norm[kBest];
1185          for (int i = 0; i < n; i++) {
1186            double value = element[i];
1187            element2[column[i]] = value;
1188          }
1189          for (int j = nIn; j < nOld; j++) {
1190            k = which[j];
1191            const OsiRowCut *thisCut2 = cs.rowCutPtr(k + numberRowCutsBefore);
1192            int nB = thisCut2->row().getNumElements();
1193            const int *columnB = thisCut2->row().getIndices();
1194            const double *elementB = thisCut2->row().getElements();
1195            assert(nB);
1196            double normB = norm[k];
1197            double product = 0.0;
1198            for (int i = 0; i < nB; i++) {
1199              double value = elementB[i];
1200              product += value * element2[columnB[i]];
1201            }
1202            double orthoScore = 1.0 - product * normNew * normB;
1203            if (orthoScore >= testValue) {
1204              ortho[k] = CoinMin(orthoScore, ortho[k]);
1205              double test = score[k] + ortho[k];
1206              if (test > best) {
1207                best = score[k];
1208                iBest = nPossible;
1209              }
1210              which[nPossible++] = k;
1211            } else {
1212              which[nOut++] = k;
1213            }
1214          }
1215          for (int i = 0; i < n; i++) {
1216            element2[column[i]] = 0.0;
1217          }
1218        }
1219        delete[] score;
1220        delete[] ortho;
1221        std::sort(which + nCuts, which + nOut);
1222        k = nOut - 1;
1223        for (; k >= nCuts; k--) {
1224          cs.eraseRowCut(which[k] + numberRowCutsBefore);
1225        }
1226        delete[] norm;
1227        delete[] which;
1228        numberRowCutsAfter = cs.sizeRowCuts();
1229#endif
1230      }
1231    }
1232#ifdef CBC_DEBUG
1233    {
1234      int numberRowCutsAfter = cs.sizeRowCuts();
1235      int k;
1236      int nBad = 0;
1237      for (k = numberRowCutsBefore; k < numberRowCutsAfter; k++) {
1238        OsiRowCut thisCut = cs.rowCut(k);
1239        if (thisCut.lb() > thisCut.ub() || thisCut.lb() > 1.0e8 || thisCut.ub() < -1.0e8)
1240          printf("cut from %s has bounds %g and %g!\n",
1241            generatorName_, thisCut.lb(), thisCut.ub());
1242        if (thisCut.lb() <= thisCut.ub()) {
1243          /* check size of elements.
1244                       We can allow smaller but this helps debug generators as it
1245                       is unsafe to have small elements */
1246          int n = thisCut.row().getNumElements();
1247          const int *column = thisCut.row().getIndices();
1248          const double *element = thisCut.row().getElements();
1249          assert(n);
1250          for (int i = 0; i < n; i++) {
1251            double value = element[i];
1252            if (fabs(value) <= 1.0e-12 || fabs(value) >= 1.0e20)
1253              nBad++;
1254          }
1255        }
1256        if (nBad)
1257          printf("Cut generator %s produced %d cuts of which %d had tiny or large elements\n",
1258            generatorName_, numberRowCutsAfter - numberRowCutsBefore, nBad);
1259      }
1260    }
1261#endif
1262    int numberRowCutsAfter = cs.sizeRowCuts();
1263    int numberColumnCutsAfter = cs.sizeColCuts();
1264    if (numberRowCutsBefore < numberRowCutsAfter) {
1265      for (int k = numberRowCutsBefore; k < numberRowCutsAfter; k++) {
1266        OsiRowCut thisCut = cs.rowCut(k);
1267        int n = thisCut.row().getNumElements();
1268        numberElements_ += n;
1269      }
1270#ifdef JJF_ZERO
1271      printf("generator %s generated %d row cuts\n",
1272        generatorName_, numberRowCutsAfter - numberRowCutsBefore);
1273#endif
1274      numberCuts_ += numberRowCutsAfter - numberRowCutsBefore;
1275    }
1276    if (numberColumnCutsBefore < numberColumnCutsAfter) {
1277#ifdef JJF_ZERO
1278      printf("generator %s generated %d column cuts\n",
1279        generatorName_, numberColumnCutsAfter - numberColumnCutsBefore);
1280#endif
1281      numberColumnCuts_ += numberColumnCutsAfter - numberColumnCutsBefore;
1282    }
1283    if (timing()) {
1284      // Using thread clock Id get the clock time
1285      /* TODO there should be check in configure or the Posix version C preprocessor variable
1286       * to decide whether pthread_getcpuclockid is available
1287       */
1288#if defined(_MSC_VER) || defined(__MSVCRT__) || defined(__APPLE__) || !defined(CBC_THREAD)
1289      timeInCutGenerator_ += CoinCpuTime() - time1;
1290#else
1291      if (!model_->getNumberThreads()) {
1292        timeInCutGenerator_ += CoinCpuTime() - time1;
1293      } else {
1294        clock_gettime(threadClockId, &currTime);
1295        timeInCutGenerator_ += static_cast<double>(currTime.tv_sec) + 1.0e-9*
1296          static_cast<double>(currTime.tv_nsec) - time1;
1297      }
1298#endif
1299    }
1300    // switch off if first time and no good
1301    if (node == NULL && !pass) {
1302      if (numberRowCutsAfter - numberRowCutsBefore
1303        < switchOffIfLessThan_ /*&& numberCuts_ < switchOffIfLessThan_*/) {
1304        // switch off
1305        maximumTries_ = 0;
1306        whenCutGenerator_ = -100;
1307        //whenCutGenerator_ = -100;
1308        //whenCutGeneratorInSub_ = -200;
1309      }
1310    }
1311    if (maximumTries_ > 0) {
1312      maximumTries_--;
1313      if (!maximumTries_)
1314        whenCutGenerator_ = -100;
1315    }
1316  }
1317  return returnCode;
1318}
1319void CbcCutGenerator::setHowOften(int howOften)
1320{
1321
1322  if (howOften >= 1000000) {
1323    // leave Probing every SCANCUTS_PROBING
1324    howOften = howOften % 1000000;
1325    CglProbing *generator = dynamic_cast< CglProbing * >(generator_);
1326
1327    if (generator && howOften > SCANCUTS_PROBING)
1328      howOften = SCANCUTS_PROBING + 1000000;
1329    else
1330      howOften += 1000000;
1331  }
1332  whenCutGenerator_ = howOften;
1333}
1334void CbcCutGenerator::setWhatDepth(int value)
1335{
1336  depthCutGenerator_ = value;
1337}
1338void CbcCutGenerator::setWhatDepthInSub(int value)
1339{
1340  depthCutGeneratorInSub_ = value;
1341}
1342// Add in statistics from other
1343void CbcCutGenerator::addStatistics(const CbcCutGenerator *other)
1344{
1345  // Time in cut generator
1346  timeInCutGenerator_ += other->timeInCutGenerator_;
1347  // Number times cut generator entered
1348  numberTimes_ += other->numberTimes_;
1349  // Total number of cuts added
1350  numberCuts_ += other->numberCuts_;
1351  // Total number of elements added
1352  numberElements_ += other->numberElements_;
1353  // Total number of column cuts added
1354  numberColumnCuts_ += other->numberColumnCuts_;
1355  // Total number of cuts active after (at end of n cut passes at each node)
1356  numberCutsActive_ += other->numberCutsActive_;
1357  // Number of cuts generated at root
1358  numberCutsAtRoot_ += other->numberCutsAtRoot_;
1359  // Number of cuts active at root
1360  numberActiveCutsAtRoot_ += other->numberActiveCutsAtRoot_;
1361  // Number of short cuts at root
1362  numberShortCutsAtRoot_ += other->numberShortCutsAtRoot_;
1363}
1364// Scale back statistics by factor
1365void CbcCutGenerator::scaleBackStatistics(int factor)
1366{
1367  // leave time
1368  // Number times cut generator entered
1369  numberTimes_ = (numberTimes_ + factor - 1) / factor;
1370  // Total number of cuts added
1371  numberCuts_ = (numberCuts_ + factor - 1) / factor;
1372  // Total number of elements added
1373  numberElements_ = (numberElements_ + factor - 1) / factor;
1374  // Total number of column cuts added
1375  numberColumnCuts_ = (numberColumnCuts_ + factor - 1) / factor;
1376  // Total number of cuts active after (at end of n cut passes at each node)
1377  numberCutsActive_ = (numberCutsActive_ + factor - 1) / factor;
1378  // Number of cuts generated at root
1379  numberCutsAtRoot_ = (numberCutsAtRoot_ + factor - 1) / factor;
1380  // Number of cuts active at root
1381  numberActiveCutsAtRoot_ = (numberActiveCutsAtRoot_ + factor - 1) / factor;
1382  // Number of short cuts at root
1383  numberShortCutsAtRoot_ = (numberShortCutsAtRoot_ + factor - 1) / factor;
1384}
1385// Create C++ lines to get to current state
1386void CbcCutGenerator::generateTuning(FILE *fp)
1387{
1388  fprintf(fp, "// Cbc tuning for generator %s\n", generatorName_);
1389  fprintf(fp, "   generator->setHowOften(%d);\n", whenCutGenerator_);
1390  fprintf(fp, "   generator->setSwitchOffIfLessThan(%d);\n", switchOffIfLessThan_);
1391  fprintf(fp, "   generator->setWhatDepth(%d);\n", depthCutGenerator_);
1392  fprintf(fp, "   generator->setInaccuracy(%d);\n", inaccuracy_);
1393  if (timing())
1394    fprintf(fp, "   generator->setTiming(true);\n");
1395  if (normal())
1396    fprintf(fp, "   generator->setNormal(true);\n");
1397  if (atSolution())
1398    fprintf(fp, "   generator->setAtSolution(true);\n");
1399  if (whenInfeasible())
1400    fprintf(fp, "   generator->setWhenInfeasible(true);\n");
1401  if (needsOptimalBasis())
1402    fprintf(fp, "   generator->setNeedsOptimalBasis(true);\n");
1403  if (mustCallAgain())
1404    fprintf(fp, "   generator->setMustCallAgain(true);\n");
1405  if (whetherToUse())
1406    fprintf(fp, "   generator->setWhetherToUse(true);\n");
1407}
1408
1409/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
1410*/
Note: See TracBrowser for help on using the repository browser.