Ignore:
Timestamp:
Aug 21, 2009 12:19:13 PM (10 years ago)
Author:
forrest
Message:

fixes

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Cbc/src/unitTestClp.cpp

    r1200 r1212  
    7979
    8080//#############################################################################
    81 
     81// testSwitch -2 unitTest, -1 normal (==2)
    8282int CbcClpUnitTest (const CbcModel & saveModel, std::string& dirMiplib,
    83                     bool unitTestOnly,
     83                    int testSwitch,
    8484                    double * stuff)
    8585{
     
    107107  std::vector<double> objValueC ;
    108108  std::vector<double> objValue ;
    109   std::vector<int> strategy ;
     109  std::vector<int> testSet ;
    110110  /*
    111111    And a macro to make the vector creation marginally readable.
     
    113113#define PUSH_MPS(zz_mpsName_zz,\
    114114                 zz_nRows_zz,zz_nCols_zz,zz_objValue_zz,zz_objValueC_zz, \
    115                  zz_strategy_zz) \
     115                 zz_testSet_zz) \
    116116  mpsName.push_back(zz_mpsName_zz) ; \
    117117  nRows.push_back(zz_nRows_zz) ; \
    118118  nCols.push_back(zz_nCols_zz) ; \
    119119  objValueC.push_back(zz_objValueC_zz) ; \
    120   strategy.push_back(zz_strategy_zz) ; \
     120  testSet.push_back(zz_testSet_zz) ; \
    121121  objValue.push_back(zz_objValue_zz) ;
    122 
    123   if (unitTestOnly) {
    124     PUSH_MPS("p0033",16,33,3089,2520.57,7);
    125     PUSH_MPS("p0201",133,201,7615,6875.0,7);
     122  int loSwitch=0;
     123  if (testSwitch==-2) {
     124    PUSH_MPS("p0033",16,33,3089,2520.57,0);
     125    PUSH_MPS("p0201",133,201,7615,6875.0,0);
    126126  } else {
     127    if (testSwitch==-1) {
     128      testSwitch=1;
     129    } else {
     130      loSwitch=static_cast<int>(stuff[6]);
     131      printf("Solving miplib problems in sets >= %d and <=%d\n",
     132             loSwitch,testSwitch);
     133    }
    127134    /*
    128135      Load up the problem vector. Note that the row counts here include the
     
    130137    */
    131138    // 0 for no test, 1 for some, 2 for many, 3 for all
    132     //PUSH_MPS("blend2",274,353,7.598985,6.9156751140,7);
    133     //PUSH_MPS("p2756",755,2756,3124,2688.75,7);
    134     //PUSH_MPS("seymour_1",4944,1372,410.7637014,404.35152,7);
    135     //PUSH_MPS("enigma",21,100,0.0,0.0,7);
    136     //PUSH_MPS("misc03",96,160,3360,1910.,7);
    137     //PUSH_MPS("p0201",133,201,7615,6875.0,7);
    138 #define HOWMANY 2
     139    //PUSH_MPS("blend2",274,353,7.598985,6.9156751140,0);
     140    //PUSH_MPS("p2756",755,2756,3124,2688.75,0);
     141    //PUSH_MPS("seymour_1",4944,1372,410.7637014,404.35152,0);
     142    //PUSH_MPS("enigma",21,100,0.0,0.0,0);
     143    //PUSH_MPS("misc03",96,160,3360,1910.,0);
     144    //PUSH_MPS("p0201",133,201,7615,6875.0,0);
     145#define HOWMANY 6
    139146#if HOWMANY
    140 #if HOWMANY>1
    141     PUSH_MPS("10teams",230,2025,924,917,7);
    142 #endif
    143     PUSH_MPS("air03",124,10757,340160,338864.25,7);
    144 #if HOWMANY>2
     147    PUSH_MPS("10teams",230,2025,924,917,1);
     148    PUSH_MPS("air03",124,10757,340160,338864.25,0);
    145149    PUSH_MPS("air04",823,8904,56137,55535.436,8);
    146150    PUSH_MPS("air05",426,7195,26374,25877.609,8);
    147 #endif
    148   //    PUSH_MPS("arki001",1048,1388,7580813.0459,7579599.80787,7);
    149     PUSH_MPS("bell3a",123,133,878430.32,862578.64,7);
    150 #if HOWMANY>1
    151     PUSH_MPS("bell5",91,104,8966406.49,8608417.95,7);
    152 #endif
    153     PUSH_MPS("blend2",274,353,7.598985,6.9156751140,7);
    154 #if HOWMANY>1
    155     PUSH_MPS("cap6000",2176,6000,-2451377,-2451537.325,7);
    156 #endif
    157     //    PUSH_MPS("dano3mip",3202,13873,728.1111,576.23162474,7);
    158     //PUSH_MPS("danoint",664,521,65.67,62.637280418,7);
    159     PUSH_MPS("dcmulti",290,548,188182,183975.5397,7);
    160     PUSH_MPS("dsbmip",1182,1886,-305.19817501,-305.19817501,7);
    161     PUSH_MPS("egout",98,141,568.101,149.589,7);
    162     PUSH_MPS("enigma",21,100,0.0,0.0,7);
    163 #if HOWMANY>3
    164     PUSH_MPS("fast0507",507,63009,174,172.14556668,7);
    165 #endif
    166     PUSH_MPS("fiber",363,1298,405935.18000,156082.51759,7);
    167 #if HOWMANY>1
    168     PUSH_MPS("fixnet6",478,878,3983,1200.88,7);
    169 #endif
    170     PUSH_MPS("flugpl",18,18,1201500,1167185.7,7);
    171     PUSH_MPS("gen",780,870,112313,112130.0,7);
    172 #if HOWMANY>1
    173     PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,7);
    174     PUSH_MPS("gesa2_o",1248,1224,25779856.372,25476489.678,7);
    175 #endif
    176     PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,7);
    177     PUSH_MPS("gesa3_o",1224,1152,27991042.648,27833632.451,7);
    178     PUSH_MPS("gt2",29,188,21166.000,13460.233074,7);
    179 #if HOWMANY>2
    180     PUSH_MPS("harp2",112,2993,-73899798.00,-74353341.502,7);
    181 #endif
    182     PUSH_MPS("khb05250",101,1350,106940226,95919464.0,7);
    183 #if HOWMANY>1
    184     PUSH_MPS("l152lav",97,1989,4722,4656.36,7);
    185 #endif
    186     PUSH_MPS("lseu",28,89,1120,834.68,7);
    187     PUSH_MPS("misc03",96,160,3360,1910.,7);
    188     PUSH_MPS("misc06",820,1808,12850.8607,12841.6,7);
    189 #if HOWMANY>1
    190     PUSH_MPS("misc07",212,260,2810,1415.0,7);
    191     PUSH_MPS("mitre",2054,10724,115155,114740.5184,7);
    192 #endif
    193     PUSH_MPS("mod008",6,319,307,290.9,7);
    194     PUSH_MPS("mod010",146,2655,6548,6532.08,7);
    195 #if HOWMANY>2
    196     PUSH_MPS("mod011",4480,10958,-54558535,-62121982.55,7);
    197     PUSH_MPS("modglob",291,422,20740508,20430947.,7);
    198 #endif
    199 #if HOWMANY>3
    200     PUSH_MPS("noswot",182,128,-43,-43.0,7);
    201 #endif
    202 #if HOWMANY>1
    203     PUSH_MPS("nw04",36,87482,16862,16310.66667,7);
    204 #endif
    205     PUSH_MPS("p0033",16,33,3089,2520.57,7);
    206     PUSH_MPS("p0201",133,201,7615,6875.0,7);
    207     PUSH_MPS("p0282",241,282,258411,176867.50,7);
    208     PUSH_MPS("p0548",176,548,8691,315.29,7);
    209     PUSH_MPS("p2756",755,2756,3124,2688.75,7);
    210 #if HOWMANY>2
    211     PUSH_MPS("pk1",45,86,11.0,0.0,7);
    212 #endif
    213 #if HOWMANY>1
    214     PUSH_MPS("pp08a",136,240,7350.0,2748.3452381,7);
    215     PUSH_MPS("pp08aCUTS",246,240,7350.0,5480.6061563,7);
    216 #endif
    217 #if HOWMANY>2
    218     PUSH_MPS("qiu",1192,840,-132.873137,-931.638857,7);
    219 #endif
    220     PUSH_MPS("qnet1",503,1541,16029.692681,14274.102667,7);
    221     PUSH_MPS("qnet1_o",456,1541,16029.692681,12095.571667,7);
    222     PUSH_MPS("rentacar",6803,9557,30356761,28806137.644,7);
    223     PUSH_MPS("rgn",24,180,82.1999,48.7999,7);
    224 #if HOWMANY>3
    225     PUSH_MPS("rout",291,556,1077.56,981.86428571,7);
    226     PUSH_MPS("set1ch",492,712,54537.75,32007.73,7);
    227 #endif
    228     //    PUSH_MPS("seymour",4944,1372,423,403.84647413,7);
    229     PUSH_MPS("stein27",118,27,18,13.0,7);
    230 #if HOWMANY>1
    231     PUSH_MPS("stein45",331,45,30,22.0,7);
    232 #endif
    233     PUSH_MPS("vpm1",234,378,20,15.4167,7);
    234     PUSH_MPS("vpm2",234,378,13.75,9.8892645972,7);
     151    PUSH_MPS("arki001",1048,1388,7580813.0459,7579599.80787,7);
     152    PUSH_MPS("bell3a",123,133,878430.32,862578.64,0);
     153    PUSH_MPS("bell5",91,104,8966406.49,8608417.95,1);
     154    PUSH_MPS("blend2",274,353,7.598985,6.9156751140,0);
     155    PUSH_MPS("cap6000",2176,6000,-2451377,-2451537.325,1);
     156    PUSH_MPS("dano3mip",3202,13873,728.1111,576.23162474,7);
     157    PUSH_MPS("danoint",664,521,65.67,62.637280418,6);
     158    PUSH_MPS("dcmulti",290,548,188182,183975.5397,0);
     159    PUSH_MPS("dsbmip",1182,1886,-305.19817501,-305.19817501,0);
     160    PUSH_MPS("egout",98,141,568.101,149.589,0);
     161    PUSH_MPS("enigma",21,100,0.0,0.0,0);
     162    PUSH_MPS("fast0507",507,63009,174,172.14556668,5);
     163    PUSH_MPS("fiber",363,1298,405935.18000,156082.51759,0);
     164    PUSH_MPS("fixnet6",478,878,3983,1200.88,1);
     165    PUSH_MPS("flugpl",18,18,1201500,1167185.7,0);
     166    PUSH_MPS("gen",780,870,112313,112130.0,0);
     167    PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,1);
     168    PUSH_MPS("gesa2_o",1248,1224,25779856.372,25476489.678,1);
     169    PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,0);
     170    PUSH_MPS("gesa3_o",1224,1152,27991042.648,27833632.451,0);
     171    PUSH_MPS("gt2",29,188,21166.000,13460.233074,0);
     172    PUSH_MPS("harp2",112,2993,-73899798.00,-74353341.502,6);
     173    PUSH_MPS("khb05250",101,1350,106940226,95919464.0,0);
     174    PUSH_MPS("l152lav",97,1989,4722,4656.36,1);
     175    PUSH_MPS("lseu",28,89,1120,834.68,0);
     176    PUSH_MPS("mas74",13,151,11801.18573,10482.79528,3);
     177    PUSH_MPS("mas76",12,151,40005.05414,38893.9036,2);
     178    PUSH_MPS("misc03",96,160,3360,1910.,0);
     179    PUSH_MPS("misc06",820,1808,12850.8607,12841.6,0);
     180    PUSH_MPS("misc07",212,260,2810,1415.0,1);
     181    PUSH_MPS("mitre",2054,10724,115155,114740.5184,1);
     182    PUSH_MPS("mkc",3411,5325,-553.75,-611.85,7); // this is suboptimal
     183    PUSH_MPS("mod008",6,319,307,290.9,0);
     184    PUSH_MPS("mod010",146,2655,6548,6532.08,0);
     185    PUSH_MPS("mod011",4480,10958,-54558535,-62121982.55,2);
     186    PUSH_MPS("modglob",291,422,20740508,20430947.,2);
     187    PUSH_MPS("noswot",182,128,-43,-43.0,6);
     188    PUSH_MPS("nw04",36,87482,16862,16310.66667,1);
     189    PUSH_MPS("p0033",16,33,3089,2520.57,0);
     190    PUSH_MPS("p0201",133,201,7615,6875.0,0);
     191    PUSH_MPS("p0282",241,282,258411,176867.50,0);
     192    PUSH_MPS("p0548",176,548,8691,315.29,0);
     193    PUSH_MPS("p2756",755,2756,3124,2688.75,0);
     194    PUSH_MPS("pk1",45,86,11.0,0.0,2);
     195    PUSH_MPS("pp08a",136,240,7350.0,2748.3452381,1);
     196    PUSH_MPS("pp08aCUTS",246,240,7350.0,5480.6061563,1);
     197    PUSH_MPS("qiu",1192,840,-132.873137,-931.638857,3);
     198    PUSH_MPS("qnet1",503,1541,16029.692681,14274.102667,0);
     199    PUSH_MPS("qnet1_o",456,1541,16029.692681,12095.571667,0);
     200    PUSH_MPS("rentacar",6803,9557,30356761,28806137.644,0);
     201    PUSH_MPS("rgn",24,180,82.1999,48.7999,0);
     202    PUSH_MPS("rout",291,556,1077.56,981.86428571,3);
     203    PUSH_MPS("set1ch",492,712,54537.75,32007.73,5);
     204    PUSH_MPS("seymour",4944,1372,423,403.84647413,7);
     205    PUSH_MPS("seymour_1",4944,1372,410.76370,403.84647413,5);
     206    PUSH_MPS("stein27",118,27,18,13.0,0);
     207    PUSH_MPS("stein45",331,45,30,22.0,1);
     208    PUSH_MPS("swath",884,6805,497.603,334.4968581,7);
     209    PUSH_MPS("vpm1",234,378,20,15.4167,0);
     210    PUSH_MPS("vpm2",234,378,13.75,9.8892645972,0);
    235211#endif
    236212  }
     
    244220#endif
    245221  int numberFailures=0;
     222  int numberAttempts=0;
    246223 
    247224  /*
     
    249226  */
    250227  for (m = 0 ; m < mpsName.size() ; m++) {
    251     std::cout << "  processing mps file: " << mpsName[m]
    252               << " (" << m+1 << " out of " << mpsName.size() << ")\n";
    253     /*
    254       Stage 1: Read the MPS
    255       and make sure the size of the constraint matrix is correct.
    256     */
    257     CbcModel * model = new CbcModel(saveModel);
     228    int test = testSet[m];
     229    if(testSwitch>=test&&loSwitch<=test) {
     230      numberAttempts++;
     231      std::cout << "  processing mps file: " << mpsName[m]
     232                << " (" << m+1 << " out of " << mpsName.size() << ")\n";
     233      /*
     234        Stage 1: Read the MPS
     235        and make sure the size of the constraint matrix is correct.
     236      */
     237      CbcModel * model = new CbcModel(saveModel);
    258238     
    259     std::string fn = dirMiplib+mpsName[m] ;
    260     if (!CbcTestMpsFile(fn)) {
    261       std::cout << "ERROR: Cannot find MPS file " << fn << "\n";
    262       continue;
    263     }
    264     CoinDrand48(true,1234567);
    265     //printf("RAND1 %g %g\n",CoinDrand48(true,1234567),model->randomNumberGenerator()->randomDouble());
    266     //printf("RAND1 %g\n",CoinDrand48(true,1234567));
    267     model->solver()->readMps(fn.c_str(),"") ;
    268     assert(model->getNumRows() == nRows[m]) ;
    269     assert(model->getNumCols() == nCols[m]) ;
    270 
    271     /*
    272       Stage 2: Call solver to solve the problem.
    273       then check the return code and objective.
    274     */
    275 
     239      std::string fn = dirMiplib+mpsName[m] ;
     240      if (!CbcTestMpsFile(fn)) {
     241        std::cout << "ERROR: Cannot find MPS file " << fn << "\n";
     242        continue;
     243      }
     244      CoinDrand48(true,1234567);
     245      //printf("RAND1 %g %g\n",CoinDrand48(true,1234567),model->randomNumberGenerator()->randomDouble());
     246      //printf("RAND1 %g\n",CoinDrand48(true,1234567));
     247      model->solver()->readMps(fn.c_str(),"") ;
     248      assert(model->getNumRows() == nRows[m]) ;
     249      assert(model->getNumCols() == nCols[m]) ;
     250     
     251      /*
     252        Stage 2: Call solver to solve the problem.
     253        then check the return code and objective.
     254      */
     255     
    276256#ifdef CLP_FACTORIZATION_INSTRUMENT
    277     extern double factorization_instrument(int type);
    278     double facTime1=factorization_instrument(0);
    279     printf("Factorization - initial solve %g seconds\n",
    280            facTime1);
    281     timeTakenFac += facTime1;
    282 #endif
    283     double startTime = CoinCpuTime()+CoinCpuTimeJustChildren();
    284     if (model->getMaximumNodes()>200000) {
    285       model->setMaximumNodes(200000);
    286     }
    287     OsiClpSolverInterface * si =
    288       dynamic_cast<OsiClpSolverInterface *>(model->solver()) ;
    289     ClpSimplex * modelC = NULL;
    290     if (si) {
    291       // get clp itself
    292       modelC = si->getModelPtr();
    293       if (stuff&&stuff[9]) {
    294         // vector matrix!
     257      extern double factorization_instrument(int type);
     258      double facTime1=factorization_instrument(0);
     259      printf("Factorization - initial solve %g seconds\n",
     260             facTime1);
     261      timeTakenFac += facTime1;
     262#endif
     263      double startTime = CoinCpuTime()+CoinCpuTimeJustChildren();
     264      int testMaximumNodes=200000;
     265      if (testSwitch>1)
     266        testMaximumNodes=20000000;
     267      if (model->getMaximumNodes()>testMaximumNodes) {
     268        model->setMaximumNodes(testMaximumNodes);
     269      }
     270      OsiClpSolverInterface * si =
     271        dynamic_cast<OsiClpSolverInterface *>(model->solver()) ;
     272      ClpSimplex * modelC = NULL;
     273      if (si) {
     274        // get clp itself
     275        modelC = si->getModelPtr();
    295276        ClpMatrixBase * matrix = modelC->clpMatrix();
    296         if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
    297           ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
     277        ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
     278        if (stuff&&stuff[9]&&clpMatrix) {
     279          // vector matrix!
    298280          clpMatrix->makeSpecialColumnCopy();
    299281        }
    300       }
    301       model->checkModel();
    302       modelC->tightenPrimalBounds(0.0,0,true);
    303       model->initialSolve();
    304       if (modelC->dualBound()==1.0e10) {
    305         // user did not set - so modify
    306         // get largest scaled away from bound
    307         ClpSimplex temp=*modelC;
    308         temp.dual(0,7);
    309         double largestScaled=1.0e-12;
    310         double largest=1.0e-12;
    311         int numberRows = temp.numberRows();
    312         const double * rowPrimal = temp.primalRowSolution();
    313         const double * rowLower = temp.rowLower();
    314         const double * rowUpper = temp.rowUpper();
    315         const double * rowScale = temp.rowScale();
    316         int iRow;
    317         for (iRow=0;iRow<numberRows;iRow++) {
    318           double value = rowPrimal[iRow];
    319           double above = value-rowLower[iRow];
    320           double below = rowUpper[iRow]-value;
    321           if (above<1.0e12) {
    322             largest = CoinMax(largest,above);
    323           }
    324           if (below<1.0e12) {
    325             largest = CoinMax(largest,below);
    326           }
    327           if (rowScale) {
    328             double multiplier = rowScale[iRow];
    329             above *= multiplier;
    330             below *= multiplier;
    331           }
    332           if (above<1.0e12) {
    333             largestScaled = CoinMax(largestScaled,above);
    334           }
    335           if (below<1.0e12) {
    336             largestScaled = CoinMax(largestScaled,below);
    337           }
    338         }
    339        
    340         int numberColumns = temp.numberColumns();
    341         const double * columnPrimal = temp.primalColumnSolution();
    342         const double * columnLower = temp.columnLower();
    343         const double * columnUpper = temp.columnUpper();
    344         const double * columnScale = temp.columnScale();
    345         int iColumn;
    346         for (iColumn=0;iColumn<numberColumns;iColumn++) {
    347           double value = columnPrimal[iColumn];
    348           double above = value-columnLower[iColumn];
    349           double below = columnUpper[iColumn]-value;
    350           if (above<1.0e12) {
    351             largest = CoinMax(largest,above);
    352           }
    353           if (below<1.0e12) {
    354             largest = CoinMax(largest,below);
    355           }
    356           if (columnScale) {
    357             double multiplier = 1.0/columnScale[iColumn];
    358             above *= multiplier;
    359             below *= multiplier;
    360           }
    361           if (above<1.0e12) {
    362             largestScaled = CoinMax(largestScaled,above);
    363           }
    364           if (below<1.0e12) {
    365             largestScaled = CoinMax(largestScaled,below);
    366           }
    367         }
    368         std::cout<<"Largest (scaled) away from bound "<<largestScaled
    369                  <<" unscaled "<<largest<<std::endl;
    370282#if 0
    371         modelC->setDualBound(CoinMax(1.0001e8,
    372                                      CoinMin(1000.0*largestScaled,1.00001e10)));
     283        if (clpMatrix) {
     284          int numberRows = clpMatrix->getNumRows();
     285          int numberColumns = clpMatrix->getNumCols();
     286          double * elements = clpMatrix->getMutableElements();
     287          const int * row = clpMatrix->getIndices();
     288          const CoinBigIndex * columnStart = clpMatrix->getVectorStarts();
     289          const int * columnLength = clpMatrix->getVectorLengths();
     290          double * smallest = new double [numberRows];
     291          double * largest = new double [numberRows];
     292          char * flag = new char [numberRows];
     293          CoinZeroN(flag,numberRows);
     294          for (int i=0;i<numberRows;i++) {
     295            smallest[i]=COIN_DBL_MAX;
     296            largest[i]=0.0;
     297          }
     298          for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     299            bool isInteger = modelC->isInteger(iColumn);
     300            CoinBigIndex j;
     301            for (j=columnStart[iColumn];
     302                 j<columnStart[iColumn]+columnLength[iColumn];j++) {
     303              int iRow=row[j];
     304              double value = fabs(elements[j]);
     305              if (!isInteger)
     306                flag[iRow]=1;
     307              smallest[iRow]=CoinMin(smallest[iRow],value);
     308              largest[iRow]=CoinMax(largest[iRow],value);
     309            }
     310          }
     311          double * rowLower = modelC->rowLower();
     312          double * rowUpper = modelC->rowUpper();
     313          bool changed=false;
     314          for (int i=0;i<numberRows;i++) {
     315            if (flag[i]&&smallest[i]>10.0&&false) {
     316              smallest[i]=1.0/smallest[i];
     317              if (rowLower[i]>-1.0e20)
     318                rowLower[i] *= smallest[i];
     319              if (rowUpper[i]<1.0e20)
     320                rowUpper[i] *= smallest[i];
     321              changed = true;
     322            } else {
     323              smallest[i]=0.0;
     324            }
     325          }
     326          if (changed) {
     327            printf("SCALED\n");
     328            for (int iColumn=0;iColumn<numberColumns;iColumn++) {
     329              CoinBigIndex j;
     330              for (j=columnStart[iColumn];
     331                   j<columnStart[iColumn]+columnLength[iColumn];j++) {
     332                int iRow=row[j];
     333                if (smallest[iRow])
     334                  elements[j] *= smallest[iRow];
     335              }
     336            }
     337          }
     338          delete [] smallest;
     339          delete [] largest;
     340          delete [] flag;
     341        }
     342#endif
     343        model->checkModel();
     344        modelC->tightenPrimalBounds(0.0,0,true);
     345        model->initialSolve();
     346        if (modelC->dualBound()==1.0e10) {
     347          // user did not set - so modify
     348          // get largest scaled away from bound
     349          ClpSimplex temp=*modelC;
     350          temp.dual(0,7);
     351          double largestScaled=1.0e-12;
     352          double largest=1.0e-12;
     353          int numberRows = temp.numberRows();
     354          const double * rowPrimal = temp.primalRowSolution();
     355          const double * rowLower = temp.rowLower();
     356          const double * rowUpper = temp.rowUpper();
     357          const double * rowScale = temp.rowScale();
     358          int iRow;
     359          for (iRow=0;iRow<numberRows;iRow++) {
     360            double value = rowPrimal[iRow];
     361            double above = value-rowLower[iRow];
     362            double below = rowUpper[iRow]-value;
     363            if (above<1.0e12) {
     364              largest = CoinMax(largest,above);
     365            }
     366            if (below<1.0e12) {
     367              largest = CoinMax(largest,below);
     368            }
     369            if (rowScale) {
     370              double multiplier = rowScale[iRow];
     371              above *= multiplier;
     372              below *= multiplier;
     373            }
     374            if (above<1.0e12) {
     375              largestScaled = CoinMax(largestScaled,above);
     376            }
     377            if (below<1.0e12) {
     378              largestScaled = CoinMax(largestScaled,below);
     379            }
     380          }
     381         
     382          int numberColumns = temp.numberColumns();
     383          const double * columnPrimal = temp.primalColumnSolution();
     384          const double * columnLower = temp.columnLower();
     385          const double * columnUpper = temp.columnUpper();
     386          const double * columnScale = temp.columnScale();
     387          int iColumn;
     388          for (iColumn=0;iColumn<numberColumns;iColumn++) {
     389            double value = columnPrimal[iColumn];
     390            double above = value-columnLower[iColumn];
     391            double below = columnUpper[iColumn]-value;
     392            if (above<1.0e12) {
     393              largest = CoinMax(largest,above);
     394            }
     395            if (below<1.0e12) {
     396              largest = CoinMax(largest,below);
     397            }
     398            if (columnScale) {
     399              double multiplier = 1.0/columnScale[iColumn];
     400              above *= multiplier;
     401              below *= multiplier;
     402            }
     403            if (above<1.0e12) {
     404              largestScaled = CoinMax(largestScaled,above);
     405            }
     406            if (below<1.0e12) {
     407              largestScaled = CoinMax(largestScaled,below);
     408            }
     409          }
     410          std::cout<<"Largest (scaled) away from bound "<<largestScaled
     411                   <<" unscaled "<<largest<<std::endl;
     412#if 0
     413          modelC->setDualBound(CoinMax(1.0001e8,
     414                                       CoinMin(1000.0*largestScaled,1.00001e10)));
    373415#else
    374         modelC->setDualBound(CoinMax(1.0001e9,
    375                                      CoinMin(1000.0*largestScaled,1.0001e10)));
    376 #endif
    377       }
    378     }
    379     model->setMinimumDrop(CoinMin(5.0e-2,
    380                               fabs(model->getMinimizationObjValue())*1.0e-3+1.0e-4));
    381     if (CoinAbs(model->getMaximumCutPassesAtRoot())<=100) {
    382       if (model->getNumCols()<500) {
    383         model->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
    384       } else if (model->getNumCols()<5000) {
    385         model->setMaximumCutPassesAtRoot(100); // use minimum drop
    386       } else {
    387         model->setMaximumCutPassesAtRoot(20);
    388       }
    389     }
    390     // If defaults then increase trust for small models
    391     if (model->numberStrong()==5&&model->numberBeforeTrust()==10) {
    392       int numberColumns = model->getNumCols();
    393       if (numberColumns<=50) {
    394         model->setNumberBeforeTrust(1000);
    395       } else if (numberColumns<=100) {
    396         model->setNumberBeforeTrust(100);
    397       } else if (numberColumns<=300) {
    398         model->setNumberBeforeTrust(50);
    399       }
    400     }
    401     //if (model->getNumCols()>=500) {
     416          modelC->setDualBound(CoinMax(1.0001e9,
     417                                       CoinMin(1000.0*largestScaled,1.0001e10)));
     418#endif
     419        }
     420      }
     421      model->setMinimumDrop(CoinMin(5.0e-2,
     422                                    fabs(model->getMinimizationObjValue())*1.0e-3+1.0e-4));
     423      if (CoinAbs(model->getMaximumCutPassesAtRoot())<=100) {
     424        if (model->getNumCols()<500) {
     425          model->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
     426        } else if (model->getNumCols()<5000) {
     427          model->setMaximumCutPassesAtRoot(100); // use minimum drop
     428        } else {
     429          model->setMaximumCutPassesAtRoot(20);
     430        }
     431      }
     432      // If defaults then increase trust for small models
     433      if (model->numberStrong()==5&&model->numberBeforeTrust()==10) {
     434        int numberColumns = model->getNumCols();
     435        if (numberColumns<=50) {
     436          model->setNumberBeforeTrust(1000);
     437        } else if (numberColumns<=100) {
     438          model->setNumberBeforeTrust(100);
     439        } else if (numberColumns<=300) {
     440          model->setNumberBeforeTrust(50);
     441        }
     442      }
     443      //if (model->getNumCols()>=500) {
    402444      // switch off Clp stuff
    403445      //model->setFastNodeDepth(-1);
    404     //}
    405     if (model->getNumCols()==-2756) {
    406       // p2756
    407       std::string problemName ;
    408       model->solver()->getStrParam(OsiProbName,problemName) ;
    409       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     446      //}
     447      if (model->getNumCols()==-2756) {
     448        // p2756
     449        std::string problemName ;
     450        model->solver()->getStrParam(OsiProbName,problemName) ;
     451        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     452      }
     453      if (model->getNumCols()==-201) {
     454        // p201
     455        std::string problemName ;
     456        model->solver()->getStrParam(OsiProbName,problemName) ;
     457        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     458      }
     459      if (model->getNumCols()==-104) {
     460        // bell5
     461        std::string problemName ;
     462        model->solver()->getStrParam(OsiProbName,problemName) ;
     463        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     464      }
     465      if (model->getNumCols()==-548&&model->getNumRows()==176) {
     466        // p0548
     467        std::string problemName ;
     468        model->solver()->getStrParam(OsiProbName,problemName) ;
     469        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     470      }
     471      if (model->getNumCols()==-160) {
     472        // misc03
     473        std::string problemName ;
     474        model->solver()->getStrParam(OsiProbName,problemName) ;
     475        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     476      }
     477      if (model->getNumCols()==-353) {
     478        // blend2
     479        std::string problemName ;
     480        model->solver()->getStrParam(OsiProbName,problemName) ;
     481        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     482      }
     483      if (model->getNumCols()==-100&&model->getNumRows()==21) {
     484        // enigma
     485        std::string problemName ;
     486        model->solver()->getStrParam(OsiProbName,problemName) ;
     487        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     488      }
     489      if (model->getNumCols()==-1541) {
     490        // qnet1
     491        std::string problemName ;
     492        model->solver()->getStrParam(OsiProbName,problemName) ;
     493        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     494      }
     495      if (model->getNumCols()==-10724) {
     496        // mitre
     497        std::string problemName ;
     498        model->solver()->getStrParam(OsiProbName,problemName) ;
     499        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     500      }
     501      if (model->getNumCols()==-1224) {
     502        //PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,7);
     503        // gesa2
     504        std::string problemName ;
     505        model->solver()->getStrParam(OsiProbName,problemName) ;
     506        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     507      }
     508      if (model->getNumCols()==-1152&&model->getNumRows()==1368) {
     509        //PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,7);
     510        // gesa3
     511        std::string problemName ;
     512        model->solver()->getStrParam(OsiProbName,problemName) ;
     513        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     514      }
     515      if (model->getNumCols()==-1152&&model->getNumRows()==1224) {
     516        //PUSH_MPS("gesa3_o",1224,1152,27991042.648,27833632.451,7);
     517        // gesa3
     518        std::string problemName ;
     519        model->solver()->getStrParam(OsiProbName,problemName) ;
     520        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     521      }
     522      if (model->getNumCols()==-282) {
     523        //PUSH_MPS("p0282",241,282,258411,176867.50,7);
     524        // p0282
     525        std::string problemName ;
     526        model->solver()->getStrParam(OsiProbName,problemName) ;
     527        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     528      }
     529      if (model->getNumCols()==-141) {
     530        // egout
     531        std::string problemName ;
     532        model->solver()->getStrParam(OsiProbName,problemName) ;
     533        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     534      }
     535      if (model->getNumCols()==-378) {
     536        // vpm2
     537        std::string problemName ;
     538        model->solver()->getStrParam(OsiProbName,problemName) ;
     539        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     540      }
     541      if (model->getNumCols()==-240&&model->getNumRows()==246) {
     542        // pp08aCUTS
     543        std::string problemName ;
     544        model->solver()->getStrParam(OsiProbName,problemName) ;
     545        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     546      }
     547      if (model->getNumCols()==-240&&model->getNumRows()==136) {
     548        // pp08a
     549        std::string problemName ;
     550        model->solver()->getStrParam(OsiProbName,problemName) ;
     551        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     552      }
     553      if (model->getNumCols()==-1372&&model->getNumRows()==4944) {
     554        // seymour1
     555        std::string problemName ;
     556        model->solver()->getStrParam(OsiProbName,problemName) ;
     557        model->solver()->activateRowCutDebugger(problemName.c_str()) ;
     558      }
     559      setCutAndHeuristicOptions(*model);
     560      if (si) {
     561#ifdef CLP_MULTIPLE_FACTORIZATIONS
     562        if (!modelC->factorization()->isDenseOrSmall()) {
     563          int denseCode = stuff ? static_cast<int> (stuff[4]) : -1;
     564          int smallCode = stuff ? static_cast<int> (stuff[10]) : -1;
     565          if (stuff&&stuff[8]>=1) {
     566            if (denseCode<0)
     567              denseCode=40;
     568            if (smallCode<0)
     569              smallCode=40;
     570          }
     571          if (denseCode>0)
     572            modelC->factorization()->setGoDenseThreshold(denseCode);
     573          if (smallCode>0)
     574            modelC->factorization()->setGoSmallThreshold(smallCode);
     575          if (denseCode>=modelC->numberRows()) {
     576            //printf("problem going dense\n");
     577            //modelC->factorization()->goDenseOrSmall(modelC->numberRows());
     578          }
     579        }
     580#endif
     581        if (stuff&&stuff[8]>=1) {
     582          if (modelC->numberColumns()+modelC->numberRows()<=10000&&
     583              model->fastNodeDepth()==-1)
     584            model->setFastNodeDepth(-9);
     585        }
     586      }
     587      //OsiObject * obj = new CbcBranchToFixLots(model,0.3,0.0,3,3000003);
     588      //model->addObjects(1,&obj);
     589      //delete obj;
     590      model->branchAndBound();
     591#ifdef CLP_FACTORIZATION_INSTRUMENT
     592      double facTime=factorization_instrument(0);
     593      printf("Factorization %g seconds\n",
     594             facTime);
     595      timeTakenFac += facTime;
     596#endif
     597     
     598      double timeOfSolution = CoinCpuTime()+CoinCpuTimeJustChildren()-startTime;
     599      // Print more statistics
     600      std::cout<<"Cuts at root node changed objective from "<<model->getContinuousObjective()
     601               <<" to "<<model->rootObjectiveAfterCuts()<<std::endl;
     602      int numberGenerators = model->numberCutGenerators();
     603      for (int iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
     604        CbcCutGenerator * generator = model->cutGenerator(iGenerator);
     605#ifndef CLP_INVESTIGATE
     606        CglImplication * implication = dynamic_cast<CglImplication*>(generator->generator());
     607        if (implication)
     608          continue;
     609#endif
     610        std::cout<<generator->cutGeneratorName()<<" was tried "
     611                 <<generator->numberTimesEntered()<<" times and created "
     612                 <<generator->numberCutsInTotal()<<" cuts of which "
     613                 <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
     614        if (generator->timing())
     615          std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
     616        else
     617          std::cout<<std::endl;
     618      }
     619      if (!model->status()) {
     620        double soln = model->getObjValue(); 
     621        double tolerance = CoinMax(1.0e-5,1.0e-5*CoinMin(fabs(soln),fabs(objValue[m])));
     622        //CoinRelFltEq eq(1.0e-3) ;
     623        if (fabs(soln-objValue[m])<tolerance) {
     624          std::cout
     625            <<"cbc_clp"<<" "
     626            << soln << " = " << objValue[m] << " ; okay";
     627          numProbSolved++;
     628        } else  {
     629          std::cout <<"cbc_clp" <<" " <<soln << " != " <<objValue[m]
     630                    << "; error=" << fabs(objValue[m] - soln);
     631          numberFailures++;
     632          //#ifdef COIN_DEVELOP
     633          //abort();
     634          //#endif
     635        }
     636      } else {
     637        std::cout << "cbc_clp error; too many nodes" ;
     638      }
     639      timeTaken += timeOfSolution;
     640      std::cout<<" - took " <<timeOfSolution<<" seconds.("<<
     641        model->getNodeCount()<<" / "<<model->getIterationCount()<<
     642        " ) subtotal "<<timeTaken
     643               <<" ("<<mpsName[m]<<")"<<std::endl;
     644      delete model;
    410645    }
    411     if (model->getNumCols()==-201) {
    412       // p201
    413       std::string problemName ;
    414       model->solver()->getStrParam(OsiProbName,problemName) ;
    415       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    416     }
    417     if (model->getNumCols()==-104) {
    418       // bell5
    419       std::string problemName ;
    420       model->solver()->getStrParam(OsiProbName,problemName) ;
    421       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    422     }
    423     if (model->getNumCols()==-548&&model->getNumRows()==176) {
    424       // p0548
    425       std::string problemName ;
    426       model->solver()->getStrParam(OsiProbName,problemName) ;
    427       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    428     }
    429     if (model->getNumCols()==-160) {
    430       // misc03
    431       std::string problemName ;
    432       model->solver()->getStrParam(OsiProbName,problemName) ;
    433       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    434     }
    435     if (model->getNumCols()==-353) {
    436       // blend2
    437       std::string problemName ;
    438       model->solver()->getStrParam(OsiProbName,problemName) ;
    439       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    440     }
    441     if (model->getNumCols()==-100&&model->getNumRows()==21) {
    442       // enigma
    443       std::string problemName ;
    444       model->solver()->getStrParam(OsiProbName,problemName) ;
    445       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    446     }
    447     if (model->getNumCols()==-1541) {
    448       // qnet1
    449       std::string problemName ;
    450       model->solver()->getStrParam(OsiProbName,problemName) ;
    451       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    452     }
    453     if (model->getNumCols()==-10724) {
    454       // mitre
    455       std::string problemName ;
    456       model->solver()->getStrParam(OsiProbName,problemName) ;
    457       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    458     }
    459     if (model->getNumCols()==-1224) {
    460       //PUSH_MPS("gesa2",1392,1224,25779856.372,25476489.678,7);
    461       // gesa2
    462       std::string problemName ;
    463       model->solver()->getStrParam(OsiProbName,problemName) ;
    464       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    465     }
    466     if (model->getNumCols()==-1152) {
    467       //PUSH_MPS("gesa3",1368,1152,27991042.648,27833632.451,7);
    468       // gesa3
    469       std::string problemName ;
    470       model->solver()->getStrParam(OsiProbName,problemName) ;
    471       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    472     }
    473     if (model->getNumCols()==-282) {
    474       //PUSH_MPS("p0282",241,282,258411,176867.50,7);
    475       // p0282
    476       std::string problemName ;
    477       model->solver()->getStrParam(OsiProbName,problemName) ;
    478       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    479     }
    480     if (model->getNumCols()==-141) {
    481       // egout
    482       std::string problemName ;
    483       model->solver()->getStrParam(OsiProbName,problemName) ;
    484       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    485     }
    486     if (model->getNumCols()==-378) {
    487       // vpm2
    488       std::string problemName ;
    489       model->solver()->getStrParam(OsiProbName,problemName) ;
    490       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    491     }
    492     if (model->getNumCols()==-240&&model->getNumRows()==246) {
    493       // pp08aCUTS
    494       std::string problemName ;
    495       model->solver()->getStrParam(OsiProbName,problemName) ;
    496       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    497     }
    498     if (model->getNumCols()==-240&&model->getNumRows()==136) {
    499       // pp08a
    500       std::string problemName ;
    501       model->solver()->getStrParam(OsiProbName,problemName) ;
    502       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    503     }
    504     if (model->getNumCols()==-1372&&model->getNumRows()==4944) {
    505       // seymour1
    506       std::string problemName ;
    507       model->solver()->getStrParam(OsiProbName,problemName) ;
    508       model->solver()->activateRowCutDebugger(problemName.c_str()) ;
    509     }
    510     setCutAndHeuristicOptions(*model);
    511     if (si) {
    512 #ifdef CLP_MULTIPLE_FACTORIZATIONS
    513       if (!modelC->factorization()->isDenseOrSmall()) {
    514         int denseCode = stuff ? static_cast<int> (stuff[4]) : -1;
    515         int smallCode = stuff ? static_cast<int> (stuff[10]) : -1;
    516         if (stuff&&stuff[8]>=1) {
    517           if (denseCode<0)
    518             denseCode=40;
    519           if (smallCode<0)
    520             smallCode=40;
    521         }
    522         if (denseCode>0)
    523           modelC->factorization()->setGoDenseThreshold(denseCode);
    524         if (smallCode>0)
    525           modelC->factorization()->setGoSmallThreshold(smallCode);
    526         if (denseCode>=modelC->numberRows()) {
    527           //printf("problem going dense\n");
    528           //modelC->factorization()->goDenseOrSmall(modelC->numberRows());
    529         }
    530       }
    531 #endif
    532       if (stuff&&stuff[8]>=1) {
    533         if (modelC->numberColumns()+modelC->numberRows()<=100000&&
    534             model->fastNodeDepth()==-1)
    535           model->setFastNodeDepth(-9);
    536       }
    537     }
    538     //OsiObject * obj = new CbcBranchToFixLots(model,0.3,0.0,3,3000003);
    539     //model->addObjects(1,&obj);
    540     //delete obj;
    541     model->branchAndBound();
    542 #ifdef CLP_FACTORIZATION_INSTRUMENT
    543     double facTime=factorization_instrument(0);
    544     printf("Factorization %g seconds\n",
    545            facTime);
    546     timeTakenFac += facTime;
    547 #endif
    548      
    549     double timeOfSolution = CoinCpuTime()+CoinCpuTimeJustChildren()-startTime;
    550     // Print more statistics
    551     std::cout<<"Cuts at root node changed objective from "<<model->getContinuousObjective()
    552              <<" to "<<model->rootObjectiveAfterCuts()<<std::endl;
    553     int numberGenerators = model->numberCutGenerators();
    554     for (int iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
    555       CbcCutGenerator * generator = model->cutGenerator(iGenerator);
    556 #ifndef CLP_INVESTIGATE
    557       CglImplication * implication = dynamic_cast<CglImplication*>(generator->generator());
    558       if (implication)
    559         continue;
    560 #endif
    561       std::cout<<generator->cutGeneratorName()<<" was tried "
    562                <<generator->numberTimesEntered()<<" times and created "
    563                <<generator->numberCutsInTotal()<<" cuts of which "
    564                <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
    565       if (generator->timing())
    566         std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
    567       else
    568         std::cout<<std::endl;
    569     }
    570     if (!model->status()) {
    571       double soln = model->getObjValue(); 
    572       double tolerance = CoinMax(1.0e-5,1.0e-5*CoinMin(fabs(soln),fabs(objValue[m])));
    573       //CoinRelFltEq eq(1.0e-3) ;
    574       if (fabs(soln-objValue[m])<tolerance) {
    575         std::cout
    576           <<"cbc_clp"<<" "
    577           << soln << " = " << objValue[m] << " ; okay";
    578         numProbSolved++;
    579       } else  {
    580         std::cout <<"cbc_clp" <<" " <<soln << " != " <<objValue[m]
    581                   << "; error=" << fabs(objValue[m] - soln);
    582         numberFailures++;
    583         //#ifdef COIN_DEVELOP
    584         //abort();
    585         //#endif
    586       }
    587     } else {
    588       std::cout << "error; too many nodes" ;
    589     }
    590     timeTaken += timeOfSolution;
    591     std::cout<<" - took " <<timeOfSolution<<" seconds.("<<
    592       model->getNodeCount()<<" / "<<model->getIterationCount()<<
    593       " ) subtotal "<<timeTaken<<std::endl;
    594     delete model;
    595646  }
    596647  int returnCode=0;
     
    600651    <<numProbSolved
    601652    <<" out of "
    602     <<objValue.size();
    603   int numberOnNodes = objValue.size()-numProbSolved-numberFailures;
     653    <<numberAttempts;
     654  int numberOnNodes = numberAttempts-numProbSolved-numberFailures;
    604655  if (numberFailures||numberOnNodes) {
    605656    if (numberOnNodes) {
     
    616667    <<" seconds."
    617668    <<std::endl;
    618   if (unitTestOnly) {
     669  if (testSwitch==-2) {
    619670    if(numberFailures||numberOnNodes) {
    620671      printf("****** Unit Test failed\n");
Note: See TracChangeset for help on using the changeset viewer.