source: branches/devel/Cbc/src/CbcCutGenerator.cpp @ 439

Last change on this file since 439 was 424, checked in by forrest, 13 years ago

many changes

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 12.1 KB
Line 
1// Copyright (C) 2003, International Business Machines
2// Corporation and others.  All Rights Reserved.
3#if defined(_MSC_VER)
4// Turn off compiler warning about long names
5#  pragma warning(disable:4786)
6#endif
7#include <cassert>
8#include <cmath>
9#include <cfloat>
10
11#include "OsiSolverInterface.hpp"
12#include "CbcModel.hpp"
13#include "CbcMessage.hpp"
14#include "CbcCutGenerator.hpp"
15#include "CbcBranchDynamic.hpp"
16#include "CglProbing.hpp"
17#include "CoinTime.hpp"
18
19// Default Constructor
20CbcCutGenerator::CbcCutGenerator ()
21  : model_(NULL),
22    generator_(NULL),
23    whenCutGenerator_(-1),
24    whenCutGeneratorInSub_(-100),
25    switchOffIfLessThan_(0),
26    depthCutGenerator_(-1),
27    depthCutGeneratorInSub_(-1),
28    generatorName_(NULL),
29    normal_(true),
30    atSolution_(false),
31    whenInfeasible_(false),
32    mustCallAgain_(false),
33    switchedOff_(false),
34    timing_(false),
35    timeInCutGenerator_(0.0),
36    numberTimes_(0),
37    numberCuts_(0),
38    numberColumnCuts_(0),
39    numberCutsActive_(0)
40{
41}
42// Normal constructor
43CbcCutGenerator::CbcCutGenerator(CbcModel * model,CglCutGenerator * generator,
44                                 int howOften, const char * name,
45                                 bool normal, bool atSolution, 
46                                 bool infeasible, int howOftenInSub,
47                                 int whatDepth, int whatDepthInSub,
48                                 int switchOffIfLessThan)
49  : 
50    depthCutGenerator_(whatDepth),
51    depthCutGeneratorInSub_(whatDepthInSub),
52    mustCallAgain_(false),
53    switchedOff_(false),
54    timing_(false),
55    timeInCutGenerator_(0.0),
56    numberTimes_(0),
57    numberCuts_(0),
58    numberColumnCuts_(0),
59    numberCutsActive_(0)
60{
61  model_ = model;
62  generator_=generator->clone();
63  generator_->refreshSolver(model_->solver());
64  whenCutGenerator_=howOften;
65  whenCutGeneratorInSub_ = howOftenInSub;
66  switchOffIfLessThan_=switchOffIfLessThan;
67  if (name)
68    generatorName_=strdup(name);
69  else
70    generatorName_ = strdup("Unknown");
71  normal_=normal;
72  atSolution_=atSolution;
73  whenInfeasible_=infeasible;
74}
75
76// Copy constructor
77CbcCutGenerator::CbcCutGenerator ( const CbcCutGenerator & rhs)
78{
79  model_ = rhs.model_;
80  generator_=rhs.generator_->clone();
81  generator_->refreshSolver(model_->solver());
82  whenCutGenerator_=rhs.whenCutGenerator_;
83  whenCutGeneratorInSub_ = rhs.whenCutGeneratorInSub_;
84  switchOffIfLessThan_ = rhs.switchOffIfLessThan_;
85  depthCutGenerator_=rhs.depthCutGenerator_;
86  depthCutGeneratorInSub_ = rhs.depthCutGeneratorInSub_;
87  generatorName_=strdup(rhs.generatorName_);
88  normal_=rhs.normal_;
89  atSolution_=rhs.atSolution_;
90  whenInfeasible_=rhs.whenInfeasible_;
91  mustCallAgain_ = rhs.mustCallAgain_;
92  switchedOff_ = rhs.switchedOff_;
93  timing_ = rhs.timing_;
94  timeInCutGenerator_ = rhs.timeInCutGenerator_;
95  numberTimes_ = rhs.numberTimes_;
96  numberCuts_ = rhs.numberCuts_;
97  numberColumnCuts_ = rhs.numberColumnCuts_;
98  numberCutsActive_ = rhs.numberCutsActive_;
99}
100
101// Assignment operator
102CbcCutGenerator & 
103CbcCutGenerator::operator=( const CbcCutGenerator& rhs)
104{
105  if (this!=&rhs) {
106    delete generator_;
107    free(generatorName_);
108    model_ = rhs.model_;
109    generator_=rhs.generator_->clone();
110    generator_->refreshSolver(model_->solver());
111    whenCutGenerator_=rhs.whenCutGenerator_;
112    whenCutGeneratorInSub_ = rhs.whenCutGeneratorInSub_;
113    switchOffIfLessThan_ = rhs.switchOffIfLessThan_;
114    depthCutGenerator_=rhs.depthCutGenerator_;
115    depthCutGeneratorInSub_ = rhs.depthCutGeneratorInSub_;
116    generatorName_=strdup(rhs.generatorName_);
117    normal_=rhs.normal_;
118    atSolution_=rhs.atSolution_;
119    whenInfeasible_=rhs.whenInfeasible_;
120    mustCallAgain_ = rhs.mustCallAgain_;
121    switchedOff_ = rhs.switchedOff_;
122    timing_ = rhs.timing_;
123    timeInCutGenerator_ = rhs.timeInCutGenerator_;
124    numberTimes_ = rhs.numberTimes_;
125    numberCuts_ = rhs.numberCuts_;
126    numberColumnCuts_ = rhs.numberColumnCuts_;
127    numberCutsActive_ = rhs.numberCutsActive_;
128  }
129  return *this;
130}
131
132// Destructor
133CbcCutGenerator::~CbcCutGenerator ()
134{
135  free(generatorName_);
136  delete generator_;
137}
138
139/* This is used to refresh any inforamtion.
140   It also refreshes the solver in the cut generator
141   in case generator wants to do some work
142*/
143void 
144CbcCutGenerator::refreshModel(CbcModel * model)
145{
146  model_=model;
147  generator_->refreshSolver(model_->solver());
148}
149/* Generate cuts for the model data contained in si.
150   The generated cuts are inserted into and returned in the
151   collection of cuts cs.
152*/
153bool
154CbcCutGenerator::generateCuts( OsiCuts & cs , bool fullScan, CbcNode * node)
155{
156  int howOften = whenCutGenerator_;
157  if (howOften==-100)
158    return false;
159  if (howOften>0)
160    howOften = howOften % 1000000;
161  else 
162    howOften=1;
163  if (!howOften)
164    howOften=1;
165  bool returnCode=false;
166  OsiSolverInterface * solver = model_->solver();
167  int depth;
168  if (node)
169    depth=node->depth();
170  else
171    depth=0;
172  int pass=model_->getCurrentPassNumber()-1;
173  bool doThis=(model_->getNodeCount()%howOften)==0;
174  if (depthCutGenerator_>0) {
175    doThis = (depth % depthCutGenerator_) ==0;
176    if (depth<depthCutGenerator_)
177      doThis=true; // and also at top of tree
178  }
179  // But turn off if 100
180  if (howOften==100)
181    doThis=false;
182  // Switch off if special setting
183  if (whenCutGeneratorInSub_==-200) {
184    fullScan=false;
185    doThis=false;
186  }
187  if (fullScan||doThis) {
188    double time1=0.0;
189    if (timing_)
190      time1 = CoinCpuTime();
191    //#define CBC_DEBUG
192#ifdef CBC_DEBUG
193    int numberRowCutsBefore = cs.sizeRowCuts() ;
194#endif
195    int cutsBefore = cs.sizeCuts();
196    CglTreeInfo info;
197    info.level = depth;
198    info.pass = pass;
199    info.formulation_rows = model_->numberRowsAtContinuous();
200    info.inTree = node!=NULL;
201    incrementNumberTimesEntered();
202    CglProbing* generator =
203      dynamic_cast<CglProbing*>(generator_);
204    if (!generator) {
205      // Pass across model information in case it could be useful
206      //OsiSolverInterface * solver = model_->solver();
207      //void * saveData = solver->getApplicationData();
208      //solver->setApplicationData(model_);
209      generator_->generateCuts(*solver,cs,info);
210      //solver->setApplicationData(saveData);
211    } else {
212      // Probing - return tight column bounds
213      CglTreeProbingInfo * info2 = model_->probingInfo();
214      if (info2) {
215        info2->level = depth;
216        info2->pass = pass;
217        info2->formulation_rows = model_->numberRowsAtContinuous();
218        info2->inTree = node!=NULL;
219        generator->generateCutsAndModify(*solver,cs,info2);
220      } else {
221        generator->generateCutsAndModify(*solver,cs,&info);
222      }
223      const double * tightLower = generator->tightLower();
224      const double * lower = solver->getColLower();
225      const double * tightUpper = generator->tightUpper();
226      const double * upper = solver->getColUpper();
227      const double * solution = solver->getColSolution();
228      int j;
229      int numberColumns = solver->getNumCols();
230      double primalTolerance = 1.0e-8;
231#if 0
232      int numberChanged=0,ifCut=0;
233      CoinPackedVector lbs;
234      CoinPackedVector ubs;
235      for (j=0;j<numberColumns;j++) {
236        if (solver->isInteger(j)) {
237          if (tightUpper[j]<upper[j]) {
238            numberChanged++;
239            assert (tightUpper[j]==floor(tightUpper[j]+0.5));
240            ubs.insert(j,tightUpper[j]);
241            if (tightUpper[j]<solution[j]-primalTolerance)
242              ifCut=1;
243          }
244          if (tightLower[j]>lower[j]) {
245            numberChanged++;
246            assert (tightLower[j]==floor(tightLower[j]+0.5));
247            lbs.insert(j,tightLower[j]);
248            if (tightLower[j]>solution[j]+primalTolerance)
249              ifCut=1;
250          }
251        } else {
252          if (tightUpper[j]==tightLower[j]&&
253              upper[j]>lower[j]) {
254            // fix
255            //solver->setColLower(j,tightLower[j]);
256            //solver->setColUpper(j,tightUpper[j]);
257            double value = tightUpper[j];
258            numberChanged++;
259            if (value<upper[j])
260              ubs.insert(j,value);
261            if (value>lower[j])
262              lbs.insert(j,value);
263          }
264        }
265      }
266      if (numberChanged) {
267        OsiColCut cc;
268        cc.setUbs(ubs);
269        cc.setLbs(lbs);
270        if (ifCut) {
271          cc.setEffectiveness(100.0);
272        } else {
273          cc.setEffectiveness(1.0e-5);
274        }
275        cs.insert(cc);
276      }
277      // need to resolve if some bounds changed
278      returnCode = !solver->basisIsAvailable();
279      assert (!returnCode);
280#else
281      const char * tightenBounds = generator->tightenBounds();
282      for (j=0;j<numberColumns;j++) {
283        if (solver->isInteger(j)) {
284          if (tightUpper[j]<upper[j]) {
285            double nearest = floor(tightUpper[j]+0.5);
286            //assert (fabs(tightUpper[j]-nearest)<1.0e-5); may be infeasible
287            solver->setColUpper(j,nearest);
288            if (nearest<solution[j]-primalTolerance)
289              returnCode=true;
290          }
291          if (tightLower[j]>lower[j]) {
292            double nearest = floor(tightLower[j]+0.5);
293            //assert (fabs(tightLower[j]-nearest)<1.0e-5); may be infeasible
294            solver->setColLower(j,nearest);
295            if (nearest>solution[j]+primalTolerance)
296              returnCode=true;
297          }
298        } else {
299          if (upper[j]>lower[j]) {
300            if (tightUpper[j]==tightLower[j]) {
301              // fix
302              solver->setColLower(j,tightLower[j]);
303              solver->setColUpper(j,tightUpper[j]);
304              if (tightLower[j]>solution[j]+primalTolerance||
305                  tightUpper[j]<solution[j]-primalTolerance)
306                returnCode=true;
307            } else if (tightenBounds&&tightenBounds[j]) {
308              solver->setColLower(j,CoinMax(tightLower[j],lower[j]));
309              solver->setColUpper(j,CoinMin(tightUpper[j],upper[j]));
310              if (tightLower[j]>solution[j]+primalTolerance||
311                  tightUpper[j]<solution[j]-primalTolerance)
312                returnCode=true;
313            }
314          }
315        }
316      }
317      //if (!solver->basisIsAvailable())
318      //returnCode=true;
319#endif
320#if 0
321      // Pass across info to pseudocosts
322      char * mark = new char[numberColumns];
323      memset(mark,0,numberColumns);
324      int nLook = generator->numberThisTime();
325      const int * lookedAt = generator->lookedAt();
326      const int * fixedDown = generator->fixedDown();
327      const int * fixedUp = generator->fixedUp();
328      for (j=0;j<nLook;j++)
329        mark[lookedAt[j]]=1;
330      int numberObjects = model_->numberObjects();
331      for (int i=0;i<numberObjects;i++) {
332        CbcSimpleIntegerDynamicPseudoCost * obj1 =
333          dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(model_->modifiableObject(i)) ;
334        if (obj1) {
335          int iColumn = obj1->columnNumber();
336          if (mark[iColumn])
337            obj1->setProbingInformation(fixedDown[iColumn],fixedUp[iColumn]);
338        }
339      }
340      delete [] mark;
341#endif
342    }
343#ifdef CBC_DEBUG
344    {
345      int numberRowCutsAfter = cs.sizeRowCuts() ;
346      int k ;
347      int nBad=0;
348      for (k = numberRowCutsBefore;k<numberRowCutsAfter;k++) {
349        OsiRowCut thisCut = cs.rowCut(k) ;
350        if (thisCut.lb()<=thisCut.ub()) {
351          /* check size of elements.
352             We can allow smaller but this helps debug generators as it
353             is unsafe to have small elements */
354          int n=thisCut.row().getNumElements();
355          const int * column = thisCut.row().getIndices();
356          const double * element = thisCut.row().getElements();
357          assert (n);
358          for (int i=0;i<n;i++) {
359            double value = element[i];
360            if (fabs(value)<=1.0e-12||fabs(value)>=1.0e20)
361              nBad++;
362          }
363        }
364        if (nBad) 
365          printf("Cut generator %s produced %d cuts of which %d had tiny or large elements\n",
366                 generatorName_,numberRowCutsAfter-numberRowCutsBefore,nBad);
367      }
368    }
369#endif
370    if (timing_)
371      timeInCutGenerator_ += CoinCpuTime()-time1;
372    // switch off if first time and no good
373    if (node==NULL&&!pass) {
374      if (cs.sizeCuts()-cutsBefore<CoinAbs(switchOffIfLessThan_)) {
375        whenCutGenerator_=-99;
376        whenCutGeneratorInSub_ = -200;
377      }
378    }
379  }
380  return returnCode;
381}
382void 
383CbcCutGenerator::setHowOften(int howOften) 
384{
385 
386  if (howOften>=1000000) {
387    // leave Probing every 10
388    howOften = howOften % 1000000;
389    CglProbing* generator =
390      dynamic_cast<CglProbing*>(generator_);
391   
392    if (generator&&howOften>10) 
393      howOften=10+1000000;
394    else
395      howOften += 1000000;
396  }
397  whenCutGenerator_ = howOften;
398}
399void 
400CbcCutGenerator::setWhatDepth(int value) 
401{
402  depthCutGenerator_ = value;
403}
404void 
405CbcCutGenerator::setWhatDepthInSub(int value) 
406{
407  depthCutGeneratorInSub_ = value;
408}
Note: See TracBrowser for help on using the repository browser.