1 // $Id$
2 // Copyright (C) 2006, 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 // edwin 12/5/09 carved out of CbcHeuristicRINS
7 
8 #if defined(_MSC_VER)
9 // Turn off compiler warning about long names
10 #pragma warning(disable : 4786)
11 #endif
12 #include <cassert>
13 #include <cstdlib>
14 #include <cmath>
15 #include <cfloat>
16 
17 #include "OsiSolverInterface.hpp"
18 #include "CbcModel.hpp"
19 #include "CbcMessage.hpp"
20 #include "CbcHeuristicVND.hpp"
21 #include "CbcBranchActual.hpp"
22 #include "CbcStrategy.hpp"
23 #include "CglPreProcess.hpp"
24 
25 // Default Constructor
CbcHeuristicVND()26 CbcHeuristicVND::CbcHeuristicVND()
27   : CbcHeuristic()
28 {
29   numberSolutions_ = 0;
30   numberSuccesses_ = 0;
31   numberTries_ = 0;
32   lastNode_ = -999999;
33   howOften_ = 100;
34   decayFactor_ = 0.5;
35   baseSolution_ = NULL;
36   whereFrom_ = 1 + 8 + 255 * 256;
37   stepSize_ = 0;
38   k_ = 0;
39   kmax_ = 0;
40   nDifferent_ = 0;
41 }
42 
43 // Constructor with model - assumed before cuts
44 
CbcHeuristicVND(CbcModel & model)45 CbcHeuristicVND::CbcHeuristicVND(CbcModel &model)
46   : CbcHeuristic(model)
47 {
48   numberSolutions_ = 0;
49   numberSuccesses_ = 0;
50   numberTries_ = 0;
51   lastNode_ = -999999;
52   howOften_ = 100;
53   decayFactor_ = 0.5;
54   assert(model.solver());
55   int numberColumns = model.solver()->getNumCols();
56   baseSolution_ = new double[numberColumns];
57   memset(baseSolution_, 0, numberColumns * sizeof(double));
58   whereFrom_ = 1 + 8 + 255 * 256;
59   stepSize_ = 0;
60   k_ = 0;
61   kmax_ = 0;
62   nDifferent_ = 0;
63 }
64 
65 // Destructor
~CbcHeuristicVND()66 CbcHeuristicVND::~CbcHeuristicVND()
67 {
68   delete[] baseSolution_;
69 }
70 
71 // Clone
72 CbcHeuristic *
clone() const73 CbcHeuristicVND::clone() const
74 {
75   return new CbcHeuristicVND(*this);
76 }
77 
78 // Assignment operator
79 CbcHeuristicVND &
operator =(const CbcHeuristicVND & rhs)80 CbcHeuristicVND::operator=(const CbcHeuristicVND &rhs)
81 {
82   if (this != &rhs) {
83     CbcHeuristic::operator=(rhs);
84     numberSolutions_ = rhs.numberSolutions_;
85     howOften_ = rhs.howOften_;
86     numberSuccesses_ = rhs.numberSuccesses_;
87     numberTries_ = rhs.numberTries_;
88     lastNode_ = rhs.lastNode_;
89     delete[] baseSolution_;
90     if (model_ && rhs.baseSolution_) {
91       int numberColumns = model_->solver()->getNumCols();
92       baseSolution_ = new double[numberColumns];
93       memcpy(baseSolution_, rhs.baseSolution_, numberColumns * sizeof(double));
94     } else {
95       baseSolution_ = NULL;
96     }
97     stepSize_ = rhs.stepSize_;
98     k_ = rhs.k_;
99     kmax_ = rhs.kmax_;
100     nDifferent_ = rhs.nDifferent_;
101   }
102   return *this;
103 }
104 
105 // Create C++ lines to get to current state
generateCpp(FILE * fp)106 void CbcHeuristicVND::generateCpp(FILE *fp)
107 {
108   CbcHeuristicVND other;
109   fprintf(fp, "0#include \"CbcHeuristicVND.hpp\"\n");
110   fprintf(fp, "3  CbcHeuristicVND heuristicVND(*cbcModel);\n");
111   CbcHeuristic::generateCpp(fp, "heuristicVND");
112   if (howOften_ != other.howOften_)
113     fprintf(fp, "3  heuristicVND.setHowOften(%d);\n", howOften_);
114   else
115     fprintf(fp, "4  heuristicVND.setHowOften(%d);\n", howOften_);
116   fprintf(fp, "3  cbcModel->addHeuristic(&heuristicVND);\n");
117 }
118 
119 // Copy constructor
CbcHeuristicVND(const CbcHeuristicVND & rhs)120 CbcHeuristicVND::CbcHeuristicVND(const CbcHeuristicVND &rhs)
121   : CbcHeuristic(rhs)
122   , numberSolutions_(rhs.numberSolutions_)
123   , howOften_(rhs.howOften_)
124   , numberSuccesses_(rhs.numberSuccesses_)
125   , numberTries_(rhs.numberTries_)
126   , lastNode_(rhs.lastNode_)
127 {
128   if (model_ && rhs.baseSolution_) {
129     int numberColumns = model_->solver()->getNumCols();
130     baseSolution_ = new double[numberColumns];
131     memcpy(baseSolution_, rhs.baseSolution_, numberColumns * sizeof(double));
132   } else {
133     baseSolution_ = NULL;
134   }
135   stepSize_ = rhs.stepSize_;
136   k_ = rhs.k_;
137   kmax_ = rhs.kmax_;
138   nDifferent_ = rhs.nDifferent_;
139 }
140 // Resets stuff if model changes
resetModel(CbcModel *)141 void CbcHeuristicVND::resetModel(CbcModel * /*model*/)
142 {
143   //CbcHeuristic::resetModel(model);
144   delete[] baseSolution_;
145   if (model_ && baseSolution_) {
146     int numberColumns = model_->solver()->getNumCols();
147     baseSolution_ = new double[numberColumns];
148     memset(baseSolution_, 0, numberColumns * sizeof(double));
149   } else {
150     baseSolution_ = NULL;
151   }
152 }
153 /*
154   First tries setting a variable to better value.  If feasible then
155   tries setting others.  If not feasible then tries swaps
156   Returns 1 if solution, 0 if not */
solution(double & solutionValue,double * betterSolution)157 int CbcHeuristicVND::solution(double &solutionValue,
158   double *betterSolution)
159 {
160   numCouldRun_++;
161   int returnCode = 0;
162   const double *bestSolution = model_->bestSolution();
163   if (!bestSolution)
164     return 0; // No solution found yet
165 #ifdef HEURISTIC_INFORM
166   printf("Entering heuristic %s - nRuns %d numCould %d when %d\n",
167     heuristicName(), numRuns_, numCouldRun_, when_);
168 #endif
169   if (numberSolutions_ < model_->getSolutionCount()) {
170     // new solution - add info
171     numberSolutions_ = model_->getSolutionCount();
172 
173     int numberIntegers = model_->numberIntegers();
174     const int *integerVariable = model_->integerVariable();
175 
176     int i;
177     for (i = 0; i < numberIntegers; i++) {
178       int iColumn = integerVariable[i];
179       const OsiObject *object = model_->object(i);
180       // get original bounds
181       double originalLower;
182       double originalUpper;
183       getIntegerInformation(object, originalLower, originalUpper);
184       double value = bestSolution[iColumn];
185       if (value < originalLower) {
186         value = originalLower;
187       } else if (value > originalUpper) {
188         value = originalUpper;
189       }
190     }
191   }
192   int numberNodes = model_->getNodeCount();
193   if (howOften_ == 100) {
194     if (numberNodes < lastNode_ + 12)
195       return 0;
196     // Do at 50 and 100
197     if ((numberNodes > 40 && numberNodes <= 50) || (numberNodes > 90 && numberNodes < 100))
198       numberNodes = howOften_;
199   }
200   if ((numberNodes % howOften_) == 0 && (model_->getCurrentPassNumber() <= 1 || model_->getCurrentPassNumber() == 999999)) {
201     lastNode_ = model_->getNodeCount();
202     OsiSolverInterface *solver = model_->solver();
203 
204     int numberIntegers = model_->numberIntegers();
205     const int *integerVariable = model_->integerVariable();
206 
207     const double *currentSolution = solver->getColSolution();
208     OsiSolverInterface *newSolver = cloneBut(3); // was model_->continuousSolver()->clone();
209     //const double * colLower = newSolver->getColLower();
210     //const double * colUpper = newSolver->getColUpper();
211 
212     double primalTolerance;
213     solver->getDblParam(OsiPrimalTolerance, primalTolerance);
214 
215     // Sort on distance
216     double *distance = new double[numberIntegers];
217     int *which = new int[numberIntegers];
218 
219     int i;
220     int nFix = 0;
221     double tolerance = 10.0 * primalTolerance;
222     for (i = 0; i < numberIntegers; i++) {
223       int iColumn = integerVariable[i];
224       const OsiObject *object = model_->object(i);
225       // get original bounds
226       double originalLower;
227       double originalUpper;
228       getIntegerInformation(object, originalLower, originalUpper);
229       double valueInt = bestSolution[iColumn];
230       if (valueInt < originalLower) {
231         valueInt = originalLower;
232       } else if (valueInt > originalUpper) {
233         valueInt = originalUpper;
234       }
235       baseSolution_[iColumn] = currentSolution[iColumn];
236       distance[i] = fabs(currentSolution[iColumn] - valueInt);
237       which[i] = i;
238       if (fabs(currentSolution[iColumn] - valueInt) < tolerance)
239         nFix++;
240     }
241     CoinSort_2(distance, distance + numberIntegers, which);
242     nDifferent_ = numberIntegers - nFix;
243     stepSize_ = nDifferent_ / 10;
244     k_ = stepSize_;
245     //nFix = numberIntegers-stepSize_;
246     for (i = 0; i < nFix; i++) {
247       int j = which[i];
248       int iColumn = integerVariable[j];
249       const OsiObject *object = model_->object(i);
250       // get original bounds
251       double originalLower;
252       double originalUpper;
253       getIntegerInformation(object, originalLower, originalUpper);
254       double valueInt = bestSolution[iColumn];
255       if (valueInt < originalLower) {
256         valueInt = originalLower;
257       } else if (valueInt > originalUpper) {
258         valueInt = originalUpper;
259       }
260       double nearest = floor(valueInt + 0.5);
261       newSolver->setColLower(iColumn, nearest);
262       newSolver->setColUpper(iColumn, nearest);
263     }
264     delete[] distance;
265     delete[] which;
266     if (nFix > numberIntegers / 5) {
267       //printf("%d integers have samish value\n",nFix);
268       returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue,
269         model_->getCutoff(), "CbcHeuristicVND");
270       if (returnCode < 0)
271         returnCode = 0; // returned on size
272       else
273         numRuns_++;
274       if ((returnCode & 1) != 0)
275         numberSuccesses_++;
276       //printf("return code %d",returnCode);
277       if ((returnCode & 2) != 0) {
278         // could add cut
279         returnCode &= ~2;
280         //printf("could add cut with %d elements (if all 0-1)\n",nFix);
281       } else {
282         //printf("\n");
283       }
284       numberTries_++;
285       if ((numberTries_ % 10) == 0 && numberSuccesses_ * 3 < numberTries_)
286         howOften_ += static_cast< int >(howOften_ * decayFactor_);
287     }
288 
289     delete newSolver;
290   }
291   return returnCode;
292 }
293 // update model
setModel(CbcModel * model)294 void CbcHeuristicVND::setModel(CbcModel *model)
295 {
296   model_ = model;
297   // Get a copy of original matrix
298   assert(model_->solver());
299   delete[] baseSolution_;
300   int numberColumns = model->solver()->getNumCols();
301   baseSolution_ = new double[numberColumns];
302   memset(baseSolution_, 0, numberColumns * sizeof(double));
303 }
304 
305 /* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
306 */
307