source: trunk/Cbc/src/CbcHeuristicDINS.cpp @ 1461

Last change on this file since 1461 was 1432, checked in by bjarni, 10 years ago

Added extra return at end of each source file where needed, to remove possible linefeed conflicts (NightlyBuild? errors)

File size: 15.4 KB
RevLine 
[1368]1// edwin 12/5/09 carved out of CbcHeuristicRINS
2#if defined(_MSC_VER)
3// Turn off compiler warning about long names
4#  pragma warning(disable:4786)
5#endif
6#include <cassert>
7#include <cstdlib>
8#include <cmath>
9#include <cfloat>
10
11#include "OsiSolverInterface.hpp"
12#include "CbcModel.hpp"
13#include "CbcMessage.hpp"
14#include "CbcHeuristicDINS.hpp"
15#include "CbcBranchActual.hpp"
16#include "CbcStrategy.hpp"
17#include "CglPreProcess.hpp"
18
19// Default Constructor
20CbcHeuristicDINS::CbcHeuristicDINS()
21        : CbcHeuristic()
22{
23    numberSolutions_ = 0;
24    numberSuccesses_ = 0;
25    numberTries_ = 0;
26    howOften_ = 100;
27    decayFactor_ = 0.5;
28    maximumKeepSolutions_ = 5;
29    numberKeptSolutions_ = 0;
30    numberIntegers_ = -1;
31    localSpace_ = 10;
32    values_ = NULL;
33}
34
35// Constructor with model - assumed before cuts
36
37CbcHeuristicDINS::CbcHeuristicDINS(CbcModel & model)
38        : CbcHeuristic(model)
39{
40    numberSolutions_ = 0;
41    numberSuccesses_ = 0;
42    numberTries_ = 0;
43    howOften_ = 100;
44    decayFactor_ = 0.5;
45    assert(model.solver());
46    maximumKeepSolutions_ = 5;
47    numberKeptSolutions_ = 0;
48    numberIntegers_ = -1;
49    localSpace_ = 10;
50    values_ = NULL;
51}
52
53// Destructor
54CbcHeuristicDINS::~CbcHeuristicDINS ()
55{
56    for (int i = 0; i < numberKeptSolutions_; i++)
57        delete [] values_[i];
58    delete [] values_;
59}
60
61// Clone
62CbcHeuristic *
63CbcHeuristicDINS::clone() const
64{
65    return new CbcHeuristicDINS(*this);
66}
67
68// Assignment operator
69CbcHeuristicDINS &
70CbcHeuristicDINS::operator=( const CbcHeuristicDINS & rhs)
71{
72    if (this != &rhs) {
73        CbcHeuristic::operator=(rhs);
74        numberSolutions_ = rhs.numberSolutions_;
75        howOften_ = rhs.howOften_;
76        numberSuccesses_ = rhs.numberSuccesses_;
77        numberTries_ = rhs.numberTries_;
78        for (int i = 0; i < numberKeptSolutions_; i++)
79            delete [] values_[i];
80        delete [] values_;
81        maximumKeepSolutions_ = rhs.maximumKeepSolutions_;
82        numberKeptSolutions_ = rhs.numberKeptSolutions_;
83        numberIntegers_ = rhs.numberIntegers_;
84        localSpace_ = rhs.localSpace_;
85        if (model_ && rhs.values_) {
86            assert (numberIntegers_ >= 0);
87            values_ = new int * [maximumKeepSolutions_];
88            for (int i = 0; i < maximumKeepSolutions_; i++)
89                values_[i] = CoinCopyOfArray(rhs.values_[i], numberIntegers_);
90        } else {
91            values_ = NULL;
92        }
93    }
94    return *this;
95}
96
97// Create C++ lines to get to current state
98void
99CbcHeuristicDINS::generateCpp( FILE * fp)
100{
101    CbcHeuristicDINS other;
102    fprintf(fp, "0#include \"CbcHeuristicDINS.hpp\"\n");
103    fprintf(fp, "3  CbcHeuristicDINS heuristicDINS(*cbcModel);\n");
104    CbcHeuristic::generateCpp(fp, "heuristicDINS");
105    if (howOften_ != other.howOften_)
106        fprintf(fp, "3  heuristicDINS.setHowOften(%d);\n", howOften_);
107    else
108        fprintf(fp, "4  heuristicDINS.setHowOften(%d);\n", howOften_);
109    if (maximumKeepSolutions_ != other.maximumKeepSolutions_)
110        fprintf(fp, "3  heuristicDINS.setMaximumKeep(%d);\n", maximumKeepSolutions_);
111    else
112        fprintf(fp, "4  heuristicDINS.setMaximumKeep(%d);\n", maximumKeepSolutions_);
113    fprintf(fp, "3  cbcModel->addHeuristic(&heuristicDINS);\n");
114}
115
116// Copy constructor
117CbcHeuristicDINS::CbcHeuristicDINS(const CbcHeuristicDINS & rhs)
118        :
119        CbcHeuristic(rhs),
120        numberSolutions_(rhs.numberSolutions_),
121        howOften_(rhs.howOften_),
122        numberSuccesses_(rhs.numberSuccesses_),
123        numberTries_(rhs.numberTries_),
124        maximumKeepSolutions_(rhs.maximumKeepSolutions_),
125        numberKeptSolutions_(rhs.numberKeptSolutions_),
126        numberIntegers_(rhs.numberIntegers_),
127        localSpace_(rhs.localSpace_)
128{
129    if (model_ && rhs.values_) {
130        assert (numberIntegers_ >= 0);
131        values_ = new int * [maximumKeepSolutions_];
132        for (int i = 0; i < maximumKeepSolutions_; i++)
133            values_[i] = CoinCopyOfArray(rhs.values_[i], numberIntegers_);
134    } else {
135        values_ = NULL;
136    }
137}
138// Resets stuff if model changes
139void
140CbcHeuristicDINS::resetModel(CbcModel * )
141{
142    //CbcHeuristic::resetModel(model);
143    for (int i = 0; i < numberKeptSolutions_; i++)
144        delete [] values_[i];
145    delete [] values_;
146    numberKeptSolutions_ = 0;
147    numberIntegers_ = -1;
148    numberSolutions_ = 0;
149    values_ = NULL;
150}
151/*
152  First tries setting a variable to better value.  If feasible then
153  tries setting others.  If not feasible then tries swaps
154  Returns 1 if solution, 0 if not */
155int
156CbcHeuristicDINS::solution(double & solutionValue,
157                           double * betterSolution)
158{
159    numCouldRun_++;
160    int returnCode = 0;
161    const double * bestSolution = model_->bestSolution();
162    if (!bestSolution)
163        return 0; // No solution found yet
164    if (numberSolutions_ < model_->getSolutionCount()) {
165        // new solution - add info
166        numberSolutions_ = model_->getSolutionCount();
167
168        int numberIntegers = model_->numberIntegers();
169        const int * integerVariable = model_->integerVariable();
170        if (numberIntegers_ < 0) {
171            numberIntegers_ = numberIntegers;
172            assert (!values_);
173            values_ = new int * [maximumKeepSolutions_];
174            for (int i = 0; i < maximumKeepSolutions_; i++)
175                values_[i] = NULL;
176        } else {
177            assert (numberIntegers == numberIntegers_);
178        }
179        // move solutions (0 will be most recent)
180        {
181            int * temp = values_[maximumKeepSolutions_-1];
182            for (int i = maximumKeepSolutions_ - 1; i > 0; i--)
183                values_[i] = values_[i-1];
184            if (!temp)
185                temp = new int [numberIntegers_];
186            values_[0] = temp;
187        }
188        int i;
189        for (i = 0; i < numberIntegers; i++) {
190            int iColumn = integerVariable[i];
191            double value = bestSolution[iColumn];
192            double nearest = floor(value + 0.5);
193            values_[0][i] = static_cast<int> (nearest);
194        }
195        numberKeptSolutions_ = CoinMin(numberKeptSolutions_ + 1, maximumKeepSolutions_);
196    }
197    int finalReturnCode = 0;
198    if (((model_->getNodeCount() % howOften_) == howOften_ / 2 || !model_->getNodeCount()) && (model_->getCurrentPassNumber() == 1 || model_->getCurrentPassNumber() == 999999)) {
199        OsiSolverInterface * solver = model_->solver();
200
201        int numberIntegers = model_->numberIntegers();
202        const int * integerVariable = model_->integerVariable();
203
204        const double * currentSolution = solver->getColSolution();
205        int localSpace = localSpace_;
206        // 0 means finished but no solution, 1 solution, 2 node limit
207        int status = -1;
208        double cutoff = model_->getCutoff();
209        while (status) {
210            status = 0;
211            OsiSolverInterface * newSolver = cloneBut(3); // was model_->continuousSolver()->clone();
212            const double * colLower = solver->getColLower();
213            const double * colUpper = solver->getColUpper();
214
215            double primalTolerance;
216            solver->getDblParam(OsiPrimalTolerance, primalTolerance);
217            const double * continuousSolution = newSolver->getColSolution();
218            // Space for added constraint
219            double * element = new double [numberIntegers];
220            int * column = new int [numberIntegers];
221            int i;
222            int nFix = 0;
223            int nCouldFix = 0;
224            int nCouldFix2 = 0;
225            int nBound = 0;
226            int nEl = 0;
227            double bias = localSpace;
228            int okSame = numberKeptSolutions_ - 1;
229            for (i = 0; i < numberIntegers; i++) {
230                int iColumn = integerVariable[i];
231                const OsiObject * object = model_->object(i);
232                // get original bounds
233                double originalLower;
234                double originalUpper;
235                getIntegerInformation( object, originalLower, originalUpper);
236                double valueInt = bestSolution[iColumn];
237                if (valueInt < originalLower) {
238                    valueInt = originalLower;
239                } else if (valueInt > originalUpper) {
240                    valueInt = originalUpper;
241                }
242                int intValue = static_cast<int> (floor(valueInt + 0.5));
243                double currentValue = currentSolution[iColumn];
244                double currentLower = colLower[iColumn];
245                double currentUpper = colUpper[iColumn];
246                if (fabs(valueInt - currentValue) >= 0.5) {
247                    // Re-bound
248                    nBound++;
249                    if (intValue >= currentValue) {
250                        currentLower = CoinMax(currentLower, ceil(2 * currentValue - intValue));
251                        currentUpper = intValue;
252                    } else {
253                        currentLower = intValue;
254                        currentUpper = CoinMin(currentUpper, floor(2 * currentValue - intValue));
255                    }
256                    newSolver->setColLower(iColumn, currentLower);
257                    newSolver->setColUpper(iColumn, currentUpper);
258                } else {
259                    // See if can fix
260                    bool canFix = false;
261                    double continuousValue = continuousSolution[iColumn];
262                    if (fabs(currentValue - valueInt) < 10.0*primalTolerance) {
263                        if (currentUpper - currentLower > 1.0) {
264                            // General integer variable
265                            canFix = true;
266                        } else if (fabs(continuousValue - valueInt) < 10.0*primalTolerance) {
267                            int nSame = 1;
268                            //assert (intValue==values_[0][i]);
269                            for (int k = 1; k < numberKeptSolutions_; k++) {
270                                if (intValue == values_[k][i])
271                                    nSame++;
272                            }
273                            if (nSame >= okSame) {
274                                // can fix
275                                canFix = true;
276                            } else {
277                                nCouldFix++;
278                            }
279                        } else {
280                            nCouldFix2++;
281                        }
282                    }
283                    if (canFix) {
284                        newSolver->setColLower(iColumn, intValue);
285                        newSolver->setColUpper(iColumn, intValue);
286                        nFix++;
287                    } else {
288                        if (currentUpper - currentLower > 1.0) {
289                            // General integer variable
290                            currentLower = floor(currentValue);
291                            if (intValue >= currentLower && intValue <= currentLower + 1) {
292                                newSolver->setColLower(iColumn, currentLower);
293                                newSolver->setColUpper(iColumn, currentLower + 1.0);
294                            } else {
295                                // fix
296                                double value;
297                                if (intValue < currentLower)
298                                    value = currentLower;
299                                else
300                                    value = currentLower + 1;
301                                newSolver->setColLower(iColumn, value);
302                                newSolver->setColUpper(iColumn, value);
303                                nFix++;
304                            }
305                        } else {
306                            // 0-1 (ish)
307                            column[nEl] = iColumn;
308                            if (intValue == currentLower) {
309                                bias += currentLower;
310                                element[nEl++] = 1.0;
311                            } else if (intValue == currentUpper) {
312                                bias += currentUpper;
313                                element[nEl++] = -1.0;
314                            } else {
315                                printf("bad DINS logic\n");
316                                abort();
317                            }
318                        }
319                    }
320                }
321            }
322            char generalPrint[200];
323            sprintf(generalPrint,
324                    "%d fixed, %d same as cont/int, %d same as int - %d bounded %d in cut\n",
325                    nFix, nCouldFix, nCouldFix2, nBound, nEl);
326            model_->messageHandler()->message(CBC_FPUMP2, model_->messages())
327            << generalPrint
328            << CoinMessageEol;
329            if (nFix > numberIntegers / 10) {
[1393]330#ifdef JJF_ZERO
[1368]331                newSolver->initialSolve();
332                printf("obj %g\n", newSolver->getObjValue());
333                for (i = 0; i < numberIntegers; i++) {
334                    int iColumn = integerVariable[i];
335                    printf("%d new bounds %g %g - solutions %g %g\n",
336                           iColumn, newSolver->getColLower()[iColumn],
337                           newSolver->getColUpper()[iColumn],
338                           bestSolution[iColumn],
339                           currentSolution[iColumn]);
340                }
341#endif
342                if (nEl > 0)
343                    newSolver->addRow(nEl, column, element, -COIN_DBL_MAX, bias);
344                //printf("%d integers have same value\n",nFix);
345                returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue,
346                                                 cutoff, "CbcHeuristicDINS");
347                if (returnCode < 0) {
348                    returnCode = 0; // returned on size
349                    status = 0;
350                } else {
351                    numRuns_++;
352                    if ((returnCode&2) != 0) {
353                        // could add cut as complete search
354                        returnCode &= ~2;
355                        if ((returnCode&1) != 0) {
356                            numberSuccesses_++;
357                            status = 1;
358                        } else {
359                            // no solution
360                            status = 0;
361                        }
362                    } else {
363                        if ((returnCode&1) != 0) {
364                            numberSuccesses_++;
365                            status = 1;
366                        } else {
367                            // no solution but node limit
368                            status = 2;
369                            if (nEl)
370                                localSpace -= 5;
371                            else
372                                localSpace = -1;
373                            if (localSpace < 0)
374                                status = 0;
375                        }
376                    }
377                    if ((returnCode&1) != 0) {
378                        cutoff = CoinMin(cutoff, solutionValue - model_->getCutoffIncrement());
379                        finalReturnCode = 1;
380                    }
381                }
382            }
383            delete [] element;
384            delete [] column;
385            delete newSolver;
386        }
387        numberTries_++;
388        if ((numberTries_ % 10) == 0 && numberSuccesses_*3 < numberTries_)
389            howOften_ += static_cast<int> (howOften_ * decayFactor_);
390    }
391    return finalReturnCode;
392}
393// update model
394void CbcHeuristicDINS::setModel(CbcModel * model)
395{
396    model_ = model;
397    // Get a copy of original matrix
398    assert(model_->solver());
399    for (int i = 0; i < numberKeptSolutions_; i++)
400        delete [] values_[i];
401    delete [] values_;
402    numberKeptSolutions_ = 0;
403    numberIntegers_ = -1;
404    numberSolutions_ = 0;
405    values_ = NULL;
406}
[1432]407
Note: See TracBrowser for help on using the repository browser.