source: stable/2.5/Cbc/src/CbcHeuristicRENS.cpp @ 1510

Last change on this file since 1510 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: 10.2 KB
Line 
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 "CbcHeuristicRENS.hpp"
15#include "CbcBranchActual.hpp"
16#include "CbcStrategy.hpp"
17#include "CglPreProcess.hpp"
18
19// Default Constructor
20CbcHeuristicRENS::CbcHeuristicRENS()
21        : CbcHeuristic()
22{
23    numberTries_ = 0;
24    whereFrom_ = 256 + 1;
25}
26
27// Constructor with model - assumed before cuts
28
29CbcHeuristicRENS::CbcHeuristicRENS(CbcModel & model)
30        : CbcHeuristic(model)
31{
32    numberTries_ = 0;
33    whereFrom_ = 256 + 1;
34}
35
36// Destructor
37CbcHeuristicRENS::~CbcHeuristicRENS ()
38{
39}
40
41// Clone
42CbcHeuristic *
43CbcHeuristicRENS::clone() const
44{
45    return new CbcHeuristicRENS(*this);
46}
47
48// Assignment operator
49CbcHeuristicRENS &
50CbcHeuristicRENS::operator=( const CbcHeuristicRENS & rhs)
51{
52    if (this != &rhs) {
53        CbcHeuristic::operator=(rhs);
54        numberTries_ = rhs.numberTries_;
55    }
56    return *this;
57}
58
59// Copy constructor
60CbcHeuristicRENS::CbcHeuristicRENS(const CbcHeuristicRENS & rhs)
61        :
62        CbcHeuristic(rhs),
63        numberTries_(rhs.numberTries_)
64{
65}
66// Resets stuff if model changes
67void
68CbcHeuristicRENS::resetModel(CbcModel * )
69{
70}
71int
72CbcHeuristicRENS::solution(double & solutionValue,
73                           double * betterSolution)
74{
75    int returnCode = 0;
76    const double * bestSolution = model_->bestSolution();
77    if (numberTries_ || (when() < 2 && bestSolution))
78        return 0;
79    numberTries_++;
80    OsiSolverInterface * solver = model_->solver();
81
82    int numberIntegers = model_->numberIntegers();
83    const int * integerVariable = model_->integerVariable();
84
85    const double * currentSolution = solver->getColSolution();
86    OsiSolverInterface * newSolver = cloneBut(3); // was model_->continuousSolver()->clone();
87    const double * colLower = newSolver->getColLower();
88    const double * colUpper = newSolver->getColUpper();
89
90    double primalTolerance;
91    solver->getDblParam(OsiPrimalTolerance, primalTolerance);
92
93    int i;
94    int numberFixed = 0;
95    int numberTightened = 0;
96    int numberAtBound = 0;
97    int numberColumns = newSolver->getNumCols();
98    int numberContinuous = numberColumns - numberIntegers;
99
100    for (i = 0; i < numberIntegers; i++) {
101        int iColumn = integerVariable[i];
102        double value = currentSolution[iColumn];
103        double lower = colLower[iColumn];
104        double upper = colUpper[iColumn];
105        value = CoinMax(value, lower);
106        value = CoinMin(value, upper);
107#define RENS_FIX_ONLY_LOWER
108#ifndef RENS_FIX_ONLY_LOWER
109        if (fabs(value - floor(value + 0.5)) < 1.0e-8) {
110            value = floor(value + 0.5);
111            if (value == lower || value == upper)
112                numberAtBound++;
113            newSolver->setColLower(iColumn, value);
114            newSolver->setColUpper(iColumn, value);
115            numberFixed++;
116        } else if (colUpper[iColumn] - colLower[iColumn] >= 2.0) {
117            numberTightened++;
118            newSolver->setColLower(iColumn, floor(value));
119            newSolver->setColUpper(iColumn, ceil(value));
120        }
121#else
122        if (fabs(value - floor(value + 0.5)) < 1.0e-8 &&
123                floor(value + 0.5) == lower) {
124            value = floor(value + 0.5);
125            numberAtBound++;
126            newSolver->setColLower(iColumn, value);
127            newSolver->setColUpper(iColumn, value);
128            numberFixed++;
129        } else if (colUpper[iColumn] - colLower[iColumn] >= 2.0) {
130            numberTightened++;
131            if (fabs(value - floor(value + 0.5)) < 1.0e-8) {
132                value = floor(value + 0.5);
133                if (value < upper) {
134                    newSolver->setColLower(iColumn, CoinMax(value - 1.0, lower));
135                    newSolver->setColUpper(iColumn, CoinMin(value + 1.0, upper));
136                } else {
137                    newSolver->setColLower(iColumn, upper - 1.0);
138                }
139            } else {
140                newSolver->setColLower(iColumn, floor(value));
141                newSolver->setColUpper(iColumn, ceil(value));
142            }
143        }
144#endif
145    }
146    if (numberFixed > numberIntegers / 5) {
147        if (numberContinuous > numberIntegers && numberFixed < numberColumns / 5) {
148#define RENS_FIX_CONTINUOUS
149#ifdef RENS_FIX_CONTINUOUS
150            const double * colLower = newSolver->getColLower();
151            //const double * colUpper = newSolver->getColUpper();
152            int nAtLb = 0;
153            double sumDj = 0.0;
154            const double * dj = newSolver->getReducedCost();
155            double direction = newSolver->getObjSense();
156            for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
157                if (!newSolver->isInteger(iColumn)) {
158                    double value = currentSolution[iColumn];
159                    if (value < colLower[iColumn] + 1.0e-8) {
160                        double djValue = dj[iColumn] * direction;
161                        nAtLb++;
162                        sumDj += djValue;
163                    }
164                }
165            }
166            if (nAtLb) {
167                // fix some continuous
168                double * sort = new double[nAtLb];
169                int * which = new int [nAtLb];
170                double threshold = CoinMax((0.01 * sumDj) / static_cast<double>(nAtLb), 1.0e-6);
171                int nFix2 = 0;
172                for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
173                    if (!newSolver->isInteger(iColumn)) {
174                        double value = currentSolution[iColumn];
175                        if (value < colLower[iColumn] + 1.0e-8) {
176                            double djValue = dj[iColumn] * direction;
177                            if (djValue > threshold) {
178                                sort[nFix2] = -djValue;
179                                which[nFix2++] = iColumn;
180                            }
181                        }
182                    }
183                }
184                CoinSort_2(sort, sort + nFix2, which);
185                nFix2 = CoinMin(nFix2, (numberColumns - numberFixed) / 2);
186                for (int i = 0; i < nFix2; i++) {
187                    int iColumn = which[i];
188                    newSolver->setColUpper(iColumn, colLower[iColumn]);
189                }
190                delete [] sort;
191                delete [] which;
192#ifdef CLP_INVESTIGATE2
193                printf("%d integers fixed (%d tightened) (%d at bound), and %d continuous fixed at lb\n",
194                       numberFixed, numberTightened, numberAtBound, nFix2);
195#endif
196            }
197#endif
198        }
199#ifdef COIN_DEVELOP
200        printf("%d integers fixed and %d tightened\n", numberFixed, numberTightened);
201#endif
202        returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue,
203                                         model_->getCutoff(), "CbcHeuristicRENS");
204        if (returnCode < 0) {
205            returnCode = 0; // returned on size
206#ifdef RENS_FIX_CONTINUOUS
207            if (numberContinuous > numberIntegers && numberFixed >= numberColumns / 5) {
208                const double * colLower = newSolver->getColLower();
209                //const double * colUpper = newSolver->getColUpper();
210                int nAtLb = 0;
211                double sumDj = 0.0;
212                const double * dj = newSolver->getReducedCost();
213                double direction = newSolver->getObjSense();
214                for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
215                    if (!newSolver->isInteger(iColumn)) {
216                        double value = currentSolution[iColumn];
217                        if (value < colLower[iColumn] + 1.0e-8) {
218                            double djValue = dj[iColumn] * direction;
219                            nAtLb++;
220                            sumDj += djValue;
221                        }
222                    }
223                }
224                if (nAtLb) {
225                    // fix some continuous
226                    double * sort = new double[nAtLb];
227                    int * which = new int [nAtLb];
228                    double threshold = CoinMax((0.01 * sumDj) / static_cast<double>(nAtLb), 1.0e-6);
229                    int nFix2 = 0;
230                    for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
231                        if (!newSolver->isInteger(iColumn)) {
232                            double value = currentSolution[iColumn];
233                            if (value < colLower[iColumn] + 1.0e-8) {
234                                double djValue = dj[iColumn] * direction;
235                                if (djValue > threshold) {
236                                    sort[nFix2] = -djValue;
237                                    which[nFix2++] = iColumn;
238                                }
239                            }
240                        }
241                    }
242                    CoinSort_2(sort, sort + nFix2, which);
243                    nFix2 = CoinMin(nFix2, (numberColumns - numberFixed) / 2);
244                    for (int i = 0; i < nFix2; i++) {
245                        int iColumn = which[i];
246                        newSolver->setColUpper(iColumn, colLower[iColumn]);
247                    }
248                    delete [] sort;
249                    delete [] which;
250#ifdef CLP_INVESTIGATE2
251                    printf("%d integers fixed (%d tightened) (%d at bound), and %d continuous fixed at lb\n",
252                           numberFixed, numberTightened, numberAtBound, nFix2);
253#endif
254                }
255                returnCode = smallBranchAndBound(newSolver, numberNodes_, betterSolution, solutionValue,
256                                                 model_->getCutoff(), "CbcHeuristicRENS");
257#endif
258            }
259        }
260        //printf("return code %d",returnCode);
261        if ((returnCode&2) != 0) {
262            // could add cut
263            returnCode &= ~2;
264#ifdef COIN_DEVELOP
265            if (!numberTightened && numberFixed == numberAtBound)
266                printf("could add cut with %d elements\n", numberFixed);
267#endif
268        } else {
269            //printf("\n");
270        }
271    }
272
273    delete newSolver;
274    return returnCode;
275}
276// update model
277void CbcHeuristicRENS::setModel(CbcModel * model)
278{
279    model_ = model;
280}
281
Note: See TracBrowser for help on using the repository browser.