source: trunk/Cbc/src/CbcCutGenerator.cpp @ 2344

Last change on this file since 2344 was 2344, checked in by forrest, 2 years ago

change int to CoinBigIndex?

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