source: trunk/Cbc/src/CbcSolverHeuristics.cpp @ 1955

Last change on this file since 1955 was 1955, checked in by forrest, 8 years ago

add more twiddles for advanced use of heuristics

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 76.9 KB
Line 
1/* $Id: CbcSolverHeuristics.cpp 1955 2013-08-12 19:33:18Z forrest $ */
2// Copyright (C) 2007, 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
7/*! \file CbcSolverHeuristics.cpp
8    \brief Second level routines for the cbc stand-alone solver.
9*/
10
11#include "CbcConfig.h"
12#include "CoinPragma.hpp"
13
14
15#include "CoinTime.hpp"
16
17#include "OsiClpSolverInterface.hpp"
18
19#include "ClpPresolve.hpp"
20
21#include "CbcOrClpParam.hpp"
22
23#include "CbcModel.hpp"
24
25#include "CbcHeuristicLocal.hpp"
26#include "CbcHeuristicPivotAndFix.hpp"
27//#include "CbcHeuristicPivotAndComplement.hpp"
28#include "CbcHeuristicRandRound.hpp"
29#include "CbcHeuristicGreedy.hpp"
30#include "CbcHeuristicFPump.hpp"
31#include "CbcHeuristicRINS.hpp"
32#include "CbcHeuristicDW.hpp"
33
34#include "CbcHeuristicDiveCoefficient.hpp"
35#include "CbcHeuristicDiveFractional.hpp"
36#include "CbcHeuristicDiveGuided.hpp"
37#include "CbcHeuristicDiveVectorLength.hpp"
38#include "CbcHeuristicDivePseudoCost.hpp"
39#include "CbcHeuristicDiveLineSearch.hpp"
40
41#include "CbcStrategy.hpp"
42#include "OsiAuxInfo.hpp"
43
44#include "ClpSimplexOther.hpp"
45
46// Crunch down model
47void
48crunchIt(ClpSimplex * model)
49{
50#ifdef JJF_ZERO
51    model->dual();
52#else
53    int numberColumns = model->numberColumns();
54    int numberRows = model->numberRows();
55    // Use dual region
56    double * rhs = model->dualRowSolution();
57    int * whichRow = new int[3*numberRows];
58    int * whichColumn = new int[2*numberColumns];
59    int nBound;
60    ClpSimplex * small = static_cast<ClpSimplexOther *> (model)->crunch(rhs, whichRow, whichColumn,
61                         nBound, false, false);
62    if (small) {
63        small->dual();
64        if (small->problemStatus() == 0) {
65            model->setProblemStatus(0);
66            static_cast<ClpSimplexOther *> (model)->afterCrunch(*small, whichRow, whichColumn, nBound);
67        } else if (small->problemStatus() != 3) {
68            model->setProblemStatus(1);
69        } else {
70            if (small->problemStatus() == 3) {
71                // may be problems
72                small->computeObjectiveValue();
73                model->setObjectiveValue(small->objectiveValue());
74                model->setProblemStatus(3);
75            } else {
76                model->setProblemStatus(3);
77            }
78        }
79        delete small;
80    } else {
81        model->setProblemStatus(1);
82    }
83    delete [] whichRow;
84    delete [] whichColumn;
85#endif
86}
87/*
88  On input
89  doAction - 0 just fix in original and return NULL
90             1 return fixed non-presolved solver
91             2 as one but use presolve Inside this
92             3 use presolve and fix ones with large cost
93             ? do heuristics and set best solution
94             ? do BAB and just set best solution
95             10+ then use lastSolution and relax a few
96             -2 cleanup afterwards if using 2
97  On output - number fixed
98*/
99OsiClpSolverInterface *
100fixVubs(CbcModel & model, int skipZero2,
101        int & doAction,
102        CoinMessageHandler * /*generalMessageHandler*/,
103        const double * lastSolution, double dextra[6],
104        int extra[5])
105{
106    if (doAction == 11 && !lastSolution)
107        lastSolution = model.bestSolution();
108    assert (((doAction >= 0 && doAction <= 3) && !lastSolution) || (doAction == 11 && lastSolution));
109    double fractionIntFixed = dextra[3];
110    double fractionFixed = dextra[4];
111    double fixAbove = dextra[2];
112    double fixAboveValue = (dextra[5] > 0.0) ? dextra[5] : 1.0;
113#ifdef COIN_DETAIL
114    double time1 = CoinCpuTime();
115#endif
116    int leaveIntFree = extra[1];
117    OsiSolverInterface * originalSolver = model.solver();
118    OsiClpSolverInterface * originalClpSolver = dynamic_cast< OsiClpSolverInterface*> (originalSolver);
119    ClpSimplex * originalLpSolver = originalClpSolver->getModelPtr();
120    int * originalColumns = NULL;
121    OsiClpSolverInterface * clpSolver;
122    ClpSimplex * lpSolver;
123    ClpPresolve pinfo;
124    assert(originalSolver->getObjSense() > 0);
125    if (doAction == 2 || doAction == 3) {
126        double * saveLB = NULL;
127        double * saveUB = NULL;
128        int numberColumns = originalLpSolver->numberColumns();
129        if (fixAbove > 0.0) {
130#ifdef COIN_DETAIL
131            double time1 = CoinCpuTime();
132#endif
133            originalClpSolver->initialSolve();
134            COIN_DETAIL_PRINT(printf("first solve took %g seconds\n", CoinCpuTime() - time1));
135            double * columnLower = originalLpSolver->columnLower() ;
136            double * columnUpper = originalLpSolver->columnUpper() ;
137            const double * solution = originalLpSolver->primalColumnSolution();
138            saveLB = CoinCopyOfArray(columnLower, numberColumns);
139            saveUB = CoinCopyOfArray(columnUpper, numberColumns);
140            const double * objective = originalLpSolver->getObjCoefficients() ;
141            int iColumn;
142            int nFix = 0;
143            int nArt = 0;
144            for (iColumn = 0; iColumn < numberColumns; iColumn++) {
145                if (objective[iColumn] > fixAbove) {
146                    if (solution[iColumn] < columnLower[iColumn] + 1.0e-8) {
147                        columnUpper[iColumn] = columnLower[iColumn];
148                        nFix++;
149                    } else {
150                        nArt++;
151                    }
152                } else if (objective[iColumn] < -fixAbove) {
153                    if (solution[iColumn] > columnUpper[iColumn] - 1.0e-8) {
154                        columnLower[iColumn] = columnUpper[iColumn];
155                        nFix++;
156                    } else {
157                        nArt++;
158                    }
159                }
160            }
161            COIN_DETAIL_PRINT(printf("%d artificials fixed, %d left as in solution\n", nFix, nArt));
162            lpSolver = pinfo.presolvedModel(*originalLpSolver, 1.0e-8, true, 10);
163            if (!lpSolver || doAction == 2) {
164                // take off fixing in original
165                memcpy(columnLower, saveLB, numberColumns*sizeof(double));
166                memcpy(columnUpper, saveUB, numberColumns*sizeof(double));
167            }
168            delete [] saveLB;
169            delete [] saveUB;
170            if (!lpSolver) {
171                // try again
172                pinfo.destroyPresolve();
173                lpSolver = pinfo.presolvedModel(*originalLpSolver, 1.0e-8, true, 10);
174                assert (lpSolver);
175            }
176        } else {
177            lpSolver = pinfo.presolvedModel(*originalLpSolver, 1.0e-8, true, 10);
178            assert (lpSolver);
179        }
180        clpSolver = new OsiClpSolverInterface(lpSolver, true);
181        assert(lpSolver == clpSolver->getModelPtr());
182        numberColumns = lpSolver->numberColumns();
183        originalColumns = CoinCopyOfArray(pinfo.originalColumns(), numberColumns);
184        doAction = 1;
185    } else {
186        OsiSolverInterface * solver = originalSolver->clone();
187        clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
188        lpSolver = clpSolver->getModelPtr();
189    }
190    // Tighten bounds
191    lpSolver->tightenPrimalBounds(0.0, 11, true);
192    int numberColumns = clpSolver->getNumCols() ;
193    double * saveColumnLower = CoinCopyOfArray(lpSolver->columnLower(), numberColumns);
194    double * saveColumnUpper = CoinCopyOfArray(lpSolver->columnUpper(), numberColumns);
195    //char generalPrint[200];
196    const double *objective = lpSolver->getObjCoefficients() ;
197    double *columnLower = lpSolver->columnLower() ;
198    double *columnUpper = lpSolver->columnUpper() ;
199    int numberRows = clpSolver->getNumRows();
200    int iRow, iColumn;
201
202    // Row copy
203    CoinPackedMatrix matrixByRow(*clpSolver->getMatrixByRow());
204    const double * elementByRow = matrixByRow.getElements();
205    const int * column = matrixByRow.getIndices();
206    const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
207    const int * rowLength = matrixByRow.getVectorLengths();
208
209    // Column copy
210    CoinPackedMatrix  matrixByCol(*clpSolver->getMatrixByCol());
211    //const double * element = matrixByCol.getElements();
212    const int * row = matrixByCol.getIndices();
213    const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
214    const int * columnLength = matrixByCol.getVectorLengths();
215
216    const double * rowLower = clpSolver->getRowLower();
217    const double * rowUpper = clpSolver->getRowUpper();
218
219    // Get maximum size of VUB tree
220    // otherColumn is one fixed to 0 if this one zero
221    int nEl = matrixByCol.getNumElements();
222    CoinBigIndex * fixColumn = new CoinBigIndex [numberColumns+1];
223    int * otherColumn = new int [nEl];
224    int * fix = new int[numberColumns];
225    char * mark = new char [numberColumns];
226    memset(mark, 0, numberColumns);
227    int numberInteger = 0;
228    int numberOther = 0;
229    fixColumn[0] = 0;
230    double large = lpSolver->largeValue(); // treat bounds > this as infinite
231#ifndef NDEBUG
232    double large2 = 1.0e10 * large;
233#endif
234    double tolerance = lpSolver->primalTolerance();
235    int * check = new int[numberRows];
236    for (iRow = 0; iRow < numberRows; iRow++) {
237        check[iRow] = -2; // don't check
238        if (rowLower[iRow] < -1.0e6 && rowUpper[iRow] > 1.0e6)
239            continue;// unlikely
240        // possible row
241        int numberPositive = 0;
242        int iPositive = -1;
243        int numberNegative = 0;
244        int iNegative = -1;
245        CoinBigIndex rStart = rowStart[iRow];
246        CoinBigIndex rEnd = rowStart[iRow] + rowLength[iRow];
247        CoinBigIndex j;
248        int kColumn;
249        for (j = rStart; j < rEnd; ++j) {
250            double value = elementByRow[j];
251            kColumn = column[j];
252            if (columnUpper[kColumn] > columnLower[kColumn]) {
253                if (value > 0.0) {
254                    numberPositive++;
255                    iPositive = kColumn;
256                } else {
257                    numberNegative++;
258                    iNegative = kColumn;
259                }
260            }
261        }
262        if (numberPositive == 1 && numberNegative == 1)
263            check[iRow] = -1; // try both
264        if (numberPositive == 1 && rowLower[iRow] > -1.0e20)
265            check[iRow] = iPositive;
266        else if (numberNegative == 1 && rowUpper[iRow] < 1.0e20)
267            check[iRow] = iNegative;
268    }
269    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
270        fix[iColumn] = -1;
271        if (columnUpper[iColumn] > columnLower[iColumn] + 1.0e-8) {
272            if (clpSolver->isInteger(iColumn))
273                numberInteger++;
274            if (columnLower[iColumn] == 0.0) {
275                bool infeasible = false;
276                fix[iColumn] = 0;
277                // fake upper bound
278                double saveUpper = columnUpper[iColumn];
279                columnUpper[iColumn] = 0.0;
280                for (CoinBigIndex i = columnStart[iColumn];
281                        i < columnStart[iColumn] + columnLength[iColumn]; i++) {
282                    iRow = row[i];
283                    if (check[iRow] != -1 && check[iRow] != iColumn)
284                        continue; // unlikely
285                    // possible row
286                    int infiniteUpper = 0;
287                    int infiniteLower = 0;
288                    double maximumUp = 0.0;
289                    double maximumDown = 0.0;
290                    double newBound;
291                    CoinBigIndex rStart = rowStart[iRow];
292                    CoinBigIndex rEnd = rowStart[iRow] + rowLength[iRow];
293                    CoinBigIndex j;
294                    int kColumn;
295                    // Compute possible lower and upper ranges
296                    for (j = rStart; j < rEnd; ++j) {
297                        double value = elementByRow[j];
298                        kColumn = column[j];
299                        if (value > 0.0) {
300                            if (columnUpper[kColumn] >= large) {
301                                ++infiniteUpper;
302                            } else {
303                                maximumUp += columnUpper[kColumn] * value;
304                            }
305                            if (columnLower[kColumn] <= -large) {
306                                ++infiniteLower;
307                            } else {
308                                maximumDown += columnLower[kColumn] * value;
309                            }
310                        } else if (value < 0.0) {
311                            if (columnUpper[kColumn] >= large) {
312                                ++infiniteLower;
313                            } else {
314                                maximumDown += columnUpper[kColumn] * value;
315                            }
316                            if (columnLower[kColumn] <= -large) {
317                                ++infiniteUpper;
318                            } else {
319                                maximumUp += columnLower[kColumn] * value;
320                            }
321                        }
322                    }
323                    // Build in a margin of error
324                    maximumUp += 1.0e-8 * fabs(maximumUp);
325                    maximumDown -= 1.0e-8 * fabs(maximumDown);
326                    double maxUp = maximumUp + infiniteUpper * 1.0e31;
327                    double maxDown = maximumDown - infiniteLower * 1.0e31;
328                    if (maxUp <= rowUpper[iRow] + tolerance &&
329                            maxDown >= rowLower[iRow] - tolerance) {
330                        //printf("Redundant row in vubs %d\n",iRow);
331                    } else {
332                        if (maxUp < rowLower[iRow] - 100.0*tolerance ||
333                                maxDown > rowUpper[iRow] + 100.0*tolerance) {
334                            infeasible = true;
335                            break;
336                        }
337                        double lower = rowLower[iRow];
338                        double upper = rowUpper[iRow];
339                        for (j = rStart; j < rEnd; ++j) {
340                            double value = elementByRow[j];
341                            kColumn = column[j];
342                            double nowLower = columnLower[kColumn];
343                            double nowUpper = columnUpper[kColumn];
344                            if (value > 0.0) {
345                                // positive value
346                                if (lower > -large) {
347                                    if (!infiniteUpper) {
348                                        assert(nowUpper < large2);
349                                        newBound = nowUpper +
350                                                   (lower - maximumUp) / value;
351                                        // relax if original was large
352                                        if (fabs(maximumUp) > 1.0e8)
353                                            newBound -= 1.0e-12 * fabs(maximumUp);
354                                    } else if (infiniteUpper == 1 && nowUpper > large) {
355                                        newBound = (lower - maximumUp) / value;
356                                        // relax if original was large
357                                        if (fabs(maximumUp) > 1.0e8)
358                                            newBound -= 1.0e-12 * fabs(maximumUp);
359                                    } else {
360                                        newBound = -COIN_DBL_MAX;
361                                    }
362                                    if (newBound > nowLower + 1.0e-12 && newBound > -large) {
363                                        // Tighten the lower bound
364                                        // check infeasible (relaxed)
365                                        if (nowUpper < newBound) {
366                                            if (nowUpper - newBound <
367                                                    -100.0*tolerance) {
368                                                infeasible = true;
369                                                break;
370                                            }
371                                        }
372                                    }
373                                }
374                                if (upper < large) {
375                                    if (!infiniteLower) {
376                                        assert(nowLower > - large2);
377                                        newBound = nowLower +
378                                                   (upper - maximumDown) / value;
379                                        // relax if original was large
380                                        if (fabs(maximumDown) > 1.0e8)
381                                            newBound += 1.0e-12 * fabs(maximumDown);
382                                    } else if (infiniteLower == 1 && nowLower < -large) {
383                                        newBound =   (upper - maximumDown) / value;
384                                        // relax if original was large
385                                        if (fabs(maximumDown) > 1.0e8)
386                                            newBound += 1.0e-12 * fabs(maximumDown);
387                                    } else {
388                                        newBound = COIN_DBL_MAX;
389                                    }
390                                    if (newBound < nowUpper - 1.0e-12 && newBound < large) {
391                                        // Tighten the upper bound
392                                        // check infeasible (relaxed)
393                                        if (nowLower > newBound) {
394                                            if (newBound - nowLower <
395                                                    -100.0*tolerance) {
396                                                infeasible = true;
397                                                break;
398                                            } else {
399                                                newBound = nowLower;
400                                            }
401                                        }
402                                        if (!newBound || (clpSolver->isInteger(kColumn) && newBound < 0.999)) {
403                                            // fix to zero
404                                            if (!mark[kColumn]) {
405                                                otherColumn[numberOther++] = kColumn;
406                                                mark[kColumn] = 1;
407                                                if (check[iRow] == -1)
408                                                    check[iRow] = iColumn;
409                                                else
410                                                    assert(check[iRow] == iColumn);
411                                            }
412                                        }
413                                    }
414                                }
415                            } else {
416                                // negative value
417                                if (lower > -large) {
418                                    if (!infiniteUpper) {
419                                        assert(nowLower < large2);
420                                        newBound = nowLower +
421                                                   (lower - maximumUp) / value;
422                                        // relax if original was large
423                                        if (fabs(maximumUp) > 1.0e8)
424                                            newBound += 1.0e-12 * fabs(maximumUp);
425                                    } else if (infiniteUpper == 1 && nowLower < -large) {
426                                        newBound = (lower - maximumUp) / value;
427                                        // relax if original was large
428                                        if (fabs(maximumUp) > 1.0e8)
429                                            newBound += 1.0e-12 * fabs(maximumUp);
430                                    } else {
431                                        newBound = COIN_DBL_MAX;
432                                    }
433                                    if (newBound < nowUpper - 1.0e-12 && newBound < large) {
434                                        // Tighten the upper bound
435                                        // check infeasible (relaxed)
436                                        if (nowLower > newBound) {
437                                            if (newBound - nowLower <
438                                                    -100.0*tolerance) {
439                                                infeasible = true;
440                                                break;
441                                            } else {
442                                                newBound = nowLower;
443                                            }
444                                        }
445                                        if (!newBound || (clpSolver->isInteger(kColumn) && newBound < 0.999)) {
446                                            // fix to zero
447                                            if (!mark[kColumn]) {
448                                                otherColumn[numberOther++] = kColumn;
449                                                mark[kColumn] = 1;
450                                                if (check[iRow] == -1)
451                                                    check[iRow] = iColumn;
452                                                else
453                                                    assert(check[iRow] == iColumn);
454                                            }
455                                        }
456                                    }
457                                }
458                                if (upper < large) {
459                                    if (!infiniteLower) {
460                                        assert(nowUpper < large2);
461                                        newBound = nowUpper +
462                                                   (upper - maximumDown) / value;
463                                        // relax if original was large
464                                        if (fabs(maximumDown) > 1.0e8)
465                                            newBound -= 1.0e-12 * fabs(maximumDown);
466                                    } else if (infiniteLower == 1 && nowUpper > large) {
467                                        newBound =   (upper - maximumDown) / value;
468                                        // relax if original was large
469                                        if (fabs(maximumDown) > 1.0e8)
470                                            newBound -= 1.0e-12 * fabs(maximumDown);
471                                    } else {
472                                        newBound = -COIN_DBL_MAX;
473                                    }
474                                    if (newBound > nowLower + 1.0e-12 && newBound > -large) {
475                                        // Tighten the lower bound
476                                        // check infeasible (relaxed)
477                                        if (nowUpper < newBound) {
478                                            if (nowUpper - newBound <
479                                                    -100.0*tolerance) {
480                                                infeasible = true;
481                                                break;
482                                            }
483                                        }
484                                    }
485                                }
486                            }
487                        }
488                    }
489                }
490                for (int i = fixColumn[iColumn]; i < numberOther; i++)
491                    mark[otherColumn[i]] = 0;
492                // reset bound unless infeasible
493                if (!infeasible || !clpSolver->isInteger(iColumn))
494                    columnUpper[iColumn] = saveUpper;
495                else if (clpSolver->isInteger(iColumn))
496                    columnLower[iColumn] = 1.0;
497            }
498        }
499        fixColumn[iColumn+1] = numberOther;
500    }
501    delete [] check;
502    delete [] mark;
503    // Now do reverse way
504    int * counts = new int [numberColumns];
505    CoinZeroN(counts, numberColumns);
506    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
507        for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++)
508            counts[otherColumn[i]]++;
509    }
510    numberOther = 0;
511    CoinBigIndex * fixColumn2 = new CoinBigIndex [numberColumns+1];
512    int * otherColumn2 = new int [fixColumn[numberColumns]];
513    fixColumn2[0] = 0;
514    for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
515        numberOther += counts[iColumn];
516        counts[iColumn] = 0;
517        fixColumn2[iColumn+1] = numberOther;
518    }
519    // Create other way
520    for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
521        for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) {
522            int jColumn = otherColumn[i];
523            CoinBigIndex put = fixColumn2[jColumn] + counts[jColumn];
524            counts[jColumn]++;
525            otherColumn2[put] = iColumn;
526        }
527    }
528    // get top layer i.e. those which are not fixed by any other
529    int kLayer = 0;
530    while (true) {
531        int numberLayered = 0;
532        for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
533            if (fix[iColumn] == kLayer) {
534                for (int i = fixColumn2[iColumn]; i < fixColumn2[iColumn+1]; i++) {
535                    int jColumn = otherColumn2[i];
536                    if (fix[jColumn] == kLayer) {
537                        fix[iColumn] = kLayer + 100;
538                    }
539                }
540            }
541            if (fix[iColumn] == kLayer) {
542                numberLayered++;
543            }
544        }
545        if (numberLayered) {
546            kLayer += 100;
547        } else {
548            break;
549        }
550    }
551    for (int iPass = 0; iPass < 2; iPass++) {
552        for (int jLayer = 0; jLayer < kLayer; jLayer++) {
553            int check[] = { -1, 0, 1, 2, 3, 4, 5, 10, 50, 100, 500, 1000, 5000, 10000, COIN_INT_MAX};
554            int nCheck = static_cast<int> (sizeof(check) / sizeof(int));
555            int countsI[20];
556            int countsC[20];
557            assert (nCheck <= 20);
558            memset(countsI, 0, nCheck*sizeof(int));
559            memset(countsC, 0, nCheck*sizeof(int));
560            check[nCheck-1] = numberColumns;
561            int numberLayered = 0;
562            int numberInteger = 0;
563            for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
564                if (fix[iColumn] == jLayer) {
565                    numberLayered++;
566                    int nFix = fixColumn[iColumn+1] - fixColumn[iColumn];
567                    if (iPass) {
568                        // just integers
569                        nFix = 0;
570                        for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) {
571                            int jColumn = otherColumn[i];
572                            if (clpSolver->isInteger(jColumn))
573                                nFix++;
574                        }
575                    }
576                    int iFix;
577                    for (iFix = 0; iFix < nCheck; iFix++) {
578                        if (nFix <= check[iFix])
579                            break;
580                    }
581                    assert (iFix < nCheck);
582                    if (clpSolver->isInteger(iColumn)) {
583                        numberInteger++;
584                        countsI[iFix]++;
585                    } else {
586                        countsC[iFix]++;
587                    }
588                }
589            }
590#ifdef COIN_DETAIL
591            if (numberLayered) {
592                printf("%d (%d integer) at priority %d\n", numberLayered, numberInteger, 1 + (jLayer / 100));
593                char buffer[50];
594                for (int i = 1; i < nCheck; i++) {
595                    if (countsI[i] || countsC[i]) {
596                        if (i == 1)
597                            sprintf(buffer, " ==    zero            ");
598                        else if (i < nCheck - 1)
599                            sprintf(buffer, "> %6d and <= %6d ", check[i-1], check[i]);
600                        else
601                            sprintf(buffer, "> %6d                ", check[i-1]);
602                        printf("%s %8d integers and %8d continuous\n", buffer, countsI[i], countsC[i]);
603                    }
604                }
605            }
606#endif
607        }
608    }
609    delete [] counts;
610    // Now do fixing
611    {
612        // switch off presolve and up weight
613        ClpSolve solveOptions;
614        //solveOptions.setPresolveType(ClpSolve::presolveOff,0);
615        solveOptions.setSolveType(ClpSolve::usePrimalorSprint);
616        //solveOptions.setSpecialOption(1,3,30); // sprint
617        int numberColumns = lpSolver->numberColumns();
618        int iColumn;
619        bool allSlack = true;
620        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
621            if (lpSolver->getColumnStatus(iColumn) == ClpSimplex::basic) {
622                allSlack = false;
623                break;
624            }
625        }
626        if (allSlack)
627            solveOptions.setSpecialOption(1, 2, 50); // idiot
628        lpSolver->setInfeasibilityCost(1.0e11);
629        lpSolver->defaultFactorizationFrequency();
630        if (doAction != 11)
631            lpSolver->initialSolve(solveOptions);
632        double * columnLower = lpSolver->columnLower();
633        double * columnUpper = lpSolver->columnUpper();
634        double * fullSolution = lpSolver->primalColumnSolution();
635        const double * dj = lpSolver->dualColumnSolution();
636        int iPass = 0;
637#define MAXPROB 2
638        ClpSimplex models[MAXPROB];
639        int kPass = -1;
640        int kLayer = 0;
641        int skipZero = 0;
642        if (skipZero2 == -1)
643            skipZero2 = 40; //-1;
644        /* 0 fixed to 0 by choice
645           1 lb of 1 by choice
646           2 fixed to 0 by another
647           3 as 2 but this go
648           -1 free
649        */
650        char * state = new char [numberColumns];
651        for (iColumn = 0; iColumn < numberColumns; iColumn++)
652            state[iColumn] = -1;
653        while (true) {
654            double largest = -0.1;
655            double smallest = 1.1;
656            int iLargest = -1;
657            int iSmallest = -1;
658            int atZero = 0;
659            int atOne = 0;
660            int toZero = 0;
661            int toOne = 0;
662            int numberFree = 0;
663            int numberGreater = 0;
664            columnLower = lpSolver->columnLower();
665            columnUpper = lpSolver->columnUpper();
666            fullSolution = lpSolver->primalColumnSolution();
667            if (doAction == 11) {
668                {
669                    double * columnLower = lpSolver->columnLower();
670                    double * columnUpper = lpSolver->columnUpper();
671                    //    lpSolver->dual();
672                    memcpy(columnLower, saveColumnLower, numberColumns*sizeof(double));
673                    memcpy(columnUpper, saveColumnUpper, numberColumns*sizeof(double));
674                    //    lpSolver->dual();
675                    int iColumn;
676                    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
677                        if (columnUpper[iColumn] > columnLower[iColumn] + 1.0e-8) {
678                            if (clpSolver->isInteger(iColumn)) {
679                                double value = lastSolution[iColumn];
680                                int iValue = static_cast<int> (value + 0.5);
681                                assert (fabs(value - static_cast<double> (iValue)) < 1.0e-3);
682                                assert (iValue >= columnLower[iColumn] &&
683                                        iValue <= columnUpper[iColumn]);
684                                columnLower[iColumn] = iValue;
685                                columnUpper[iColumn] = iValue;
686                            }
687                        }
688                    }
689                    lpSolver->initialSolve(solveOptions);
690                    memcpy(columnLower, saveColumnLower, numberColumns*sizeof(double));
691                    memcpy(columnUpper, saveColumnUpper, numberColumns*sizeof(double));
692                }
693                for (iColumn = 0; iColumn < numberColumns; iColumn++) {
694                    if (columnUpper[iColumn] > columnLower[iColumn] + 1.0e-8) {
695                        if (clpSolver->isInteger(iColumn)) {
696                            double value = lastSolution[iColumn];
697                            int iValue = static_cast<int> (value + 0.5);
698                            assert (fabs(value - static_cast<double> (iValue)) < 1.0e-3);
699                            assert (iValue >= columnLower[iColumn] &&
700                                    iValue <= columnUpper[iColumn]);
701                            if (!fix[iColumn]) {
702                                if (iValue == 0) {
703                                    state[iColumn] = 0;
704                                    assert (!columnLower[iColumn]);
705                                    columnUpper[iColumn] = 0.0;
706                                } else if (iValue == 1) {
707                                    state[iColumn] = 1;
708                                    columnLower[iColumn] = 1.0;
709                                } else {
710                                    // leave fixed
711                                    columnLower[iColumn] = iValue;
712                                    columnUpper[iColumn] = iValue;
713                                }
714                            } else if (iValue == 0) {
715                                state[iColumn] = 10;
716                                columnUpper[iColumn] = 0.0;
717                            } else {
718                                // leave fixed
719                                columnLower[iColumn] = iValue;
720                                columnUpper[iColumn] = iValue;
721                            }
722                        }
723                    }
724                }
725                int jLayer = 0;
726                int nFixed = -1;
727                int nTotalFixed = 0;
728                while (nFixed) {
729                    nFixed = 0;
730                    for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
731                        if (columnUpper[iColumn] == 0.0 && fix[iColumn] == jLayer) {
732                            for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) {
733                                int jColumn = otherColumn[i];
734                                if (columnUpper[jColumn]) {
735                                    bool canFix = true;
736                                    for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) {
737                                        int kColumn = otherColumn2[k];
738                                        if (state[kColumn] == 1) {
739                                            canFix = false;
740                                            break;
741                                        }
742                                    }
743                                    if (canFix) {
744                                        columnUpper[jColumn] = 0.0;
745                                        nFixed++;
746                                    }
747                                }
748                            }
749                        }
750                    }
751                    nTotalFixed += nFixed;
752                    jLayer += 100;
753                }
754                COIN_DETAIL_PRINT(printf("This fixes %d variables in lower priorities\n", nTotalFixed));
755                break;
756            }
757            for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
758                if (!clpSolver->isInteger(iColumn) || fix[iColumn] > kLayer)
759                    continue;
760                // skip if fixes nothing
761                if (fixColumn[iColumn+1] - fixColumn[iColumn] <= skipZero2)
762                    continue;
763                double value = fullSolution[iColumn];
764                if (value > 1.00001) {
765                    numberGreater++;
766                    continue;
767                }
768                double lower = columnLower[iColumn];
769                double upper = columnUpper[iColumn];
770                if (lower == upper) {
771                    if (lower)
772                        atOne++;
773                    else
774                        atZero++;
775                    continue;
776                }
777                if (value < 1.0e-7) {
778                    toZero++;
779                    columnUpper[iColumn] = 0.0;
780                    state[iColumn] = 10;
781                    continue;
782                }
783                if (value > 1.0 - 1.0e-7) {
784                    toOne++;
785                    columnLower[iColumn] = 1.0;
786                    state[iColumn] = 1;
787                    continue;
788                }
789                numberFree++;
790                // skip if fixes nothing
791                if (fixColumn[iColumn+1] - fixColumn[iColumn] <= skipZero)
792                    continue;
793                if (value < smallest) {
794                    smallest = value;
795                    iSmallest = iColumn;
796                }
797                if (value > largest) {
798                    largest = value;
799                    iLargest = iColumn;
800                }
801            }
802            if (toZero || toOne)
803              COIN_DETAIL_PRINT(printf("%d at 0 fixed and %d at one fixed\n", toZero, toOne));
804            COIN_DETAIL_PRINT(printf("%d variables free, %d fixed to 0, %d to 1 - smallest %g, largest %g\n",
805                                     numberFree, atZero, atOne, smallest, largest));
806            if (numberGreater && !iPass)
807              COIN_DETAIL_PRINT(printf("%d variables have value > 1.0\n", numberGreater));
808            //skipZero2=0; // leave 0 fixing
809            int jLayer = 0;
810            int nFixed = -1;
811            int nTotalFixed = 0;
812            while (nFixed) {
813                nFixed = 0;
814                for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
815                    if (columnUpper[iColumn] == 0.0 && fix[iColumn] == jLayer) {
816                        for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) {
817                            int jColumn = otherColumn[i];
818                            if (columnUpper[jColumn]) {
819                                bool canFix = true;
820                                for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) {
821                                    int kColumn = otherColumn2[k];
822                                    if (state[kColumn] == 1) {
823                                        canFix = false;
824                                        break;
825                                    }
826                                }
827                                if (canFix) {
828                                    columnUpper[jColumn] = 0.0;
829                                    nFixed++;
830                                }
831                            }
832                        }
833                    }
834                }
835                nTotalFixed += nFixed;
836                jLayer += 100;
837            }
838            COIN_DETAIL_PRINT(printf("This fixes %d variables in lower priorities\n", nTotalFixed));
839            if (iLargest < 0 || numberFree <= leaveIntFree)
840                break;
841            double movement;
842            int way;
843            if (smallest <= 1.0 - largest && smallest < 0.2 && largest < fixAboveValue) {
844                columnUpper[iSmallest] = 0.0;
845                state[iSmallest] = 0;
846                movement = smallest;
847                way = -1;
848            } else {
849                columnLower[iLargest] = 1.0;
850                state[iLargest] = 1;
851                movement = 1.0 - largest;
852                way = 1;
853            }
854            double saveObj = lpSolver->objectiveValue();
855            iPass++;
856            kPass = iPass % MAXPROB;
857            models[kPass] = *lpSolver;
858            if (way == -1) {
859                // fix others
860                for (int i = fixColumn[iSmallest]; i < fixColumn[iSmallest+1]; i++) {
861                    int jColumn = otherColumn[i];
862                    if (state[jColumn] == -1) {
863                        columnUpper[jColumn] = 0.0;
864                        state[jColumn] = 3;
865                    }
866                }
867            }
868            double maxCostUp = COIN_DBL_MAX;
869            objective = lpSolver->getObjCoefficients() ;
870            if (way == -1)
871                maxCostUp = (1.0 - movement) * objective[iSmallest];
872            lpSolver->setDualObjectiveLimit(saveObj + maxCostUp);
873            crunchIt(lpSolver);
874            double moveObj = lpSolver->objectiveValue() - saveObj;
875            COIN_DETAIL_PRINT(printf("movement %s was %g costing %g\n",
876                                     (way == -1) ? "down" : "up", movement, moveObj));
877            if (way == -1 && (moveObj >= maxCostUp || lpSolver->status())) {
878                // go up
879                columnLower = models[kPass].columnLower();
880                columnUpper = models[kPass].columnUpper();
881                columnLower[iSmallest] = 1.0;
882                columnUpper[iSmallest] = saveColumnUpper[iSmallest];
883                *lpSolver = models[kPass];
884                columnLower = lpSolver->columnLower();
885                columnUpper = lpSolver->columnUpper();
886                fullSolution = lpSolver->primalColumnSolution();
887                dj = lpSolver->dualColumnSolution();
888                columnLower[iSmallest] = 1.0;
889                columnUpper[iSmallest] = saveColumnUpper[iSmallest];
890                state[iSmallest] = 1;
891                // unfix others
892                for (int i = fixColumn[iSmallest]; i < fixColumn[iSmallest+1]; i++) {
893                    int jColumn = otherColumn[i];
894                    if (state[jColumn] == 3) {
895                        columnUpper[jColumn] = saveColumnUpper[jColumn];
896                        state[jColumn] = -1;
897                    }
898                }
899                crunchIt(lpSolver);
900            }
901            models[kPass] = *lpSolver;
902        }
903        lpSolver->dual();
904        COIN_DETAIL_PRINT(printf("Fixing took %g seconds\n", CoinCpuTime() - time1));
905        columnLower = lpSolver->columnLower();
906        columnUpper = lpSolver->columnUpper();
907        fullSolution = lpSolver->primalColumnSolution();
908        dj = lpSolver->dualColumnSolution();
909        int * sort = new int[numberColumns];
910        double * dsort = new double[numberColumns];
911        int chunk = 20;
912        int iRelax = 0;
913        //double fractionFixed=6.0/8.0;
914        // relax while lots fixed
915        while (true) {
916            if (skipZero2 > 10 && doAction < 10)
917                break;
918            iRelax++;
919            int n = 0;
920            double sum0 = 0.0;
921            double sum00 = 0.0;
922            double sum1 = 0.0;
923            for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
924                if (!clpSolver->isInteger(iColumn) || fix[iColumn] > kLayer)
925                    continue;
926                // skip if fixes nothing
927                if (fixColumn[iColumn+1] - fixColumn[iColumn] == 0 && doAction < 10)
928                    continue;
929                double djValue = dj[iColumn];
930                if (state[iColumn] == 1) {
931                    assert (columnLower[iColumn]);
932                    assert (fullSolution[iColumn] > 0.1);
933                    if (djValue > 0.0) {
934                        //printf("YY dj of %d at %g is %g\n",iColumn,value,djValue);
935                        sum1 += djValue;
936                        sort[n] = iColumn;
937                        dsort[n++] = -djValue;
938                    } else {
939                        //printf("dj of %d at %g is %g\n",iColumn,value,djValue);
940                    }
941                } else if (state[iColumn] == 0 || state[iColumn] == 10) {
942                    assert (fullSolution[iColumn] < 0.1);
943                    assert (!columnUpper[iColumn]);
944                    double otherValue = 0.0;
945                    int nn = 0;
946                    for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) {
947                        int jColumn = otherColumn[i];
948                        if (columnUpper[jColumn] == 0.0) {
949                            if (dj[jColumn] < -1.0e-5) {
950                                nn++;
951                                otherValue += dj[jColumn]; // really need to look at rest
952                            }
953                        }
954                    }
955                    if (djValue < -1.0e-2 || otherValue < -1.0e-2) {
956                        //printf("XX dj of %d at %g is %g - %d out of %d contribute %g\n",iColumn,value,djValue,
957                        // nn,fixColumn[iColumn+1]-fixColumn[iColumn],otherValue);
958                        if (djValue < 1.0e-8) {
959                            sum0 -= djValue;
960                            sum00 -= otherValue;
961                            sort[n] = iColumn;
962                            if (djValue < -1.0e-2)
963                                dsort[n++] = djValue + otherValue;
964                            else
965                                dsort[n++] = djValue + 0.001 * otherValue;
966                        }
967                    } else {
968                        //printf("dj of %d at %g is %g - no contribution from %d\n",iColumn,value,djValue,
969                        //   fixColumn[iColumn+1]-fixColumn[iColumn]);
970                    }
971                }
972            }
973            CoinSort_2(dsort, dsort + n, sort);
974            double * originalColumnLower = saveColumnLower;
975            double * originalColumnUpper = saveColumnUpper;
976            double * lo = CoinCopyOfArray(columnLower, numberColumns);
977            double * up = CoinCopyOfArray(columnUpper, numberColumns);
978            for (int k = 0; k < CoinMin(chunk, n); k++) {
979                iColumn = sort[k];
980                state[iColumn] = -2;
981            }
982            memcpy(columnLower, originalColumnLower, numberColumns*sizeof(double));
983            memcpy(columnUpper, originalColumnUpper, numberColumns*sizeof(double));
984            int nFixed = 0;
985            int nFixed0 = 0;
986            int nFixed1 = 0;
987            for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
988                if (state[iColumn] == 0 || state[iColumn] == 10) {
989                    columnUpper[iColumn] = 0.0;
990                    assert (lo[iColumn] == 0.0);
991                    nFixed++;
992                    nFixed0++;
993                    for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) {
994                        int jColumn = otherColumn[i];
995                        if (columnUpper[jColumn]) {
996                            bool canFix = true;
997                            for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) {
998                                int kColumn = otherColumn2[k];
999                                if (state[kColumn] == 1 || state[kColumn] == -2) {
1000                                    canFix = false;
1001                                    break;
1002                                }
1003                            }
1004                            if (canFix) {
1005                                columnUpper[jColumn] = 0.0;
1006                                assert (lo[jColumn] == 0.0);
1007                                nFixed++;
1008                            }
1009                        }
1010                    }
1011                } else if (state[iColumn] == 1) {
1012                    columnLower[iColumn] = 1.0;
1013                    nFixed1++;
1014                }
1015            }
1016            COIN_DETAIL_PRINT(printf("%d fixed %d orig 0 %d 1\n", nFixed, nFixed0, nFixed1));
1017            int jLayer = 0;
1018            nFixed = -1;
1019            int nTotalFixed = 0;
1020            while (nFixed) {
1021                nFixed = 0;
1022                for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
1023                    if (columnUpper[iColumn] == 0.0 && fix[iColumn] == jLayer) {
1024                        for (int i = fixColumn[iColumn]; i < fixColumn[iColumn+1]; i++) {
1025                            int jColumn = otherColumn[i];
1026                            if (columnUpper[jColumn]) {
1027                                bool canFix = true;
1028                                for (int k = fixColumn2[jColumn]; k < fixColumn2[jColumn+1]; k++) {
1029                                    int kColumn = otherColumn2[k];
1030                                    if (state[kColumn] == 1 || state[kColumn] == -2) {
1031                                        canFix = false;
1032                                        break;
1033                                    }
1034                                }
1035                                if (canFix) {
1036                                    columnUpper[jColumn] = 0.0;
1037                                    assert (lo[jColumn] == 0.0);
1038                                    nFixed++;
1039                                }
1040                            }
1041                        }
1042                    }
1043                }
1044                nTotalFixed += nFixed;
1045                jLayer += 100;
1046            }
1047            nFixed = 0;
1048            int nFixedI = 0;
1049            for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
1050                if (columnLower[iColumn] == columnUpper[iColumn]) {
1051                    if (clpSolver->isInteger(iColumn))
1052                        nFixedI++;
1053                    nFixed++;
1054                }
1055            }
1056            COIN_DETAIL_PRINT(printf("This fixes %d variables in lower priorities - total %d (%d integer) - all target %d, int target %d\n",
1057                                     nTotalFixed, nFixed, nFixedI, static_cast<int>(fractionFixed*numberColumns), static_cast<int> (fractionIntFixed*numberInteger)));
1058            int nBad = 0;
1059            int nRelax = 0;
1060            for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
1061                if (lo[iColumn] < columnLower[iColumn] ||
1062                        up[iColumn] > columnUpper[iColumn]) {
1063                    COIN_DETAIL_PRINT(printf("bad %d old %g %g, new %g %g\n", iColumn, lo[iColumn], up[iColumn],
1064                                             columnLower[iColumn], columnUpper[iColumn]));
1065                    nBad++;
1066                }
1067                if (lo[iColumn] > columnLower[iColumn] ||
1068                        up[iColumn] < columnUpper[iColumn]) {
1069                    nRelax++;
1070                }
1071            }
1072            COIN_DETAIL_PRINT(printf("%d relaxed\n", nRelax));
1073            if (iRelax > 20 && nRelax == chunk)
1074                nRelax = 0;
1075            if (iRelax > 50)
1076                nRelax = 0;
1077            assert (!nBad);
1078            delete [] lo;
1079            delete [] up;
1080            lpSolver->primal(1);
1081            if (nFixed < fractionFixed*numberColumns || nFixedI < fractionIntFixed*numberInteger || !nRelax)
1082                break;
1083        }
1084        delete [] state;
1085        delete [] sort;
1086        delete [] dsort;
1087    }
1088    delete [] fix;
1089    delete [] fixColumn;
1090    delete [] otherColumn;
1091    delete [] otherColumn2;
1092    delete [] fixColumn2;
1093    // See if was presolved
1094    if (originalColumns) {
1095        columnLower = lpSolver->columnLower();
1096        columnUpper = lpSolver->columnUpper();
1097        for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
1098            saveColumnLower[iColumn] = columnLower[iColumn];
1099            saveColumnUpper[iColumn] = columnUpper[iColumn];
1100        }
1101        pinfo.postsolve(true);
1102        columnLower = originalLpSolver->columnLower();
1103        columnUpper = originalLpSolver->columnUpper();
1104        double * newColumnLower = lpSolver->columnLower();
1105        double * newColumnUpper = lpSolver->columnUpper();
1106        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1107            int jColumn = originalColumns[iColumn];
1108            columnLower[jColumn] = CoinMax(columnLower[jColumn], newColumnLower[iColumn]);
1109            columnUpper[jColumn] = CoinMin(columnUpper[jColumn], newColumnUpper[iColumn]);
1110        }
1111        numberColumns = originalLpSolver->numberColumns();
1112        delete [] originalColumns;
1113    }
1114    delete [] saveColumnLower;
1115    delete [] saveColumnUpper;
1116    if (!originalColumns) {
1117        // Basis
1118        memcpy(originalLpSolver->statusArray(), lpSolver->statusArray(), numberRows + numberColumns);
1119        memcpy(originalLpSolver->primalColumnSolution(), lpSolver->primalColumnSolution(), numberColumns*sizeof(double));
1120        memcpy(originalLpSolver->primalRowSolution(), lpSolver->primalRowSolution(), numberRows*sizeof(double));
1121        // Fix in solver
1122        columnLower = lpSolver->columnLower();
1123        columnUpper = lpSolver->columnUpper();
1124    }
1125    double * originalColumnLower = originalLpSolver->columnLower();
1126    double * originalColumnUpper = originalLpSolver->columnUpper();
1127    // number fixed
1128    doAction = 0;
1129    for ( iColumn = 0; iColumn < numberColumns; iColumn++) {
1130        originalColumnLower[iColumn] = columnLower[iColumn];
1131        originalColumnUpper[iColumn] = columnUpper[iColumn];
1132        if (columnLower[iColumn] == columnUpper[iColumn])
1133            doAction++;
1134    }
1135    COIN_DETAIL_PRINT(printf("%d fixed by vub preprocessing\n", doAction));
1136    if (originalColumns) {
1137        originalLpSolver->initialSolve();
1138    }
1139    delete clpSolver;
1140    return NULL;
1141}
1142
1143int doHeuristics(CbcModel * model, int type, CbcOrClpParam* parameters_,
1144                 int numberParameters_,int noPrinting_,int initialPumpTune)
1145{
1146#ifdef JJF_ZERO //NEW_STYLE_SOLVER==0
1147    CbcOrClpParam * parameters_ = parameters;
1148    int numberParameters_ = numberParameters;
1149    bool noPrinting_ = noPrinting_;
1150#endif
1151    char generalPrint[10000];
1152    CoinMessages generalMessages = model->messages();
1153    CoinMessageHandler * generalMessageHandler = model->messageHandler();
1154    //generalMessageHandler->setPrefix(false);
1155    bool anyToDo = false;
1156    int logLevel = parameters_[whichParam(CLP_PARAM_INT_LOGLEVEL, numberParameters_, parameters_)].intValue();
1157    int useFpump = parameters_[whichParam(CBC_PARAM_STR_FPUMP, numberParameters_, parameters_)].currentOptionAsInteger();
1158    int useRounding = parameters_[whichParam(CBC_PARAM_STR_ROUNDING, numberParameters_, parameters_)].currentOptionAsInteger();
1159    int useGreedy = parameters_[whichParam(CBC_PARAM_STR_GREEDY, numberParameters_, parameters_)].currentOptionAsInteger();
1160    int useCombine = parameters_[whichParam(CBC_PARAM_STR_COMBINE, numberParameters_, parameters_)].currentOptionAsInteger();
1161    int useProximity = parameters_[whichParam(CBC_PARAM_STR_PROXIMITY, numberParameters_, parameters_)].currentOptionAsInteger();
1162    int useCrossover = parameters_[whichParam(CBC_PARAM_STR_CROSSOVER2, numberParameters_, parameters_)].currentOptionAsInteger();
1163    //int usePivotC = parameters_[whichParam(CBC_PARAM_STR_PIVOTANDCOMPLEMENT, numberParameters_, parameters_)].currentOptionAsInteger();
1164    int usePivotF = parameters_[whichParam(CBC_PARAM_STR_PIVOTANDFIX, numberParameters_, parameters_)].currentOptionAsInteger();
1165    int useRand = parameters_[whichParam(CBC_PARAM_STR_RANDROUND, numberParameters_, parameters_)].currentOptionAsInteger();
1166    int useRINS = parameters_[whichParam(CBC_PARAM_STR_RINS, numberParameters_, parameters_)].currentOptionAsInteger();
1167    int useRENS = parameters_[whichParam(CBC_PARAM_STR_RENS, numberParameters_, parameters_)].currentOptionAsInteger();
1168    int useDINS = parameters_[whichParam(CBC_PARAM_STR_DINS, numberParameters_, parameters_)].currentOptionAsInteger();
1169    int useDIVING2 = parameters_[whichParam(CBC_PARAM_STR_DIVINGS, numberParameters_, parameters_)].currentOptionAsInteger();
1170    int useNaive = parameters_[whichParam(CBC_PARAM_STR_NAIVE, numberParameters_, parameters_)].currentOptionAsInteger();
1171    int useDW = parameters_[whichParam(CBC_PARAM_STR_DW, numberParameters_, parameters_)].currentOptionAsInteger();
1172    int kType = (type < 10) ? type : 1;
1173    assert (kType == 1 || kType == 2);
1174    // FPump done first as it only works if no solution
1175    if (useFpump >= kType && useFpump <= kType + 1) {
1176        anyToDo = true;
1177        CbcHeuristicFPump heuristic4(*model);
1178        double dextra3 = parameters_[whichParam(CBC_PARAM_DBL_SMALLBAB, numberParameters_, parameters_)].doubleValue();
1179        heuristic4.setFractionSmall(dextra3);
1180        double dextra1 = parameters_[whichParam(CBC_PARAM_DBL_ARTIFICIALCOST, numberParameters_, parameters_)].doubleValue();
1181        if (dextra1)
1182            heuristic4.setArtificialCost(dextra1);
1183        heuristic4.setMaximumPasses(parameters_[whichParam(CBC_PARAM_INT_FPUMPITS, numberParameters_, parameters_)].intValue());
1184        if (parameters_[whichParam(CBC_PARAM_INT_FPUMPITS, numberParameters_, parameters_)].intValue() == 21)
1185            heuristic4.setIterationRatio(1.0);
1186        int pumpTune = parameters_[whichParam(CBC_PARAM_INT_FPUMPTUNE, numberParameters_, parameters_)].intValue();
1187        int pumpTune2 = parameters_[whichParam(CBC_PARAM_INT_FPUMPTUNE2, numberParameters_, parameters_)].intValue();
1188        if (pumpTune > 0) {
1189            bool printStuff = (pumpTune != initialPumpTune || logLevel > 1 || pumpTune2 > 0)
1190                              && !noPrinting_;
1191            if (printStuff) {
1192                generalMessageHandler->message(CBC_GENERAL, generalMessages)
1193                << "Options for feasibility pump - "
1194                << CoinMessageEol;
1195            }
1196            /*
1197            >=10000000 for using obj
1198            >=1000000 use as accumulate switch
1199            >=1000 use index+1 as number of large loops
1200            >=100 use dextra1 as cutoff
1201            %100 == 10,20 etc for experimentation
1202            1 == fix ints at bounds, 2 fix all integral ints, 3 and continuous at bounds
1203            4 and static continuous, 5 as 3 but no internal integers
1204            6 as 3 but all slack basis!
1205            */
1206            double value = model->solver()->getObjSense() * model->solver()->getObjValue();
1207            int w = pumpTune / 10;
1208            int i = w % 10;
1209            w /= 10;
1210            int c = w % 10;
1211            w /= 10;
1212            int r = w;
1213            int accumulate = r / 1000;
1214            r -= 1000 * accumulate;
1215            if (accumulate >= 10) {
1216                int which = accumulate / 10;
1217                accumulate -= 10 * which;
1218                which--;
1219                // weights and factors
1220                double weight[] = {0.01, 0.01, 0.1, 0.1, 0.5, 0.5, 1.0, 1.0, 5.0, 5.0};
1221                double factor[] = {0.1, 0.5, 0.1, 0.5, 0.1, 0.5, 0.1, 0.5, 0.1, 0.5};
1222                heuristic4.setInitialWeight(weight[which]);
1223                heuristic4.setWeightFactor(factor[which]);
1224                if (printStuff) {
1225                    sprintf(generalPrint, "Initial weight for objective %g, decay factor %g",
1226                            weight[which], factor[which]);
1227                    generalMessageHandler->message(CBC_GENERAL, generalMessages)
1228                    << generalPrint
1229                    << CoinMessageEol;
1230                }
1231
1232            }
1233            // fake cutoff
1234            if (c) {
1235                double cutoff;
1236                model->solver()->getDblParam(OsiDualObjectiveLimit, cutoff);
1237                cutoff = CoinMin(cutoff, value + 0.05 * fabs(value) * c);
1238                double fakeCutoff = parameters_[whichParam(CBC_PARAM_DBL_FAKECUTOFF, numberParameters_, parameters_)].doubleValue();
1239                if (fakeCutoff)
1240                    cutoff = fakeCutoff;
1241                heuristic4.setFakeCutoff(cutoff);
1242                if (printStuff) {
1243                    sprintf(generalPrint, "Fake cutoff of %g", cutoff);
1244                    generalMessageHandler->message(CBC_GENERAL, generalMessages)
1245                    << generalPrint
1246                    << CoinMessageEol;
1247                }
1248            }
1249            int offRandomEtc = 0;
1250            if (pumpTune2) {
1251                if ((pumpTune2 / 1000) != 0) {
1252                    offRandomEtc = 1000000 * (pumpTune2 / 1000);
1253                    if (printStuff) {
1254                        generalMessageHandler->message(CBC_GENERAL, generalMessages)
1255                        << "Feasibility pump may run twice"
1256                        << CoinMessageEol;
1257                    }
1258                    pumpTune2 = pumpTune2 % 1000;
1259                }
1260                if ((pumpTune2 / 100) != 0) {
1261                    offRandomEtc += 100 * (pumpTune2 / 100);
1262                    if (printStuff) {
1263                        generalMessageHandler->message(CBC_GENERAL, generalMessages)
1264                        << "Not using randomized objective"
1265                        << CoinMessageEol;
1266                    }
1267                }
1268                int maxAllowed = pumpTune2 % 100;
1269                if (maxAllowed) {
1270                    offRandomEtc += 1000 * maxAllowed;
1271                    if (printStuff) {
1272                        sprintf(generalPrint, "Fixing if same for %d passes",
1273                                maxAllowed);
1274                        generalMessageHandler->message(CBC_GENERAL, generalMessages)
1275                        << generalPrint
1276                        << CoinMessageEol;
1277                    }
1278                }
1279            }
1280            if (accumulate) {
1281                heuristic4.setAccumulate(accumulate);
1282                if (printStuff) {
1283                    if (accumulate) {
1284                        sprintf(generalPrint, "Accumulate of %d", accumulate);
1285                        generalMessageHandler->message(CBC_GENERAL, generalMessages)
1286                        << generalPrint
1287                        << CoinMessageEol;
1288                    }
1289                }
1290            }
1291            if (r) {
1292                // also set increment
1293                //double increment = (0.01*i+0.005)*(fabs(value)+1.0e-12);
1294                double increment = 0.0;
1295                double fakeIncrement = parameters_[whichParam(CBC_PARAM_DBL_FAKEINCREMENT, numberParameters_, parameters_)].doubleValue();
1296                if (fakeIncrement)
1297                    increment = fakeIncrement;
1298                if (increment>=0.0)
1299                  heuristic4.setAbsoluteIncrement(increment);
1300                else
1301                  heuristic4.setRelativeIncrement(-increment);
1302                heuristic4.setMaximumRetries(r + 1);
1303                if (printStuff) {
1304                    if (increment) {
1305                      if (increment>0.0)
1306                        sprintf(generalPrint, "Absolute increment of %g", increment);
1307                      else
1308                        sprintf(generalPrint, "Relative increment of %g", -increment);
1309                        generalMessageHandler->message(CBC_GENERAL, generalMessages)
1310                        << generalPrint
1311                        << CoinMessageEol;
1312                    }
1313                    sprintf(generalPrint, "%d retries", r + 1);
1314                    generalMessageHandler->message(CBC_GENERAL, generalMessages)
1315                    << generalPrint
1316                    << CoinMessageEol;
1317                }
1318            }
1319            if (i + offRandomEtc) {
1320                heuristic4.setFeasibilityPumpOptions(i*10 + offRandomEtc);
1321                if (printStuff) {
1322                    sprintf(generalPrint, "Feasibility pump options of %d",
1323                            i*10 + offRandomEtc);
1324                    generalMessageHandler->message(CBC_GENERAL, generalMessages)
1325                    << generalPrint
1326                    << CoinMessageEol;
1327                }
1328            }
1329            pumpTune = pumpTune % 100;
1330            if (pumpTune == 6)
1331                pumpTune = 13;
1332            heuristic4.setWhen((pumpTune % 10) + 10);
1333            if (printStuff) {
1334                sprintf(generalPrint, "Tuning (fixing) %d", pumpTune % 10);
1335                generalMessageHandler->message(CBC_GENERAL, generalMessages)
1336                << generalPrint
1337                << CoinMessageEol;
1338            }
1339        }
1340        heuristic4.setHeuristicName("feasibility pump");
1341        //#define ROLF
1342#ifdef ROLF
1343        CbcHeuristicFPump pump(*model);
1344        pump.setMaximumTime(60);
1345        pump.setMaximumPasses(100);
1346        pump.setMaximumRetries(1);
1347        pump.setFixOnReducedCosts(0);
1348        pump.setHeuristicName("Feasibility pump");
1349        pump.setFractionSmall(1.0);
1350        pump.setWhen(13);
1351        model->addHeuristic(&pump);
1352#else
1353        model->addHeuristic(&heuristic4);
1354#endif
1355    }
1356    if (useRounding >= type && useRounding >= kType && useRounding <= kType + 1) {
1357        CbcRounding heuristic1(*model);
1358        heuristic1.setHeuristicName("rounding");
1359        model->addHeuristic(&heuristic1) ;
1360        anyToDo = true;
1361    }
1362    if (useGreedy >= type && useGreedy >= kType && useGreedy <= kType + 1) {
1363        CbcHeuristicGreedyCover heuristic3(*model);
1364        heuristic3.setHeuristicName("greedy cover");
1365        CbcHeuristicGreedyEquality heuristic3a(*model);
1366        heuristic3a.setHeuristicName("greedy equality");
1367        model->addHeuristic(&heuristic3);
1368        model->addHeuristic(&heuristic3a);
1369        anyToDo = true;
1370    }
1371    if ((useRENS==7 && kType==1) || (useRENS==8 && kType==2)) {
1372        useRENS=1+2*(useRENS-7);
1373        CbcHeuristicRENS heuristic6a(*model);
1374        heuristic6a.setHeuristicName("RENSdj");
1375        heuristic6a.setFractionSmall(0.6/*3.4*/);
1376        heuristic6a.setFeasibilityPumpOptions(3);
1377        heuristic6a.setNumberNodes(10);
1378        heuristic6a.setWhereFrom(4*256+4*1);
1379        heuristic6a.setWhen(2);
1380        heuristic6a.setRensType(1+16);
1381        model->addHeuristic(&heuristic6a) ;
1382        heuristic6a.setHeuristicName("RENSub");
1383        heuristic6a.setFractionSmall(0.4);
1384        heuristic6a.setFeasibilityPumpOptions(1008003);
1385        heuristic6a.setNumberNodes(50);
1386        heuristic6a.setWhereFrom(4*256+4*1);
1387        heuristic6a.setWhen(2);
1388        heuristic6a.setRensType(2+16);
1389        model->addHeuristic(&heuristic6a) ;
1390    }
1391    if (useRENS >= kType && useRENS <= kType + 1) {
1392#ifndef JJF_ONE
1393        CbcHeuristicRENS heuristic6(*model);
1394        heuristic6.setHeuristicName("RENS");
1395        heuristic6.setFractionSmall(0.4);
1396        heuristic6.setFeasibilityPumpOptions(1008003);
1397        int nodes [] = { -2, 50, 50, 50, 200, 1000, 10000};
1398        heuristic6.setNumberNodes(nodes[useRENS]);
1399#else
1400        CbcHeuristicVND heuristic6(*model);
1401        heuristic6.setHeuristicName("VND");
1402        heuristic5.setFractionSmall(0.5);
1403        heuristic5.setDecayFactor(5.0);
1404#endif
1405        model->addHeuristic(&heuristic6) ;
1406        anyToDo = true;
1407    }
1408    if (useNaive >= kType && useNaive <= kType + 1) {
1409        CbcHeuristicNaive heuristic5b(*model);
1410        heuristic5b.setHeuristicName("Naive");
1411        heuristic5b.setFractionSmall(0.4);
1412        heuristic5b.setNumberNodes(50);
1413        model->addHeuristic(&heuristic5b) ;
1414        anyToDo = true;
1415    }
1416    int useDIVING = 0;
1417    {
1418        int useD;
1419        useD = parameters_[whichParam(CBC_PARAM_STR_DIVINGV, numberParameters_, parameters_)].currentOptionAsInteger();
1420        useDIVING |= 1 * ((useD >= kType) ? 1 : 0);
1421        useD = parameters_[whichParam(CBC_PARAM_STR_DIVINGG, numberParameters_, parameters_)].currentOptionAsInteger();
1422        useDIVING |= 2 * ((useD >= kType) ? 1 : 0);
1423        useD = parameters_[whichParam(CBC_PARAM_STR_DIVINGF, numberParameters_, parameters_)].currentOptionAsInteger();
1424        useDIVING |= 4 * ((useD >= kType) ? 1 : 0);
1425        useD = parameters_[whichParam(CBC_PARAM_STR_DIVINGC, numberParameters_, parameters_)].currentOptionAsInteger();
1426        useDIVING |= 8 * ((useD >= kType) ? 1 : 0);
1427        useD = parameters_[whichParam(CBC_PARAM_STR_DIVINGL, numberParameters_, parameters_)].currentOptionAsInteger();
1428        useDIVING |= 16 * ((useD >= kType) ? 1 : 0);
1429        useD = parameters_[whichParam(CBC_PARAM_STR_DIVINGP, numberParameters_, parameters_)].currentOptionAsInteger();
1430        useDIVING |= 32 * ((useD >= kType) ? 1 : 0);
1431    }
1432    if (useDIVING2 >= kType && useDIVING2 <= kType + 1) {
1433        int diveOptions = parameters_[whichParam(CBC_PARAM_INT_DIVEOPT, numberParameters_, parameters_)].intValue();
1434        if (diveOptions < 0 || diveOptions > 10)
1435            diveOptions = 2;
1436        CbcHeuristicJustOne heuristicJustOne(*model);
1437        heuristicJustOne.setHeuristicName("DiveAny");
1438        heuristicJustOne.setWhen(diveOptions);
1439        // add in others
1440        CbcHeuristicDiveCoefficient heuristicDC(*model);
1441        heuristicDC.setHeuristicName("DiveCoefficient");
1442        heuristicJustOne.addHeuristic(&heuristicDC, 1.0) ;
1443        CbcHeuristicDiveFractional heuristicDF(*model);
1444        heuristicDF.setHeuristicName("DiveFractional");
1445        heuristicJustOne.addHeuristic(&heuristicDF, 1.0) ;
1446        CbcHeuristicDiveGuided heuristicDG(*model);
1447        heuristicDG.setHeuristicName("DiveGuided");
1448        heuristicJustOne.addHeuristic(&heuristicDG, 1.0) ;
1449        CbcHeuristicDiveLineSearch heuristicDL(*model);
1450        heuristicDL.setHeuristicName("DiveLineSearch");
1451        heuristicJustOne.addHeuristic(&heuristicDL, 1.0) ;
1452        CbcHeuristicDivePseudoCost heuristicDP(*model);
1453        heuristicDP.setHeuristicName("DivePseudoCost");
1454        heuristicJustOne.addHeuristic(&heuristicDP, 1.0) ;
1455        CbcHeuristicDiveVectorLength heuristicDV(*model);
1456        heuristicDV.setHeuristicName("DiveVectorLength");
1457        heuristicJustOne.addHeuristic(&heuristicDV, 1.0) ;
1458        // Now normalize probabilities
1459        heuristicJustOne.normalizeProbabilities();
1460        model->addHeuristic(&heuristicJustOne) ;
1461    }
1462
1463    if (useDIVING > 0) {
1464        int majorIterations=64;
1465        int diveOptions2 = parameters_[whichParam(CBC_PARAM_INT_DIVEOPT, numberParameters_, parameters_)].intValue();
1466        int diveOptions;
1467        if (diveOptions2 > 99) {
1468            // switch on various active set stuff
1469            diveOptions = diveOptions2%100;
1470            diveOptions2 /= 100;
1471        } else {
1472            diveOptions = diveOptions2;
1473            diveOptions2 = 0;
1474        }
1475        if (diveOptions < 0 || diveOptions > 9)
1476            diveOptions = 2;
1477        if ((useDIVING&1) != 0) {
1478            CbcHeuristicDiveVectorLength heuristicDV(*model);
1479            heuristicDV.setHeuristicName("DiveVectorLength");
1480            heuristicDV.setWhen(diveOptions);
1481            if (diveOptions2) {
1482              heuristicDV.setMaxIterations(majorIterations);
1483              heuristicDV.setPercentageToFix(0.0);
1484              heuristicDV.setMaxSimplexIterations(COIN_INT_MAX);
1485              heuristicDV.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
1486            }
1487            model->addHeuristic(&heuristicDV) ;
1488        }
1489        if ((useDIVING&2) != 0) {
1490            CbcHeuristicDiveGuided heuristicDG(*model);
1491            heuristicDG.setHeuristicName("DiveGuided");
1492            heuristicDG.setWhen(diveOptions);
1493            if (diveOptions2) {
1494              heuristicDG.setMaxIterations(majorIterations);
1495              heuristicDG.setPercentageToFix(0.0);
1496              heuristicDG.setMaxSimplexIterations(COIN_INT_MAX);
1497              heuristicDG.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
1498            }
1499            model->addHeuristic(&heuristicDG) ;
1500        }
1501        if ((useDIVING&4) != 0) {
1502            CbcHeuristicDiveFractional heuristicDF(*model);
1503            heuristicDF.setHeuristicName("DiveFractional");
1504            heuristicDF.setWhen(diveOptions);
1505            if (diveOptions2) {
1506              heuristicDF.setMaxIterations(majorIterations);
1507              heuristicDF.setPercentageToFix(0.0);
1508              heuristicDF.setMaxSimplexIterations(COIN_INT_MAX);
1509              heuristicDF.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
1510            }
1511            model->addHeuristic(&heuristicDF) ;
1512        }
1513        if ((useDIVING&8) != 0) {
1514            CbcHeuristicDiveCoefficient heuristicDC(*model);
1515            heuristicDC.setHeuristicName("DiveCoefficient");
1516            heuristicDC.setWhen(diveOptions);
1517            if (diveOptions2) {
1518              heuristicDC.setMaxIterations(majorIterations);
1519              heuristicDC.setPercentageToFix(0.0);
1520              heuristicDC.setMaxSimplexIterations(COIN_INT_MAX);
1521              heuristicDC.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
1522            }
1523            model->addHeuristic(&heuristicDC) ;
1524        }
1525        if ((useDIVING&16) != 0) {
1526            CbcHeuristicDiveLineSearch heuristicDL(*model);
1527            heuristicDL.setHeuristicName("DiveLineSearch");
1528            heuristicDL.setWhen(diveOptions);
1529            if (diveOptions2) {
1530              heuristicDL.setMaxIterations(majorIterations);
1531              heuristicDL.setPercentageToFix(0.0);
1532              heuristicDL.setMaxSimplexIterations(COIN_INT_MAX);
1533              heuristicDL.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
1534            }
1535            model->addHeuristic(&heuristicDL) ;
1536        }
1537        if ((useDIVING&32) != 0) {
1538            CbcHeuristicDivePseudoCost heuristicDP(*model);
1539            heuristicDP.setHeuristicName("DivePseudoCost");
1540            heuristicDP.setWhen(diveOptions /*+ diveOptions2*/);
1541            if (diveOptions2) {
1542              heuristicDP.setMaxIterations(majorIterations);
1543              heuristicDP.setPercentageToFix(0.0);
1544              heuristicDP.setMaxSimplexIterations(COIN_INT_MAX);
1545              heuristicDP.setMaxSimplexIterationsAtRoot(COIN_INT_MAX-(diveOptions2-1));
1546            }
1547            model->addHeuristic(&heuristicDP) ;
1548        }
1549        anyToDo = true;
1550    }
1551#ifdef JJF_ZERO
1552    if (usePivotC >= type && usePivotC <= kType + 1) {
1553        CbcHeuristicPivotAndComplement heuristic(*model);
1554        heuristic.setHeuristicName("pivot and complement");
1555        heuristic.setFractionSmall(10.0); // normally 0.5
1556        model->addHeuristic(&heuristic);
1557        anyToDo = true;
1558    }
1559#endif
1560    if (usePivotF >= type && usePivotF <= kType + 1) {
1561        CbcHeuristicPivotAndFix heuristic(*model);
1562        heuristic.setHeuristicName("pivot and fix");
1563        heuristic.setFractionSmall(10.0); // normally 0.5
1564        model->addHeuristic(&heuristic);
1565        anyToDo = true;
1566    }
1567    if (useRand >= type && useRand <= kType + 1) {
1568        CbcHeuristicRandRound heuristic(*model);
1569        heuristic.setHeuristicName("randomized rounding");
1570        heuristic.setFractionSmall(10.0); // normally 0.5
1571        model->addHeuristic(&heuristic);
1572        anyToDo = true;
1573    }
1574    if (useDINS >= kType && useDINS <= kType + 1) {
1575        CbcHeuristicDINS heuristic5a(*model);
1576        heuristic5a.setHeuristicName("DINS");
1577        heuristic5a.setFractionSmall(0.6);
1578        if (useDINS < 4)
1579            heuristic5a.setDecayFactor(5.0);
1580        else
1581            heuristic5a.setDecayFactor(1.5);
1582        heuristic5a.setNumberNodes(1000);
1583        model->addHeuristic(&heuristic5a) ;
1584        anyToDo = true;
1585    }
1586    if (useRINS >= kType && useRINS <= kType + 1) {
1587        CbcHeuristicRINS heuristic5(*model);
1588        heuristic5.setHeuristicName("RINS");
1589        if (useRINS < 4) {
1590            heuristic5.setFractionSmall(0.5);
1591            heuristic5.setDecayFactor(5.0);
1592        } else {
1593            heuristic5.setFractionSmall(0.6);
1594            heuristic5.setDecayFactor(1.5);
1595        }
1596        model->addHeuristic(&heuristic5) ;
1597        anyToDo = true;
1598    }
1599    if (useDW >= kType && useDW <= kType + 1) {
1600        CbcHeuristicDW heuristic13(*model);
1601        heuristic13.setHeuristicName("Dantzig-Wolfe");
1602        heuristic13.setNumberPasses(100);
1603        int numberIntegers=0;
1604        const OsiSolverInterface * solver = model->solver();
1605        int numberColumns = solver->getNumCols();
1606        for (int i=0;i<numberColumns;i++) {
1607          if(solver->isInteger(i))
1608            numberIntegers++;
1609        }
1610        heuristic13.setNumberNeeded(CoinMin(200,numberIntegers/10));
1611        model->addHeuristic(&heuristic13);
1612        anyToDo = true;
1613    }
1614    if (useCombine >= kType && useCombine <= kType + 1) {
1615        CbcHeuristicLocal heuristic2(*model);
1616        heuristic2.setHeuristicName("combine solutions");
1617        heuristic2.setFractionSmall(0.5);
1618        heuristic2.setSearchType(1);
1619        model->addHeuristic(&heuristic2);
1620        anyToDo = true;
1621    }
1622    if ((useProximity >= kType && useProximity <= kType + 1)||
1623        (kType == 1 && useProximity >= 4)) {
1624        CbcHeuristicProximity heuristic2a(*model);
1625        heuristic2a.setHeuristicName("Proximity Search");
1626        heuristic2a.setFractionSmall(9999999.0);
1627        heuristic2a.setNumberNodes(30);
1628        heuristic2a.setFeasibilityPumpOptions(-2);
1629        if (useProximity>=4) {
1630          const int nodes[]={10,100,300};
1631          heuristic2a.setNumberNodes(nodes[useProximity-4]);
1632          // more print out and stronger feasibility pump
1633          if (useProximity==6)
1634            heuristic2a.setFeasibilityPumpOptions(-3);
1635        }
1636        model->addHeuristic(&heuristic2a);
1637        anyToDo = true;
1638    }
1639    if (useCrossover >= kType && useCrossover <= kType + 1) {
1640        CbcHeuristicCrossover heuristic2a(*model);
1641        heuristic2a.setHeuristicName("crossover");
1642        heuristic2a.setFractionSmall(0.3);
1643        // just fix at lower
1644        heuristic2a.setWhen(11);
1645        model->addHeuristic(&heuristic2a);
1646        model->setMaximumSavedSolutions(5);
1647        anyToDo = true;
1648    }
1649    int heurSwitches = parameters_[whichParam(CBC_PARAM_INT_HOPTIONS, numberParameters_, parameters_)].intValue() % 100;
1650    if (heurSwitches) {
1651        for (int iHeur = 0; iHeur < model->numberHeuristics(); iHeur++) {
1652            CbcHeuristic * heuristic = model->heuristic(iHeur);
1653            heuristic->setSwitches(heurSwitches);
1654        }
1655    }
1656    if (type == 2 && anyToDo) {
1657        // Do heuristics
1658#ifndef JJF_ONE
1659        // clean copy
1660        CbcModel model2(*model);
1661        // But get rid of heuristics in model
1662        model->doHeuristicsAtRoot(2);
1663        if (logLevel <= 1)
1664            model2.solver()->setHintParam(OsiDoReducePrint, true, OsiHintTry);
1665        OsiBabSolver defaultC;
1666        //solver_->setAuxiliaryInfo(&defaultC);
1667        model2.passInSolverCharacteristics(&defaultC);
1668        // Save bounds
1669        int numberColumns = model2.solver()->getNumCols();
1670        model2.createContinuousSolver();
1671        bool cleanModel = !model2.numberIntegers() && !model2.numberObjects();
1672        model2.findIntegers(false);
1673        int heurOptions = (parameters_[whichParam(CBC_PARAM_INT_HOPTIONS, numberParameters_, parameters_)].intValue() / 100) % 100;
1674        if (heurOptions == 0 || heurOptions == 2) {
1675            model2.doHeuristicsAtRoot(1);
1676        } else if (heurOptions == 1 || heurOptions == 3) {
1677            model2.setMaximumNodes(-1);
1678            CbcStrategyDefault strategy(0, 5, 5);
1679            strategy.setupPreProcessing(1, 0);
1680            model2.setStrategy(strategy);
1681            model2.branchAndBound();
1682        }
1683        if (cleanModel)
1684            model2.zapIntegerInformation(false);
1685        if (model2.bestSolution()) {
1686            double value = model2.getMinimizationObjValue();
1687            model->setCutoff(value);
1688            model->setBestSolution(model2.bestSolution(), numberColumns, value);
1689            model->setSolutionCount(1);
1690            model->setNumberHeuristicSolutions(1);
1691        }
1692#else
1693        if (logLevel <= 1)
1694            model->solver()->setHintParam(OsiDoReducePrint, true, OsiHintTry);
1695        OsiBabSolver defaultC;
1696        //solver_->setAuxiliaryInfo(&defaultC);
1697        model->passInSolverCharacteristics(&defaultC);
1698        // Save bounds
1699        int numberColumns = model->solver()->getNumCols();
1700        model->createContinuousSolver();
1701        bool cleanModel = !model->numberIntegers() && !model->numberObjects();
1702        model->findIntegers(false);
1703        model->doHeuristicsAtRoot(1);
1704        if (cleanModel)
1705            model->zapIntegerInformation(false);
1706#endif
1707        return 0;
1708    } else {
1709        return 0;
1710    }
1711}
1712
Note: See TracBrowser for help on using the repository browser.