Ignore:
Timestamp:
Sep 15, 2006 4:55:05 PM (13 years ago)
Author:
forrest
Message:

many changes

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/devel/Cbc/src/CoinSolve.cpp

    r417 r424  
    6666#include "CglTwomir.hpp"
    6767#include "CglDuplicateRow.hpp"
     68#include "CglStored.hpp"
     69#include "CglLandP.hpp"
    6870
    6971#include "CbcModel.hpp"
     
    8991                        std::string & check);
    9092static void generateCode(const char * fileName,int type,int preProcess);
    91 //#############################################################################
    92 
    9393#ifdef NDEBUG
    9494#undef NDEBUG
    9595#endif
     96//#############################################################################
     97//  Start any fake main program
     98//  End any fake main program
     99//#############################################################################
     100
    96101// Allow for interrupts
    97102// But is this threadsafe ? (so switched off by option)
     
    110115
    111116int mainTest (int argc, const char *argv[],int algorithm,
    112               ClpSimplex empty, bool doPresolve,int switchOff);
     117              ClpSimplex empty, bool doPresolve,int switchOff,bool doVector);
    113118void CbcClpUnitTest (const CbcModel & saveModel);
    114119int CbcOrClpRead_mode=1;
     
    462467  return numberDrop;
    463468}
     469void checkSOS(CbcModel * babModel, const OsiSolverInterface * solver)
     470{
     471#ifdef COIN_DEVELOP
     472  //const double *objective = solver->getObjCoefficients() ;
     473  const double *columnLower = solver->getColLower() ;
     474  const double * columnUpper = solver->getColUpper() ;
     475  const double * solution = solver->getColSolution();
     476  //int numberColumns = solver->getNumCols() ;
     477  //int numberRows = solver->getNumRows();
     478  //double direction = solver->getObjSense();
     479  //int iRow,iColumn;
     480
     481  // Row copy
     482  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
     483  //const double * elementByRow = matrixByRow.getElements();
     484  //const int * column = matrixByRow.getIndices();
     485  //const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
     486  const int * rowLength = matrixByRow.getVectorLengths();
     487
     488  // Column copy
     489  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
     490  const double * element = matrixByCol.getElements();
     491  const int * row = matrixByCol.getIndices();
     492  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
     493  const int * columnLength = matrixByCol.getVectorLengths();
     494
     495  const double * rowLower = solver->getRowLower();
     496  const double * rowUpper = solver->getRowUpper();
     497  CbcObject ** objects = babModel->objects();
     498  int numberObjects = babModel->numberObjects();
     499  for (int iObj = 0;iObj<numberObjects;iObj++) {
     500    CbcSOS * objSOS =
     501      dynamic_cast <CbcSOS *>(objects[iObj]) ;
     502    if (objSOS) {
     503      int n=objSOS->numberMembers();
     504      const int * which = objSOS->members();
     505      const double * weight = objSOS->weights();
     506      int type = objSOS->sosType();
     507      // convexity row?
     508      int iColumn;
     509      iColumn=which[0];
     510      int j;
     511      int convex=-1;
     512      for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
     513        int iRow = row[j];
     514        double value = element[j];
     515        if (rowLower[iRow]==1.0&&rowUpper[iRow]==1.0&&
     516            value==1.0) {
     517          // possible
     518          if (rowLength[iRow]==n) {
     519            if (convex==-1)
     520              convex=iRow;
     521            else
     522              convex=-2;
     523          }
     524        }
     525      }
     526      printf ("set %d of type %d has %d members - possible convexity row %d\n",
     527              iObj,type,n,convex);
     528      for (int i=0;i<n;i++) {
     529        iColumn = which[i];
     530        int convex2=-1;
     531        for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
     532          int iRow = row[j];
     533          if (iRow==convex) {
     534            double value = element[j];
     535            if (value==1.0) {
     536              convex2=iRow;
     537            }
     538          }
     539        }
     540        if (convex2<0&&convex>=0) {
     541          printf("odd convexity row\n");
     542          convex=-2;
     543        }
     544        printf("col %d has weight %g and value %g, bounds %g %g\n",
     545               iColumn,weight[i],solution[iColumn],columnLower[iColumn],
     546               columnUpper[iColumn]);
     547      }
     548    }
     549  }
     550#endif
     551}
    464552int main (int argc, const char *argv[])
    465553{
     
    480568    model.setNumberBeforeTrust(21);
    481569    int cutPass=-1234567;
     570    int tunePreProcess=5;
    482571    OsiSolverInterface * solver = model.solver();
    483572    OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
     
    500589#ifdef COIN_HAS_ASL
    501590    ampl_info info;
     591    memset(&info,0,sizeof(info));
    502592    if (argc>2&&!strcmp(argv[2],"-AMPL")) {
    503593      usingAmpl=true;
     
    570660    int dualize=0;
    571661    int doCrash=0;
     662    int doVector=0;
    572663    int doSprint=-1;
    573664    int doScaling=1;
    574665    // set reasonable defaults
    575666    int preSolve=5;
    576     int preProcess=4;
     667    int preProcess=1;
    577668    bool useStrategy=false;
    578669    bool preSolveFile=false;
     
    645736    parameters[whichParam(INCREMENT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcCutoffIncrement));
    646737    // Set up likely cut generators and defaults
    647     parameters[whichParam(PREPROCESS,numberParameters,parameters)].setCurrentOption("sos");
     738    parameters[whichParam(PREPROCESS,numberParameters,parameters)].setCurrentOption("on");
    648739    parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(128|64|1);
    649740    parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].setIntValue(-1);
     
    676767    probingGen.setRowCuts(3);
    677768    // set default action (0=off,1=on,2=root)
    678     int probingAction=3;
     769    int probingAction=1;
    679770    parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
    680771
     
    715806    int twomirAction=2;
    716807    parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption("root");
     808    CglLandP landpGen;
     809    // set default action (0=off,1=on,2=root)
     810    int landpAction=0;
     811    parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption("off");
     812    // Stored cuts
     813    bool storedCuts = false;
    717814
    718815    bool useRounding=true;
     
    815912        if (type==BAB&&goodModel) {
    816913          // check if any integers
    817           if (!lpSolver->integerInformation())
     914#ifdef COIN_HAS_ASL
     915          if (info.numberSos&&doSOS&&usingAmpl) {
     916            // SOS
     917            numberSOS = info.numberSos;
     918          }
     919#endif
     920          if (!lpSolver->integerInformation()&&!numberSOS&&
     921              !clpSolver->numberSOS())
    818922            type=DUALSIMPLEX;
    819923        }
     
    9951099              else if (parameters[iParam].type()==CUTPASS)
    9961100                cutPass = value;
     1101              else if (parameters[iParam].type()==PROCESSTUNE)
     1102                tunePreProcess = value;
    9971103              else if (parameters[iParam].type()==VERBOSE)
    9981104                verbose = value;
     
    11251231              doCrash=action;
    11261232              break;
     1233            case VECTOR:
     1234              doVector=action;
     1235              break;
    11271236            case MESSAGES:
    11281237              lpSolver->messageHandler()->setPrefix(action!=0);
     
    11781287              twomirAction = action;
    11791288              break;
     1289            case LANDPCUTS:
     1290              defaultSettings=false; // user knows what she is doing
     1291              landpAction = action;
     1292              break;
    11801293            case ROUNDING:
    11811294              defaultSettings=false; // user knows what she is doing
     
    11941307              mixedAction = action;
    11951308              twomirAction = action;
     1309              landpAction = action;
    11961310              parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption(action);
    11971311              parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption(action);
     
    12051319              parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption(action);
    12061320              parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption(action);
     1321              if (!action) {
     1322                landpAction = action;
     1323                parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption(action);
     1324              }
    12071325              break;
    12081326            case HEURISTICSTRATEGY:
     
    13471465              solveOptions.setSpecialOption(4,presolveOptions);
    13481466              solveOptions.setSpecialOption(5,printOptions);
     1467              if (doVector) {
     1468                ClpMatrixBase * matrix = lpSolver->clpMatrix();
     1469                if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
     1470                  ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
     1471                  clpMatrix->makeSpecialColumnCopy();
     1472                }
     1473              }
    13491474              if (method==ClpSolve::useDual) {
    13501475                // dual
     
    15371662          case OUTDUPROWS:
    15381663            if (goodModel) {
    1539               int nOut = outDupRow(clpSolver);
     1664              int numberRows = clpSolver->getNumRows();
     1665              //int nOut = outDupRow(clpSolver);
     1666              CglDuplicateRow dupcuts(clpSolver);
     1667              storedCuts = dupcuts.outDuplicates(clpSolver)!=0;
     1668              int nOut = numberRows-clpSolver->getNumRows();
    15401669              if (nOut&&!noPrinting)
    15411670                printf("%d rows eliminated\n",nOut);
     
    15661695            break;
    15671696          case MIPLIB:
    1568             // User can set options - main differenec is lack of model and CglPreProcess
     1697            // User can set options - main difference is lack of model and CglPreProcess
    15691698            goodModel=true;
    15701699/*
     
    15951724                OsiClpSolverInterface * si =
    15961725                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
    1597                 if (si->getModelPtr()->tightenPrimalBounds()!=0) {
     1726                ClpSimplex * clpSolver = si->getModelPtr();
     1727                if (clpSolver->tightenPrimalBounds()!=0) {
    15981728                  std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
    15991729                  exit(1);
    16001730                }
    1601                 si->getModelPtr()->dual();  // clean up
     1731                if (clpSolver->dualBound()==1.0e10) {
     1732                  // user did not set - so modify
     1733                  // get largest scaled away from bound
     1734                  double largest=1.0e-12;
     1735                  int numberRows = clpSolver->numberRows();
     1736                  const double * rowPrimal = clpSolver->primalRowSolution();
     1737                  const double * rowLower = clpSolver->rowLower();
     1738                  const double * rowUpper = clpSolver->rowUpper();
     1739                  const double * rowScale = clpSolver->rowScale();
     1740                  int iRow;
     1741                  for (iRow=0;iRow<numberRows;iRow++) {
     1742                    double value = rowPrimal[iRow];
     1743                    double above = value-rowLower[iRow];
     1744                    double below = rowUpper[iRow]-value;
     1745                    if (rowScale) {
     1746                      double multiplier = rowScale[iRow];
     1747                      above *= multiplier;
     1748                      below *= multiplier;
     1749                    }
     1750                    if (above<1.0e12)
     1751                      largest = CoinMax(largest,above);
     1752                    if (below<1.0e12)
     1753                      largest = CoinMax(largest,below);
     1754                  }
     1755                 
     1756                  int numberColumns = clpSolver->numberColumns();
     1757                  const double * columnPrimal = clpSolver->primalColumnSolution();
     1758                  const double * columnLower = clpSolver->columnLower();
     1759                  const double * columnUpper = clpSolver->columnUpper();
     1760                  const double * columnScale = clpSolver->columnScale();
     1761                  int iColumn;
     1762                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
     1763                    double value = columnPrimal[iColumn];
     1764                    double above = value-columnLower[iColumn];
     1765                    double below = columnUpper[iColumn]-value;
     1766                    if (columnScale) {
     1767                      double multiplier = 1.0/columnScale[iColumn];
     1768                      above *= multiplier;
     1769                      below *= multiplier;
     1770                    }
     1771                    if (above<1.0e12)
     1772                      largest = CoinMax(largest,above);
     1773                    if (below<1.0e12)
     1774                      largest = CoinMax(largest,below);
     1775                  }
     1776                  //if (!noPrinting)
     1777                  //std::cout<<"Largest (scaled) away from bound "<<largest<<std::endl;
     1778                  clpSolver->setDualBound(CoinMin(1000.0*largest,1.00001e10));
     1779                }
     1780                clpSolver->dual();  // clean up
    16021781              }
    16031782              // If user made settings then use them
     
    16871866              double timeLeft = babModel->getMaximumSeconds();
    16881867              int numberOriginalColumns = babModel->solver()->getNumCols();
    1689               if (preProcess==6) {
     1868              if (preProcess==7) {
    16901869                // use strategy instead
    16911870                preProcess=0;
     
    16931872              }
    16941873              if (preProcess&&type==BAB) {
     1874                // See if sos from mps file
     1875                if (numberSOS==0&&clpSolver->numberSOS()&&doSOS) {
     1876                  // SOS
     1877                  numberSOS = clpSolver->numberSOS();
     1878                  const CoinSet * setInfo = clpSolver->setInfo();
     1879                  sosStart = new int [numberSOS+1];
     1880                  sosType = new char [numberSOS];
     1881                  int i;
     1882                  int nTotal=0;
     1883                  sosStart[0]=0;
     1884                  for ( i=0;i<numberSOS;i++) {
     1885                    int type = setInfo[i].setType();
     1886                    int n=setInfo[i].numberEntries();
     1887                    sosType[i]=type;
     1888                    nTotal += n;
     1889                    sosStart[i+1] = nTotal;
     1890                  }
     1891                  sosIndices = new int[nTotal];
     1892                  sosReference = new double [nTotal];
     1893                  for (i=0;i<numberSOS;i++) {
     1894                    int n=setInfo[i].numberEntries();
     1895                    const int * which = setInfo[i].which();
     1896                    const double * weights = setInfo[i].weights();
     1897                    int base = sosStart[i];
     1898                    for (int j=0;j<n;j++) {
     1899                      int k=which[j];
     1900                      sosIndices[j+base]=k;
     1901                      sosReference[j+base] = weights ? weights[j] : (double) j;
     1902                    }
     1903                  }
     1904                }
    16951905                saveSolver=babModel->solver()->clone();
    16961906                /* Do not try and produce equality cliques and
     
    17101920                  // Add in generators
    17111921                  process.addCutGenerator(&generator1);
    1712                   int translate[]={9999,0,0,-1,2,3};
     1922                  int translate[]={9999,0,0,-1,2,3,-2};
    17131923                  process.messageHandler()->setLogLevel(babModel->logLevel());
    17141924#ifdef COIN_HAS_ASL
     
    17331943                    delete [] prohibited;
    17341944                  }
    1735                   solver2 = process.preProcessNonDefault(*saveSolver,translate[preProcess],10);
     1945                  solver2 = process.preProcessNonDefault(*saveSolver,translate[preProcess],10,
     1946                                                         tunePreProcess);
    17361947                  // Tell solver we are not in Branch and Cut
    17371948                  saveSolver->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
     
    17912002                std::string problemName ;
    17922003                babModel->solver()->getStrParam(OsiProbName,problemName) ;
    1793                 //babModel->solver()->activateRowCutDebugger(problemName.c_str()) ;
     2004                babModel->solver()->activateRowCutDebugger(problemName.c_str()) ;
    17942005                twomirGen.probname_=strdup(problemName.c_str());
    17952006                // checking seems odd
     
    18032014                int * priority = new int [numberColumns];
    18042015                const double * objective = babModel->getObjCoefficients();
     2016                const double * lower = babModel->getColLower() ;
     2017                const double * upper = babModel->getColUpper() ;
     2018                const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
     2019                const int * columnLength = matrix->getVectorLengths();
    18052020                int iColumn;
    18062021                int n=0;
     
    18092024                    sort[n]=n;
    18102025                    if (useCosts==1)
    1811                       dsort[n++]=-objective[iColumn];
    1812                     else
     2026                      dsort[n++]=-fabs(objective[iColumn]);
     2027                    else if (useCosts==2)
    18132028                      dsort[n++]=iColumn;
     2029                    else if (useCosts==3)
     2030                      dsort[n++]=upper[iColumn]-lower[iColumn];
     2031                    else if (useCosts==4)
     2032                      dsort[n++]=-(upper[iColumn]-lower[iColumn]);
     2033                    else if (useCosts==5)
     2034                      dsort[n++]=-columnLength[iColumn];
    18142035                  }
    18152036                }
     
    18582079              int switches[20];
    18592080              int numberGenerators=0;
    1860               if (probingAction==1) {
    1861                 babModel->addCutGenerator(&probingGen,-1,"Probing");
     2081              int translate[6]={-100,-1,-99,-98,1,1};
     2082              if (probingAction) {
     2083                if (probingAction==5)
     2084                  probingGen.setRowCuts(-3); // strengthening etc just at root
     2085                babModel->addCutGenerator(&probingGen,translate[probingAction],"Probing");
    18622086                switches[numberGenerators++]=0;
    1863               } else if (probingAction>=2) {
    1864                 babModel->addCutGenerator(&probingGen,-101+probingAction,"Probing");
     2087              }
     2088              if (gomoryAction) {
     2089                babModel->addCutGenerator(&gomoryGen,translate[gomoryAction],"Gomory");
     2090                switches[numberGenerators++]=-1;
     2091              }
     2092              if (knapsackAction) {
     2093                babModel->addCutGenerator(&knapsackGen,translate[knapsackAction],"Knapsack");
    18652094                switches[numberGenerators++]=0;
    18662095              }
    1867               if (gomoryAction==1) {
    1868                 babModel->addCutGenerator(&gomoryGen,-1,"Gomory");
     2096              if (redsplitAction) {
     2097                babModel->addCutGenerator(&redsplitGen,translate[redsplitAction],"Reduce-and-split");
    18692098                switches[numberGenerators++]=1;
    1870               } else if (gomoryAction>=2) {
    1871                 babModel->addCutGenerator(&gomoryGen,-101+gomoryAction,"Gomory");
     2099              }
     2100              if (cliqueAction) {
     2101                babModel->addCutGenerator(&cliqueGen,translate[cliqueAction],"Clique");
     2102                switches[numberGenerators++]=0;
     2103              }
     2104              if (mixedAction) {
     2105                babModel->addCutGenerator(&mixedGen,translate[mixedAction],"MixedIntegerRounding2");
    18722106                switches[numberGenerators++]=-1;
    18732107              }
    1874               if (knapsackAction==1) {
    1875                 babModel->addCutGenerator(&knapsackGen,-1,"Knapsack");
    1876                 switches[numberGenerators++]=0;
    1877               } else if (knapsackAction>=2) {
    1878                 babModel->addCutGenerator(&knapsackGen,-101+knapsackAction,"Knapsack");
    1879                 switches[numberGenerators++]=0;
    1880               }
    1881               if (redsplitAction==1) {
    1882                 babModel->addCutGenerator(&redsplitGen,-1,"Reduce-and-split");
     2108              if (flowAction) {
     2109                babModel->addCutGenerator(&flowGen,translate[flowAction],"FlowCover");
    18832110                switches[numberGenerators++]=1;
    1884               } else if (redsplitAction>=2) {
    1885                 babModel->addCutGenerator(&redsplitGen,-101+redsplitAction,"Reduce-and-split");
     2111              }
     2112              if (twomirAction) {
     2113                babModel->addCutGenerator(&twomirGen,translate[twomirAction],"TwoMirCuts");
    18862114                switches[numberGenerators++]=1;
    18872115              }
    1888               if (cliqueAction==1) {
    1889                 babModel->addCutGenerator(&cliqueGen,-1,"Clique");
     2116              if (landpAction) {
     2117                babModel->addCutGenerator(&landpGen,translate[landpAction],"LiftAndProject");
    18902118                switches[numberGenerators++]=1;
    1891               } else if (cliqueAction>=2) {
    1892                 babModel->addCutGenerator(&cliqueGen,-101+cliqueAction,"Clique");
    1893                 switches[numberGenerators++]=-1;
    1894               }
    1895               if (mixedAction==1) {
    1896                 babModel->addCutGenerator(&mixedGen,-1,"MixedIntegerRounding2");
    1897                 switches[numberGenerators++]=1;
    1898               } else if (mixedAction>=2) {
    1899                 babModel->addCutGenerator(&mixedGen,-101+mixedAction,"MixedIntegerRounding2");
    1900                 switches[numberGenerators++]=-1;
    1901               }
    1902               if (flowAction==1) {
    1903                 babModel->addCutGenerator(&flowGen,-1,"FlowCover");
    1904                 switches[numberGenerators++]=1;
    1905               } else if (flowAction>=2) {
    1906                 babModel->addCutGenerator(&flowGen,-101+flowAction,"FlowCover");
    1907                 switches[numberGenerators++]=1;
    1908               }
    1909               if (twomirAction==1) {
    1910                 babModel->addCutGenerator(&twomirGen,-1,"TwoMirCuts");
    1911                 switches[numberGenerators++]=1;
    1912               } else if (twomirAction>=2) {
    1913                 babModel->addCutGenerator(&twomirGen,-101+twomirAction,"TwoMirCuts");
    1914                 switches[numberGenerators++]=1;
    1915               }
     2119              }
     2120              if (storedCuts)
     2121                babModel->setSpecialOptions(babModel->specialOptions()|64);
    19162122              // Say we want timings
    19172123              numberGenerators = babModel->numberCutGenerators();
     
    19942200              }
    19952201              // probably faster to use a basis to get integer solutions
    1996               babModel->setSpecialOptions(2);
     2202              babModel->setSpecialOptions(babModel->specialOptions()|2);
    19972203              currentBranchModel = babModel;
    19982204              OsiSolverInterface * strengthenedModel=NULL;
    19992205              if (type==BAB||type==MIPLIB) {
    2000                 int moreMipOptions = parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].intValue();
     2206                int moreMipOptions = parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].intValue();
    20012207                if (moreMipOptions>=0) {
    20022208                  printf("more mip options %d\n",moreMipOptions);
    2003                   babModel->setSearchStrategy(moreMipOptions);
     2209                  if (((moreMipOptions+1)%1000000)!=0)
     2210                    babModel->setSearchStrategy(moreMipOptions%1000000);
     2211                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
     2212                  // go faster stripes
     2213                  if( moreMipOptions >=999999) {
     2214                    if (osiclp) {
     2215                      int save = osiclp->specialOptions();
     2216                      osiclp->setupForRepeatedUse(2,0);
     2217                      osiclp->setSpecialOptions(save|osiclp->specialOptions());
     2218                    }
     2219                  }
    20042220                }
    20052221              }
     
    22512467                  babModel->setStrategy(strategy);
    22522468                }
     2469                checkSOS(babModel, babModel->solver());
    22532470                babModel->branchAndBound(statistics);
     2471                checkSOS(babModel, babModel->solver());
    22542472              } else if (type==MIPLIB) {
    22552473                CbcStrategyDefault strategy(true,5,5);
    2256                 // Set up pre-processing to find sos if wanted
     2474                // Set up pre-processing
     2475                int translate2[]={9999,1,1,3,2,4,5};
    22572476                if (preProcess)
    2258                   strategy.setupPreProcessing(2);
     2477                  strategy.setupPreProcessing(translate2[preProcess]);
    22592478                babModel->setStrategy(strategy);
    22602479                CbcClpUnitTest(*babModel);
     
    22782497                         <<" to "<<babModel->rootObjectiveAfterCuts()<<std::endl;
    22792498               
     2499                numberGenerators = babModel->numberCutGenerators();
    22802500                for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
    22812501                  CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
     
    23082528                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
    23092529                }
     2530                checkSOS(babModel, babModel->solver());
    23102531              }
    23112532              if (type==STRENGTHEN&&strengthenedModel)
     
    24872708              }
    24882709              if (canOpen) {
    2489                 int status =lpSolver->readMps(fileName.c_str(),
     2710                int status =clpSolver->readMps(fileName.c_str(),
    24902711                                                   keepImportNames!=0,
    24912712                                                   allowImportErrors!=0);
     
    25762797                bool deleteModel2=false;
    25772798                ClpSimplex * model2 = lpSolver;
     2799#ifdef COIN_HAS_ASL
     2800                if (info.numberSos&&doSOS&&usingAmpl) {
     2801                  // SOS
     2802                  numberSOS = info.numberSos;
     2803                  sosStart = info.sosStart;
     2804                  sosIndices = info.sosIndices;
     2805                  sosReference = info.sosReference;
     2806                  preSolve=false;
     2807                  clpSolver->setSOSData(numberSOS,info.sosType,sosStart,sosIndices,sosReference);
     2808                }
     2809#endif
    25782810                if (preSolve) {
    25792811                  ClpPresolve pinfo;
     
    25992831
    26002832                  }
     2833                  model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
     2834                  if (deleteModel2)
     2835                    delete model2;
    26012836                } else {
    26022837                  printf("Saving model on %s\n",
    26032838                           fileName.c_str());
    2604                 }
    2605 #if 0
    2606                 // Convert names
    2607                 int iRow;
    2608                 int numberRows=model2->numberRows();
    2609                 int iColumn;
    2610                 int numberColumns=model2->numberColumns();
    2611 
    2612                 char ** rowNames = NULL;
    2613                 char ** columnNames = NULL;
    2614                 if (model2->lengthNames()) {
    2615                   rowNames = new char * [numberRows];
    2616                   for (iRow=0;iRow<numberRows;iRow++) {
    2617                     rowNames[iRow] =
    2618                       strdup(model2->rowName(iRow).c_str());
    2619 #ifdef STRIPBLANKS
    2620                     char * xx = rowNames[iRow];
    2621                     int i;
    2622                     int length = strlen(xx);
    2623                     int n=0;
    2624                     for (i=0;i<length;i++) {
    2625                       if (xx[i]!=' ')
    2626                         xx[n++]=xx[i];
     2839                  if (numberSOS) {
     2840                    // Convert names
     2841                    int iRow;
     2842                    int numberRows=model2->numberRows();
     2843                    int iColumn;
     2844                    int numberColumns=model2->numberColumns();
     2845                   
     2846                    char ** rowNames = NULL;
     2847                    char ** columnNames = NULL;
     2848                    if (model2->lengthNames()) {
     2849                      rowNames = new char * [numberRows];
     2850                      for (iRow=0;iRow<numberRows;iRow++) {
     2851                        rowNames[iRow] =
     2852                          strdup(model2->rowName(iRow).c_str());
     2853                      }
     2854                     
     2855                      columnNames = new char * [numberColumns];
     2856                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
     2857                        columnNames[iColumn] =
     2858                          strdup(model2->columnName(iColumn).c_str());
     2859                      }
    26272860                    }
    2628                     xx[n]='\0';
    2629 #endif
    2630                   }
    2631                  
    2632                   columnNames = new char * [numberColumns];
    2633                   for (iColumn=0;iColumn<numberColumns;iColumn++) {
    2634                     columnNames[iColumn] =
    2635                       strdup(model2->columnName(iColumn).c_str());
    2636 #ifdef STRIPBLANKS
    2637                     char * xx = columnNames[iColumn];
    2638                     int i;
    2639                     int length = strlen(xx);
    2640                     int n=0;
    2641                     for (i=0;i<length;i++) {
    2642                       if (xx[i]!=' ')
    2643                         xx[n++]=xx[i];
     2861                    clpSolver->writeMpsNative(fileName.c_str(),(const char **) rowNames,(const char **) columnNames,
     2862                                              (outputFormat-1)/2,1+((outputFormat-1)&1));
     2863                    if (rowNames) {
     2864                      for (iRow=0;iRow<numberRows;iRow++) {
     2865                        free(rowNames[iRow]);
     2866                      }
     2867                      delete [] rowNames;
     2868                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
     2869                        free(columnNames[iColumn]);
     2870                      }
     2871                      delete [] columnNames;
    26442872                    }
    2645                     xx[n]='\0';
    2646 #endif
     2873                  } else {
     2874                    model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
    26472875                  }
    26482876                }
    2649                 CoinMpsIO writer;
    2650                 writer.setMpsData(*model2->matrix(), COIN_DBL_MAX,
    2651                                   model2->getColLower(), model2->getColUpper(),
    2652                                   model2->getObjCoefficients(),
    2653                                   (const char*) 0 /*integrality*/,
    2654                                   model2->getRowLower(), model2->getRowUpper(),
    2655                                   columnNames, rowNames);
    2656                 // Pass in array saying if each variable integer
    2657                 writer.copyInIntegerInformation(model2->integerInformation());
    2658                 writer.setObjectiveOffset(model2->objectiveOffset());
    2659                 writer.writeMps(fileName.c_str(),0,1,1);
    2660                 if (rowNames) {
    2661                   for (iRow=0;iRow<numberRows;iRow++) {
    2662                     free(rowNames[iRow]);
    2663                   }
    2664                   delete [] rowNames;
    2665                   for (iColumn=0;iColumn<numberColumns;iColumn++) {
    2666                     free(columnNames[iColumn]);
    2667                   }
    2668                   delete [] columnNames;
    2669                 }
    2670 #else
    2671                 model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
    2672 #endif
    2673                 if (deleteModel2)
    2674                   delete model2;
    26752877                time2 = CoinCpuTime();
    26762878                totalTime += time2-time1;
     
    33353537              lpSolver->setSpecialOptions(0);
    33363538              mainTest(nFields,fields,algorithm,*lpSolver,
    3337                        (preSolve!=0),specialOptions);
     3539                       (preSolve!=0),specialOptions,doVector!=0);
    33383540            }
    33393541            break;
     
    33503552              }
    33513553              mainTest(nFields,fields,false,*lpSolver,(preSolve!=0),
    3352                        false);
     3554                       false,doVector!=0);
    33533555            }
    33543556            break;
     
    34093611            // Replace the sample code by whatever you want
    34103612            if (goodModel) {
     3613#if 1
    34113614              printf("Dummy user cbc code - model has %d rows and %d columns\n",
    34123615                     model.getNumRows(),model.getNumCols());
    3413   // Reduce printout
    3414   //solver1.setHintParam(OsiDoReducePrint,true,OsiHintTry);
    3415   OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (model.solver());
    3416   // go faster stripes
    3417   if (osiclp&&0) {
    3418     // Turn this off if you get problems
    3419     // Used to be automatically set
    3420     osiclp->setSpecialOptions(128);
    3421     if(osiclp->getNumRows()<300&&osiclp->getNumCols()<500) {
    3422       //osiclp->setupForRepeatedUse(2,0);
    3423       osiclp->setupForRepeatedUse(0,0);
    3424     }
    3425   }
    3426   // Allow rounding heuristic
    3427 
    3428   CbcRounding heuristic1(model);
    3429   model.addHeuristic(&heuristic1);
    3430 
    3431   // Do initial solve to continuous
    3432   ClpPrimalColumnSteepest steepest(5);
    3433   osiclp->getModelPtr()->setPrimalColumnPivotAlgorithm(steepest);
    3434   osiclp->getModelPtr()->setPerturbation(50);
    3435   osiclp->getModelPtr()->setInfeasibilityCost(1.0e9);
    3436   osiclp->setHintParam(OsiDoPresolveInInitial,false,OsiHintTry);
    3437   osiclp->setHintParam(OsiDoDualInInitial,false,OsiHintTry);
    3438   //osiclp->setHintParam(OsiDoDualInResolve,false,OsiHintTry);
    3439   model.setSpecialOptions(model.specialOptions()|4);
    3440   osiclp->getModelPtr()->defaultFactorizationFrequency();
    3441   {
    3442     ClpSimplex * clp = osiclp->getModelPtr();
    3443     // fix integers to 1
    3444     int numberColumns = clp->numberColumns();
    3445     double * lower = clp->columnLower();
    3446     int i;
    3447     for (i=0;i<numberColumns;i++) {
    3448       if (osiclp->isInteger(i))
    3449         lower[i]=1.0;
    3450     }
    3451     clp->primal();
    3452     double objValue = clp->objectiveValue();
    3453     osiclp->setDblParam(OsiDualObjectiveLimit,objValue+1.0e-4);
    3454     // unfix integers
    3455     for (i=0;i<numberColumns;i++) {
    3456       if (osiclp->isInteger(i))
    3457         lower[i]=0.0;
    3458     }
    3459     clp->primal();
    3460     //clp->dual();
    3461     int nArt=0;
    3462     int nFixed0=0,nFixed1=0;
    3463     double gap=objValue-clp->objectiveValue();
    3464     // for those at one fix anyway
    3465     double gap2=1.0;
    3466     const double * solution = clp->primalColumnSolution();
    3467     const double * dj = clp->dualColumnSolution();
    3468     const double * objective = clp->objective();
    3469     double * upper = clp->columnUpper();
    3470     for (i=0;i<numberColumns;i++) {
    3471       if (objective[i]>1.0e5&&solution[i]>1.0e-8)
    3472         nArt++;
    3473       if (osiclp->isInteger(i)) {
    3474         if(dj[i]>gap+1.0e-4) {
    3475           nFixed0++;
    3476           upper[i]=0.0;
    3477         }
    3478         if(-dj[i]>gap2+1.0e-4) {
    3479           nFixed1++;
    3480         lower[i]=1.0;
    3481         }
    3482       }
    3483     }
    3484     printf("%d artificials, %d fixed to 0, %d fixed to 1\n",nArt,nFixed0,nFixed1);
    3485     //osiclp->getModelPtr()->setPerturbation(100);
    3486     osiclp->setWarmStart(NULL); // set basis in osiclp
    3487   }
    3488   osiclp->initialSolve();
    3489 
    3490   // Switch off strong branching if wanted
    3491   // model.setNumberStrong(0);
    3492   // Do more strong branching if small
    3493   model.setNumberStrong(0);
    3494   model.setNumberBeforeTrust(0);
    3495 
    3496   // TEMP - set gap - better to see if no improvement in last few nodes
    3497   model.setAllowableGap(600.0);
    3498   // message levels
    3499   model.messageHandler()->setLogLevel(2);
    3500   model.solver()->messageHandler()->setLogLevel(2);
    3501   // Do complete search
    3502  
    3503   model.branchAndBound();
     3616              // Reduce printout
     3617              model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
     3618              // Do complete search
     3619              model.branchAndBound();
     3620              double objectiveValue=model.getMinimizationObjValue();
     3621              int iStat = model.status();
     3622              int iStat2 = model.secondaryStatus();
     3623#else
     3624              // Way of using an existing piece of code
     3625              OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model.solver());
     3626              ClpSimplex * lpSolver = clpSolver->getModelPtr();
     3627              // set time from integer model
     3628              double timeToGo = model.getMaximumSeconds();
     3629              lpSolver->setMaximumSeconds(timeToGo);
     3630              fakeMain(*lpSolver,*clpSolver,model);
     3631              // My actual usage has objective only in clpSolver
     3632              double objectiveValue=clpSolver->getObjValue();
     3633              int iStat = lpSolver->status();
     3634              int iStat2 = lpSolver->secondaryStatus();
     3635#endif
     3636              // make sure solution back in correct place
     3637              clpSolver = dynamic_cast< OsiClpSolverInterface*> (model.solver());
     3638              lpSolver = clpSolver->getModelPtr();
     3639#ifdef COIN_HAS_ASL
     3640              if (usingAmpl) {
     3641                int n = clpSolver->getNumCols();
     3642                double value = objectiveValue*lpSolver->getObjSense();
     3643                char buf[300];
     3644                int pos=0;
     3645                std::string minor[]={"","","gap","nodes","time","","solutions","user ctrl-c"};
     3646                if (iStat==0) {
     3647                  if (objectiveValue<1.0e40) {
     3648                    pos += sprintf(buf+pos,"optimal," );
     3649                  } else {
     3650                    // infeasible
     3651                    iStat=1;
     3652                    pos += sprintf(buf+pos,"infeasible,");
     3653                  }
     3654                } else if (iStat==1) {
     3655                  if (iStat2!=6)
     3656                    iStat=3;
     3657                  else
     3658                    iStat=4;
     3659                  pos += sprintf(buf+pos,"stopped on %s,",minor[iStat2].c_str());
     3660                } else if (iStat==2) {
     3661                  iStat = 7;
     3662                  pos += sprintf(buf+pos,"stopped on difficulties,");
     3663                } else if (iStat==5) {
     3664                  iStat = 3;
     3665                  pos += sprintf(buf+pos,"stopped on ctrl-c,");
     3666                } else {
     3667                  pos += sprintf(buf+pos,"status unknown,");
     3668                  iStat=6;
     3669                }
     3670                info.problemStatus=iStat;
     3671                info.objValue = value;
     3672                if (objectiveValue<1.0e40)
     3673                  pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
     3674                                 value);
     3675                sprintf(buf+pos,"\n%d nodes, %d iterations",
     3676                        model.getNodeCount(),
     3677                        model.getIterationCount());
     3678                if (objectiveValue<1.0e50) {
     3679                  free(info.primalSolution);
     3680                  info.primalSolution = (double *) malloc(n*sizeof(double));
     3681                  CoinCopyN(lpSolver->primalColumnSolution(),n,info.primalSolution);
     3682                  int numberRows = lpSolver->numberRows();
     3683                  free(info.dualSolution);
     3684                  info.dualSolution = (double *) malloc(numberRows*sizeof(double));
     3685                  CoinCopyN(lpSolver->dualRowSolution(),numberRows,info.dualSolution);
     3686                } else {
     3687                  info.primalSolution=NULL;
     3688                  info.dualSolution=NULL;
     3689                }
     3690                // put buffer into info
     3691                strcpy(info.buffer,buf);
     3692#endif
     3693              }
    35043694            }
    35053695            break;
     
    42764466    }
    42774467  } else {
    4278     int translate[]={9999,0,0,-1,2,3};
     4468    int translate[]={9999,0,0,-1,2,3,-2};
    42794469    strcpy(line[numberLines++],"5  // Hand coded preprocessing");
    42804470    strcpy(line[numberLines++],"5  CglPreProcess process;");
Note: See TracChangeset for help on using the changeset viewer.