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

Last change on this file since 2340 was 2340, checked in by forrest, 3 years ago

set info in CglTreeInfo?

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