Ignore:
Timestamp:
Nov 28, 2009 6:09:56 AM (10 years ago)
Author:
forrest
Message:

final changes before cleaning

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/sandbox/Cbc/src/CbcModel.cpp

    r1286 r1315  
    4444#include "CbcHeuristic.hpp"
    4545#include "CbcHeuristicFPump.hpp"
     46#include "CbcHeuristicRINS.hpp"
    4647#include "CbcHeuristicDive.hpp"
    4748#include "CbcModel.hpp"
     
    11711172    dblParam_[CbcLargestChange] = 0.0;
    11721173    intParam_[CbcNumberBranches] = 0;
     1174    // Double check optimization directions line up
     1175    dblParam_[CbcOptimizationDirection] = solver_->getObjSense();
    11731176    strongInfo_[0] = 0;
    11741177    strongInfo_[1] = 0;
     
    18631866    maximumDepthActual_ = 0;
    18641867    numberDJFixed_ = 0.0;
     1868    if (!parentModel_) {
     1869        if ((specialOptions_&262144) != 0) {
     1870            // create empty stored cuts
     1871            //storedRowCuts_ = new CglStored(solver_->getNumCols());
     1872        } else if ((specialOptions_&524288) != 0 && storedRowCuts_) {
     1873            // tighten and set best solution
     1874            // A) tight bounds on integer variables
     1875            const double * lower = solver_->getColLower();
     1876            const double * upper = solver_->getColUpper();
     1877            const double * tightLower = storedRowCuts_->tightLower();
     1878            const double * tightUpper = storedRowCuts_->tightUpper();
     1879            int nTightened = 0;
     1880            for (int i = 0; i < numberIntegers_; i++) {
     1881                int iColumn = integerVariable_[i];
     1882                if (tightLower[iColumn] > lower[iColumn]) {
     1883                    nTightened++;
     1884                    solver_->setColLower(iColumn, tightLower[iColumn]);
     1885                }
     1886                if (tightUpper[iColumn] < upper[iColumn]) {
     1887                    nTightened++;
     1888                    solver_->setColUpper(iColumn, tightUpper[iColumn]);
     1889                }
     1890            }
     1891            if (nTightened)
     1892                printf("%d tightened by alternate cuts\n", nTightened);
     1893            if (storedRowCuts_->bestObjective() < bestObjective_) {
     1894                // B) best solution
     1895                double objValue = storedRowCuts_->bestObjective();
     1896                setBestSolution(CBC_SOLUTION, objValue,
     1897                                storedRowCuts_->bestSolution()) ;
     1898                // Do heuristics
     1899                // Allow RINS
     1900                for (int i = 0; i < numberHeuristics_; i++) {
     1901                    CbcHeuristicRINS * rins
     1902                    = dynamic_cast<CbcHeuristicRINS *> (heuristic_[i]);
     1903                    if (rins) {
     1904                        rins->setLastNode(-100);
     1905                    }
     1906                }
     1907            }
     1908        }
     1909    }
    18651910    // Do heuristics
    18661911    doHeuristicsAtRoot();
     
    22152260            if (numberUnsatisfied)   {
    22162261                // User event
    2217                 if (!eventHappened_ && feasible)
     2262                if (!eventHappened_ && feasible) {
    22182263                    feasible = solveWithCuts(cuts, maximumCutPassesAtRoot_,
    22192264                                             NULL);
    2220                 else
     2265                    if ((specialOptions_&524288) != 0 && !parentModel_
     2266                            && storedRowCuts_) {
     2267                        if (feasible) {
     2268                            /* pick up stuff and try again
     2269                            add cuts, maybe keep around
     2270                            do best solution and if so new heuristics
     2271                            obviously tighten bounds
     2272                            */
     2273                            // A and B probably done on entry
     2274                            // A) tight bounds on integer variables
     2275                            const double * lower = solver_->getColLower();
     2276                            const double * upper = solver_->getColUpper();
     2277                            const double * tightLower = storedRowCuts_->tightLower();
     2278                            const double * tightUpper = storedRowCuts_->tightUpper();
     2279                            int nTightened = 0;
     2280                            for (int i = 0; i < numberIntegers_; i++) {
     2281                                int iColumn = integerVariable_[i];
     2282                                if (tightLower[iColumn] > lower[iColumn]) {
     2283                                    nTightened++;
     2284                                    solver_->setColLower(iColumn, tightLower[iColumn]);
     2285                                }
     2286                                if (tightUpper[iColumn] < upper[iColumn]) {
     2287                                    nTightened++;
     2288                                    solver_->setColUpper(iColumn, tightUpper[iColumn]);
     2289                                }
     2290                            }
     2291                            if (nTightened)
     2292                                printf("%d tightened by alternate cuts\n", nTightened);
     2293                            if (storedRowCuts_->bestObjective() < bestObjective_) {
     2294                                // B) best solution
     2295                                double objValue = storedRowCuts_->bestObjective();
     2296                                setBestSolution(CBC_SOLUTION, objValue,
     2297                                                storedRowCuts_->bestSolution()) ;
     2298                                // Do heuristics
     2299                                // Allow RINS
     2300                                for (int i = 0; i < numberHeuristics_; i++) {
     2301                                    CbcHeuristicRINS * rins
     2302                                    = dynamic_cast<CbcHeuristicRINS *> (heuristic_[i]);
     2303                                    if (rins) {
     2304                                        rins->setLastNode(-100);
     2305                                    }
     2306                                }
     2307                                doHeuristicsAtRoot();
     2308                            }
     2309#if 0
     2310                            int nCuts = storedRowCuts_->sizeRowCuts();
     2311                            // add to global list
     2312                            for (int i = 0; i < nCuts; i++) {
     2313                                OsiRowCut newCut(*storedRowCuts_->rowCutPointer(i));
     2314                                newCut.setGloballyValidAsInteger(2);
     2315                                newCut.mutableRow().setTestForDuplicateIndex(false);
     2316                                globalCuts_.insert(newCut) ;
     2317                            }
     2318#else
     2319                            addCutGenerator(storedRowCuts_, -99, "Stored from previous run",
     2320                                            true, false, false, -200);
     2321#endif
     2322                            // Set cuts as active
     2323                            delete [] addedCuts_ ;
     2324                            maximumNumberCuts_ = cuts.sizeRowCuts();
     2325                            if (maximumNumberCuts_) {
     2326                                addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
     2327                            } else {
     2328                                addedCuts_ = NULL;
     2329                            }
     2330                            for (int i = 0; i < maximumNumberCuts_; i++)
     2331                                addedCuts_[i] = new CbcCountRowCut(*cuts.rowCutPtr(i),
     2332                                                                   NULL, -1, -1, 2);
     2333                            printf("size %d\n", cuts.sizeRowCuts());
     2334                            cuts = OsiCuts();
     2335                            currentNumberCuts_ = maximumNumberCuts_;
     2336                            feasible = solveWithCuts(cuts, maximumCutPassesAtRoot_,
     2337                                                     NULL);
     2338                            for (int i = 0; i < maximumNumberCuts_; i++)
     2339                                delete addedCuts_[i];
     2340                        }
     2341                        delete storedRowCuts_;
     2342                        storedRowCuts_ = NULL;
     2343                    }
     2344                } else {
    22212345                    feasible = false;
     2346                }
    22222347            }   else if (solverCharacteristics_->solutionAddsCuts() ||
    22232348                       solverCharacteristics_->alwaysTryCutsAtRootNode()) {
     
    25052630            if (bestSolution_)
    25062631                heuristic.setInputSolution(bestSolution_, bestObjective_);
    2507             heuristic.setFractionSmall(0.5);
     2632            heuristic.setFractionSmall(0.9);
    25082633            heuristic.setFeasibilityPumpOptions(1008013);
    25092634            // Use numberNodes to say how many are original rows
     
    25362661            }
    25372662        }
     2663    }
     2664    if ((specialOptions_&262144) != 0 && !parentModel_) {
     2665        // Save stuff and return!
     2666        storedRowCuts_->saveStuff(bestObjective_, bestSolution_,
     2667                                  solver_->getColLower(),
     2668                                  solver_->getColUpper());
     2669        delete [] lowerBefore;
     2670        delete [] upperBefore;
     2671        delete saveSolver;
     2672        return;
    25382673    }
    25392674    /*
     
    28653000        }
    28663001#endif
    2867         if (probingInfo_) {
    2868             int number01 = probingInfo_->numberIntegers();
    2869             //const fixEntry * entry = probingInfo_->fixEntries();
    2870             const int * toZero = probingInfo_->toZero();
    2871             //const int * toOne = probingInfo_->toOne();
    2872             //const int * integerVariable = probingInfo_->integerVariable();
    2873             if (toZero[number01]) {
    2874                 if (probingInfo_->packDown()) {
    2875 #ifdef CLP_INVESTIGATE
    2876                     printf("%d implications on %d 0-1\n", toZero[number01], number01);
    2877 #endif
    2878                     CglImplication implication(probingInfo_);
    2879                     addCutGenerator(&implication, 1, "ImplicationCuts", true, false, false, -200);
    2880                 } else {
    2881                     delete probingInfo_;
    2882                     probingInfo_ = NULL;
    2883                 }
    2884             } else {
    2885                 delete probingInfo_;
    2886 
    2887                 probingInfo_ = NULL;
    2888             }
    2889         }
    28903002        newNode = new CbcNode ;
    28913003        // Set objective value (not so obvious if NLP etc)
     
    29333045            newNode = NULL ;
    29343046            feasible = false ;
     3047        }
     3048    }
     3049    if (newNode && probingInfo_) {
     3050        int number01 = probingInfo_->numberIntegers();
     3051        //const fixEntry * entry = probingInfo_->fixEntries();
     3052        const int * toZero = probingInfo_->toZero();
     3053        //const int * toOne = probingInfo_->toOne();
     3054        //const int * integerVariable = probingInfo_->integerVariable();
     3055        if (toZero[number01]) {
     3056            CglTreeProbingInfo info(*probingInfo_);
     3057#if 0
     3058            OsiSolverInterface * fake = info.analyze(*solver_, 1);
     3059            if (fake) {
     3060                fake->writeMps("fake");
     3061                CglFakeClique cliqueGen(fake);
     3062                cliqueGen.setStarCliqueReport(false);
     3063                cliqueGen.setRowCliqueReport(false);
     3064                cliqueGen.setMinViolation(0.1);
     3065                addCutGenerator(&cliqueGen, 1, "Fake cliques", true, false, false, -200);
     3066                generator_[numberCutGenerators_-1]->setTiming(true);
     3067            }
     3068#endif
     3069            if (probingInfo_->packDown()) {
     3070#ifdef CLP_INVESTIGATE
     3071                printf("%d implications on %d 0-1\n", toZero[number01], number01);
     3072#endif
     3073                CglImplication implication(probingInfo_);
     3074                addCutGenerator(&implication, 1, "ImplicationCuts", true, false, false, -200);
     3075            } else {
     3076                delete probingInfo_;
     3077                probingInfo_ = NULL;
     3078            }
     3079        } else {
     3080            delete probingInfo_;
     3081
     3082            probingInfo_ = NULL;
    29353083        }
    29363084    }
     
    32043352    int lastEvery1000 = 0;
    32053353    int lastPrintEvery = 0;
     3354    int numberConsecutiveInfeasible = 0;
    32063355    while (true) {
    32073356#ifdef CBC_THREAD
     
    34663615#endif
    34673616                    numberFixed += numberFixed2;
    3468                     if (numberFixed*10 < numberColumns)
     3617                    if (numberFixed*10 < numberColumns && numberFixed*4 < numberIntegers_)
    34693618                        tryNewSearch = false;
    34703619                }
     
    34873636                    if (bestSolution_)
    34883637                        heuristic.setInputSolution(bestSolution_, bestObjective_);
    3489                     heuristic.setFractionSmall(0.6);
     3638                    heuristic.setFractionSmall(0.8);
    34903639                    heuristic.setFeasibilityPumpOptions(1008013);
    34913640                    // Use numberNodes to say how many are original rows
     
    35403689        assert(!solverCharacteristics_->solutionAddsCuts() || solverCharacteristics_->mipFeasible());
    35413690#endif
     3691#define DIVE_WHEN 1000
     3692#define DIVE_STOP 2000
     3693        int kNode = numberNodes_ % 4000;
     3694        if (numberNodes_<100000 && kNode>DIVE_WHEN && kNode <= DIVE_STOP) {
     3695            if (!parallelMode()) {
     3696                if (kNode == DIVE_WHEN + 1 || numberConsecutiveInfeasible > 1) {
     3697                    CbcCompareDefault * compare = dynamic_cast<CbcCompareDefault *>
     3698                                                  (nodeCompare_);
     3699                    if (compare) {
     3700                        //printf("Redoing tree\n");
     3701                        compare->startDive(this);
     3702                        numberConsecutiveInfeasible = 0;
     3703                    }
     3704                }
     3705            }
     3706        }
    35423707        if (cutoff > getCutoff()) {
    35433708            double newCutoff = getCutoff();
     
    37313896            // Do main work of solving node here
    37323897            doOneNode(this, node, createdNode);
     3898#if 0
     3899            if (node) {
     3900                if (createdNode) {
     3901                    printf("Node %d depth %d, created %d depth %d\n",
     3902                           node->nodeNumber(), node->depth(),
     3903                           createdNode->nodeNumber(), createdNode->depth());
     3904                } else {
     3905                    printf("Node %d depth %d,  no created node\n",
     3906                           node->nodeNumber(), node->depth());
     3907                }
     3908            } else if (createdNode) {
     3909                printf("Node exhausted, created %d depth %d\n",
     3910                       createdNode->nodeNumber(), createdNode->depth());
     3911            } else {
     3912                printf("Node exhausted,  no created node\n");
     3913                numberConsecutiveInfeasible = 2;
     3914            }
     3915#endif
     3916            //if (createdNode)
     3917            //numberConsecutiveInfeasible=0;
     3918            //else
     3919            //numberConsecutiveInfeasible++;
    37333920#ifdef CBC_THREAD
    37343921        } else if (parallelMode() > 0) {
     
    44644651{
    44654652    assert (solver_);
     4653    // Double check optimization directions line up
     4654    dblParam_[CbcOptimizationDirection] = solver_->getObjSense();
    44664655    // Check if bounds are all integral (as may get messed up later)
    44674656    checkModel();
     
    46574846        maximumNumberUpdateItems_(0),
    46584847        updateItems_(NULL),
     4848        storedRowCuts_(NULL),
    46594849        numberThreads_(0),
    46604850        threadMode_(0)
     
    48145004        maximumNumberUpdateItems_(0),
    48155005        updateItems_(NULL),
     5006        storedRowCuts_(NULL),
    48165007        numberThreads_(0),
    48175008        threadMode_(0)
     
    50555246        maximumNumberUpdateItems_(rhs.maximumNumberUpdateItems_),
    50565247        updateItems_(NULL),
     5248        storedRowCuts_(NULL),
    50575249        numberThreads_(rhs.numberThreads_),
    50585250        threadMode_(rhs.threadMode_)
     
    68026994            CbcSimpleInteger * simpleObject =
    68036995                static_cast <CbcSimpleInteger *>(object) ;
    6804             int iObject;
     6996            int iObject = simpleObject->position();
     6997#ifndef NDEBUG
    68056998            int iColumn = simpleObject->columnNumber();
    6806             for (iObject = 0 ; iObject < numberObjects_ ; iObject++) {
     6999            int jObject;
     7000            for (jObject = 0 ; jObject < numberObjects_ ; jObject++) {
    68077001                simpleObject =
    6808                     static_cast <CbcSimpleInteger *>(object_[iObject]) ;
     7002                    static_cast <CbcSimpleInteger *>(object_[jObject]) ;
    68097003                if (simpleObject->columnNumber() == iColumn)
    68107004                    break;
    68117005            }
    6812             assert (iObject < numberObjects_);
     7006            assert (jObject < numberObjects_ && iObject == jObject);
     7007#else
     7008#ifdef TIGHTEN_BOUNDS
     7009            int iColumn = simpleObject->columnNumber();
     7010#endif
     7011#endif
    68137012            update.objectNumber_ = iObject;
    68147013            addUpdateInformation(update);
     
    71967395                if (generate) {
    71977396                    if (!node && cut_obj[CUT_HISTORY-1] != -COIN_DBL_MAX &&
    7198                             fabs(cut_obj[CUT_HISTORY-1] - cut_obj[CUT_HISTORY-2]) < 1.0e-7)
     7397                            fabs(cut_obj[CUT_HISTORY-1] - cut_obj[CUT_HISTORY-2]) < 1.0e-7 +
     7398                            1.0e-6*fabs(cut_obj[CUT_HISTORY-1]))
    71997399                        generator_[i]->setIneffectualCuts(true);
    72007400                    bool mustResolve =
     
    79518151                    addCuts[i] = &theseCuts.rowCut(i) ;
    79528152                }
     8153                if ((specialOptions_&262144) != 0 && !parentModel_) {
     8154                    //save
     8155                    for (i = 0 ; i < numberToAdd ; i++)
     8156                        storedRowCuts_->addCut(*addCuts[i]);
     8157                }
    79538158                solver_->applyRowCuts(numberToAdd, addCuts) ;
    79548159                CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
     
    81028307                    }
    81038308                }
    8104                 if (numberTries == 1 && currentDepth_ && currentPassNumber_ < 10) {
    8105                     if (thisObj - lastObjective > 10.0*minimumDrop) {
     8309                if (numberTries == 1 && currentDepth_ < 12 && currentPassNumber_ < 10) {
     8310                    double drop[12] = {1.0, 2.0, 3.0, 10.0, 10.0, 10.0, 10.0, 20.0, 100.0, 100.0, 1000.0, 1000.0};
     8311                    if (thisObj - lastObjective > drop[currentDepth_]*minimumDrop) {
    81068312                        numberTries++;
    81078313#ifdef CLP_INVESTIGATE
     
    1030610512            if (!type) {
    1030710513                obj2->setNumberBeforeTrust(numberBeforeTrust_);
    10308             } else {
     10514            } else if (type == 1) {
    1030910515                int value = obj2->numberBeforeTrust();
    1031010516                value = (value * 11) / 10 + 1;
    1031110517                value = CoinMax(numberBeforeTrust_, value);
    1031210518                obj2->setNumberBeforeTrust(value);
     10519            } else {
     10520                assert (type == 2);
     10521                int value = obj2->numberBeforeTrust();
     10522                int n = CoinMax(obj2->numberTimesDown(),
     10523                                obj2->numberTimesUp());
     10524                if (n >= value)
     10525                    obj2->setNumberBeforeTrust(n + 1);
    1031310526            }
    1031410527        }
     
    1242212635    if (currentNode_ && currentNode_->depth() >= 8)
    1242312636        stateOfSearch_ += 10;
     12637    //printf("STate %d, %d nodes - parent %c - sol %d %d\n",
     12638    // stateOfSearch_,numberNodes_,parentModel_ ? 'Y' :'N',
     12639    // numberSolutions_,numberHeuristicSolutions_);
    1242412640    int anyAction = -1 ;
    1242512641    resolved = false ;
     
    1381514031        bool checkingNode = false;
    1381614032        if (feasible) {
    13817 #if 0
     14033#ifdef FUNNY_BRANCHING
    1381814034            // Far too clever
    13819             if (numberThreads_ == -10 && node->numberBranches() == 2) {
     14035            if ((numberThreads_ == -10 || true) && node->numberBranches() == 2) {
    1382014036                // see if any parent branches redundant
    1382114037                // Look at state of "node"
    1382214038                CbcNodeInfo * nodeInfo = node->nodeInfo();
    13823                 // See if any branched variables off bounds
    13824                 const double * dj = solver_->getReducedCost();
    13825                 const double * lower = solver_->getColLower();
    13826                 const double * upper = solver_->getColUpper();
    13827                 const double * solution = solver_->getColSolution();
    13828                 double direction = solver_->getObjSense() ;
    1382914039                if (nodeInfo) {
     14040                    // See if any branched variables off bounds
     14041                    const double * dj = solver_->getReducedCost();
     14042                    const double * lower = solver_->getColLower();
     14043                    const double * upper = solver_->getColUpper();
     14044                    const double * solution = solver_->getColSolution();
     14045                    int numberColumns = solver_->getNumCols();
     14046                    double * currentLower = CoinCopyOfArray(lower, numberColumns);
     14047                    double * currentUpper = CoinCopyOfArray(upper, numberColumns);
     14048                    char * touched = new char[numberColumns];
     14049                    memset(touched, 0, numberColumns);
     14050                    double direction = solver_->getObjSense() ;
    1383014051                    bool canDelete = nodeInfo->numberBranchesLeft() > 0;
    1383114052                    //int numberBounds = nodeInfo->numberChangedBounds();
     
    1384414065                            double originalLower1 = object1->originalLowerBound();
    1384514066                            double originalUpper1 = object1->originalUpperBound();
     14067                            // Unset all bounds from parents
     14068                            CbcPartialNodeInfo * partial =
     14069                                dynamic_cast<CbcPartialNodeInfo *>(nodeInfo);
     14070                            touched[iColumn1] = 1;
     14071                            if (partial) {
     14072                                /* maybe don't do if obj hasn't changed
     14073                                as then you might get loop
     14074                                at present just 0-1
     14075                                as need to know original bound
     14076                                */
     14077                                int n = partial->numberChangedBounds();
     14078                                const int * which = partial->variables();
     14079                                const double * values = partial->newBounds();
     14080                                for (int i = 0; i < n; i++) {
     14081                                    int variable = which[i];
     14082                                    int k = variable & 0x3fffffff;
     14083                                    assert (k != iColumn1);
     14084                                    if (!touched[k]) {
     14085                                        if ((variable&0x80000000) == 0) {
     14086                                            // lower bound changing
     14087                                            assert (currentLower[k] == 1.0);
     14088                                            currentLower[k] = 0.0;
     14089                                        } else {
     14090                                            // upper bound changing
     14091                                            assert (currentUpper[k] == 0.0);
     14092                                            currentUpper[k] = 1.0;
     14093                                        }
     14094                                    }
     14095                                }
     14096                            }
    1384614097                            zeroOne1 = originalLower1 == 0.0 && originalUpper1 == 1.0;
    1384714098                            way1 = objectI->way();
    1384814099                            assert (way1 == -1 || way1 == 1);
     14100                            int kWay = way1;
    1384914101                            //way1 = -way1; // what last branch did
    1385014102                            // work out using bounds
     
    1385414106                            else
    1385514107                                way1 = 1;
     14108                            assert (kWay == -way1);
    1385614109                            if (way1 < 0) {
    1385714110                                // must have been down branch
     
    1392314176                                                previousBounds(node, nodeInfo, iColumn2, newLower, newUpper, 2);
    1392414177                                                solver_->setColUpper(iColumn2, newUpper);
     14178                                                assert (newLower == lower[iColumn2]);
    1392514179                                            } else {
    1392614180                                                printf("SKipping\n");
    1392714181                                            }
    13928                                         } else if (iColumn1 >= 0 && iColumn1 != iColumn2 && (!inBetween || true) && zeroOne1 && zeroOne2) {
     14182                                        } else if (iColumn1 >= 0 && iColumn1 != iColumn2 && (!inBetween || true) && zeroOne1 && zeroOne2 && false) {
    1392914183#if 1
    1393014184                                            if (true) {
     
    1396514219                                                previousBounds(node, nodeInfo, iColumn2, newLower, newUpper, 1);
    1396614220                                                solver_->setColLower(iColumn2, newLower);
     14221                                                assert (newUpper == upper[iColumn2]);
    1396714222                                            } else {
    1396814223                                                printf("SKipping\n");
    1396914224                                            }
    13970                                         } else if (iColumn1 >= 0 && iColumn1 != iColumn2 && (!inBetween || true) && zeroOne1 && zeroOne2) {
     14225                                        } else if (iColumn1 >= 0 && iColumn1 != iColumn2 && (!inBetween || true) && zeroOne1 && zeroOne2 && false) {
    1397114226#if 1
    1397214227                                            // add in bounds
     
    1399414249                        }
    1399514250                    }
     14251                    delete [] currentLower;
     14252                    delete [] currentUpper;
    1399614253                }
    1399714254            }
     
    1520015457        }
    1520115458    }
    15202 #ifdef CLP_INVESTIGATE
     15459#ifdef CLP_INVESTIGATE4
    1520315460    if (priority)
    1520415461        printf("Before fathom %d not trusted out of %d\n",
Note: See TracChangeset for help on using the changeset viewer.