source: trunk/Cbc/src/unitTestClp.cpp @ 1622

Last change on this file since 1622 was 1622, checked in by lou, 8 years ago

Clean up CbcClpUnitTest?. Easy activation for row cut debugger. Documentation.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 23.5 KB
Line 
1/* $Id: unitTestClp.cpp 1622 2011-03-25 17:20:10Z lou $ */
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6
7#include <cstdio>
8#include <string>
9#include <iostream>
10#include <iomanip>
11
12#include "CoinTime.hpp"
13#include "CbcModel.hpp"
14#include "CbcHeuristic.hpp"
15#include "CbcCutGenerator.hpp"
16#include "CbcBranchCut.hpp"
17#include "CglProbing.hpp"
18#include "OsiClpSolverInterface.hpp"
19#include "ClpFactorization.hpp"
20#include "OsiRowCutDebugger.hpp"
21//#############################################################################
22
23#ifdef NDEBUG
24#undef NDEBUG
25#endif
26
27//#############################################################################
28
29// Display message on stdout and stderr
30void testingMessage( const char * const msg )
31{
32    std::cout << msg;
33    //cout <<endl <<"*****************************************"
34    //     <<endl <<msg <<endl;
35}
36
37//#############################################################################
38
39static inline bool CbcTestFile(const std::string name)
40{
41    FILE *fp = fopen(name.c_str(), "r");
42    if (fp) {
43        fclose(fp);
44        return true;
45    }
46    return false;
47}
48
49//#############################################################################
50
51bool CbcTestMpsFile(std::string& fname)
52{
53    if (CbcTestFile(fname)) {
54        return true;
55    }
56    if (CbcTestFile(fname + ".mps")) {
57        fname += ".mps";
58        return true;
59    }
60    if (CbcTestFile(fname + ".MPS")) {
61        fname += ".MPS";
62        return true;
63    }
64#ifdef COIN_HAS_ZLIB
65    if (CbcTestFile(fname + ".gz")) {
66        return true;
67    }
68    if (CbcTestFile(fname + ".mps.gz")) {
69        fname += ".mps";
70        return true;
71    }
72    if (CbcTestFile(fname + ".MPS.gz")) {
73        fname += ".MPS";
74        return true;
75    }
76    if (CbcTestFile(fname + ".MPS.GZ")) {
77        fname += ".MPS";
78        return true;
79    }
80#endif
81    return false;
82}
83
84//#############################################################################
85/*
86  jjf: testSwitch -2 unitTest, -1 normal (==2)
87
88  MiplibTest might be more appropriate.
89
90  TestSwitch and stuff[6] together control how much of miplib is executed:
91    For testSwitch set to:
92      -2: solve p0033 and p0201 only (the unit test)
93      -1: solve miplib sets #0 and #1
94       0: solve nothing
95       k: execute sets j:k, where j is determined by the value of stuff[6]
96  The last parameter of PUSH_MPS specifies the test set membership.
97
98  For -miplib, -extra2 sets testSwitch, -extra3 sets stuff[6]. The command
99    cbc -extra2 -2 -miplib
100  will execute the unit test on the miplib directory.
101
102  dirMiplib should end in the directory separator character for the platform.
103
104  If you want to activate the row cut debugger for a given problem, change the
105  last parameter of the PUSH_MPS macro for the problem to true.
106
107  Returns 0 if all goes well, -1 if the Miplib directory is missing, otherwise
108     100*(number with bad objective)+(number that exceeded node limit)
109*/
110int CbcClpUnitTest (const CbcModel &saveModel, const std::string &dirMiplib,
111                    int testSwitch, const double *stuff)
112{
113  // Stop Windows popup
114  WindowsErrorPopupBlocker() ;
115  unsigned int m ;
116
117  // Do an existence check.
118  std::string test1 = dirMiplib+"p0033";
119  bool doTest = CbcTestMpsFile(test1);
120  if (!doTest) {
121    std::cout
122      << "Not doing miplib run as can't find mps files." << std::endl
123      << "Perhaps you're trying to read gzipped (.gz) files without libz?"
124      << std::endl ;
125    return (-1) ;
126  }
127  int dfltPrecision = std::cout.precision() ;
128/*
129  Set the range of problems to be tested. testSwitch = -2 is special and is
130  picked up below.
131*/
132    int loSet = 0 ;
133    int hiSet = 0 ;
134    if (testSwitch == -1) {
135      loSet = 0 ;
136      hiSet = 1 ;
137    } else if (testSwitch >= 0) {
138      loSet = static_cast<int>(stuff[6]) ;
139      hiSet = testSwitch ;
140      std::cout
141        << "Solving miplib problems in sets " << loSet
142        << ":" << hiSet << "." << std::endl ;
143    }
144/*
145  Vectors to hold test problem names and characteristics.
146*/
147    std::vector<std::string> mpsName ;
148    std::vector<int> nRows ;
149    std::vector<int> nCols ;
150    std::vector<double> objValueC ;
151    std::vector<double> objValue ;
152    std::vector<int> testSet ;
153    std::vector<bool> rowCutDebugger ;
154/*
155  A macro to make the vector creation marginally readable. Parameters are
156  name, rows, columns, integer objective, continuous objective, set ID,
157  row cut debugger
158
159  To enable the row cut debugger for a given problem, change the last
160  parameter to true. Don't forget to turn it off before committing changes!
161*/
162#define PUSH_MPS(zz_mpsName_zz,\
163                 zz_nRows_zz,zz_nCols_zz,zz_objValue_zz,zz_objValueC_zz, \
164                 zz_testSet_zz, zz_rcDbg_zz) \
165  mpsName.push_back(zz_mpsName_zz) ; \
166  nRows.push_back(zz_nRows_zz) ; \
167  nCols.push_back(zz_nCols_zz) ; \
168  objValueC.push_back(zz_objValueC_zz) ; \
169  testSet.push_back(zz_testSet_zz) ; \
170  objValue.push_back(zz_objValue_zz) ; \
171  rowCutDebugger.push_back(zz_rcDbg_zz) ;
172/*
173  Push the miplib problems. Except for -2 (unitTest), push all, even if we're
174  not going to do all of them.
175*/
176  if (testSwitch == -2) {
177      PUSH_MPS("p0033", 16, 33, 3089, 2520.57, 0, false);
178      PUSH_MPS("p0201", 133, 201, 7615, 6875.0, 0, false);
179      // PUSH_MPS("flugpl", 18, 18, 1201500, 1167185.7, 0, false);
180  } else {
181/*
182  Load up the problem vector. Note that the row counts here include the
183  objective function.
184*/
185    PUSH_MPS("10teams", 230, 2025, 924, 917, 1, false);
186    PUSH_MPS("air03", 124, 10757, 340160, 338864.25, 0, false);
187    PUSH_MPS("air04", 823, 8904, 56137, 55535.436, 2, false);
188    PUSH_MPS("air05", 426, 7195, 26374, 25877.609, 2, false);
189    PUSH_MPS("arki001", 1048, 1388, 7580813.0459, 7579599.80787, 7, false);
190    PUSH_MPS("bell3a", 123, 133, 878430.32, 862578.64, 0, false);
191    PUSH_MPS("bell5", 91, 104, 8966406.49, 8608417.95, 1, false);
192    PUSH_MPS("blend2", 274, 353, 7.598985, 6.9156751140, 0, false);
193    PUSH_MPS("cap6000", 2176, 6000, -2451377, -2451537.325, 1, false);
194    PUSH_MPS("dano3mip", 3202, 13873, 728.1111, 576.23162474, 7, false);
195    PUSH_MPS("danoint", 664, 521, 65.67, 62.637280418, 6, false);
196    PUSH_MPS("dcmulti", 290, 548, 188182, 183975.5397, 0, false);
197    PUSH_MPS("dsbmip", 1182, 1886, -305.19817501, -305.19817501, 0, false);
198    PUSH_MPS("egout", 98, 141, 568.101, 149.589, 0, false);
199    PUSH_MPS("enigma", 21, 100, 0.0, 0.0, 0, false);
200    PUSH_MPS("fast0507", 507, 63009, 174, 172.14556668, 5, false);
201    PUSH_MPS("fiber", 363, 1298, 405935.18000, 156082.51759, 0, false);
202    PUSH_MPS("fixnet6", 478, 878, 3983, 1200.88, 1, false);
203    PUSH_MPS("flugpl", 18, 18, 1201500, 1167185.7, 0, false);
204    PUSH_MPS("gen", 780, 870, 112313, 112130.0, 0, false);
205    PUSH_MPS("gesa2", 1392, 1224, 25779856.372, 25476489.678, 1, false);
206    PUSH_MPS("gesa2_o", 1248, 1224, 25779856.372, 25476489.678, 1, false);
207    PUSH_MPS("gesa3", 1368, 1152, 27991042.648, 27833632.451, 0, false);
208    PUSH_MPS("gesa3_o", 1224, 1152, 27991042.648, 27833632.451, 0, false);
209    PUSH_MPS("gt2", 29, 188, 21166.000, 13460.233074, 0, false);
210    PUSH_MPS("harp2", 112, 2993, -73899798.00, -74353341.502, 6, false);
211    PUSH_MPS("khb05250", 101, 1350, 106940226, 95919464.0, 0, false);
212    PUSH_MPS("l152lav", 97, 1989, 4722, 4656.36, 1, false);
213    PUSH_MPS("lseu", 28, 89, 1120, 834.68, 0, false);
214    PUSH_MPS("mas74", 13, 151, 11801.18573, 10482.79528, 3, false);
215    PUSH_MPS("mas76", 12, 151, 40005.05414, 38893.9036, 2, false);
216    PUSH_MPS("misc03", 96, 160, 3360, 1910., 0, false);
217    PUSH_MPS("misc06", 820, 1808, 12850.8607, 12841.6, 0, false);
218    PUSH_MPS("misc07", 212, 260, 2810, 1415.0, 1, false);
219    PUSH_MPS("mitre", 2054, 10724, 115155, 114740.5184, 1, false);
220    PUSH_MPS("mkc", 3411, 5325, -553.75, -611.85, 7, false); // suboptimal
221    PUSH_MPS("mod008", 6, 319, 307, 290.9, 0, false);
222    PUSH_MPS("mod010", 146, 2655, 6548, 6532.08, 0, false);
223    PUSH_MPS("mod011", 4480, 10958, -54558535, -62121982.55, 2, false);
224    PUSH_MPS("modglob", 291, 422, 20740508, 20430947., 2, false);
225    PUSH_MPS("noswot", 182, 128, -43, -43.0, 6, false);
226    PUSH_MPS("nw04", 36, 87482, 16862, 16310.66667, 1, false);
227    PUSH_MPS("p0033", 16, 33, 3089, 2520.57, 0, false);
228    PUSH_MPS("p0201", 133, 201, 7615, 6875.0, 0, false);
229    PUSH_MPS("p0282", 241, 282, 258411, 176867.50, 0, false);
230    PUSH_MPS("p0548", 176, 548, 8691, 315.29, 0, false);
231    PUSH_MPS("p2756", 755, 2756, 3124, 2688.75, 0, false);
232    PUSH_MPS("pk1", 45, 86, 11.0, 0.0, 2, false);
233    PUSH_MPS("pp08a", 136, 240, 7350.0, 2748.3452381, 1, false);
234    PUSH_MPS("pp08aCUTS", 246, 240, 7350.0, 5480.6061563, 1, false);
235    PUSH_MPS("qiu", 1192, 840, -132.873137, -931.638857, 3, false);
236    PUSH_MPS("qnet1", 503, 1541, 16029.692681, 14274.102667, 0, false);
237    PUSH_MPS("qnet1_o", 456, 1541, 16029.692681, 12095.571667, 0, false);
238    PUSH_MPS("rentacar", 6803, 9557, 30356761, 28806137.644, 0, false);
239    PUSH_MPS("rgn", 24, 180, 82.1999, 48.7999, 0, false);
240    PUSH_MPS("rout", 291, 556, 1077.56, 981.86428571, 3, false);
241    PUSH_MPS("set1ch", 492, 712, 54537.75, 32007.73, 5, false);
242    PUSH_MPS("seymour", 4944, 1372, 423, 403.84647413, 7, false);
243    PUSH_MPS("seymour_1", 4944, 1372, 410.76370, 403.84647413, 5, false);
244    PUSH_MPS("stein27", 118, 27, 18, 13.0, 0, false);
245    PUSH_MPS("stein45", 331, 45, 30, 22.0, 1, false);
246    PUSH_MPS("swath", 884, 6805, 497.603, 334.4968581, 7, false);
247    PUSH_MPS("vpm1", 234, 378, 20, 15.4167, 0, false);
248    PUSH_MPS("vpm2", 234, 378, 13.75, 9.8892645972, 0, false);
249  }
250#undef PUSH_MPS
251
252/*
253  Normally the problems are executed in order. Define RANDOM_ORDER below to
254  randomize.
255
256  #define RANDOM_ORDER
257*/
258  int which[100];
259  int nLoop = static_cast<int>(mpsName.size());
260  assert (nLoop <= 100);
261  for (int i = 0; i < nLoop; i++) which[i] = i;
262
263# ifdef RANDOM_ORDER
264  unsigned int iTime = static_cast<unsigned int>(CoinGetTimeOfDay()-1.256e9);
265  std::cout << "Time (seed) " << iTime << "." << std::endl ;
266  double sort[100];
267  CoinDrand48(true,iTime);
268  for (int i = 0; i < nLoop; i++) sort[i] = CoinDrand48();
269  CoinSort_2(sort,sort+nLoop,which);
270# endif
271
272  int problemCnt = 0;
273  for (m = 0 ; m < mpsName.size() ; m++) {
274    int setID = testSet[m];
275    if (loSet <= setID && setID <= hiSet) problemCnt++;
276  }
277
278  int numberFailures = 0;
279  int numberAttempts = 0;
280  int numProbSolved = 0;
281  double timeTaken = 0.0;
282
283//#define CLP_FACTORIZATION_INSTRUMENT
284# ifdef CLP_FACTORIZATION_INSTRUMENT
285  double timeTakenFac = 0.0;
286# endif
287
288/*
289  Open the main loop to step through the MPS problems.
290*/
291  for (unsigned int mw = 0 ; mw < mpsName.size() ; mw++) {
292    m = which[mw];
293    int setID = testSet[m];
294    // Skip if problem is not in specified problem set(s)
295    if (!(loSet <= setID && setID <= hiSet)) continue ;
296   
297    numberAttempts++;
298    std::cout << "  processing mps file: " << mpsName[m]
299              << " (" << numberAttempts << " out of "
300              << problemCnt << ")" << std::endl ;
301/*
302  Stage 1: Read the MPS and make sure the size of the constraint matrix
303           is correct.
304*/
305    CbcModel *model = new CbcModel(saveModel) ;
306
307    std::string fn = dirMiplib+mpsName[m] ;
308    if (!CbcTestMpsFile(fn)) {
309      std::cout << "ERROR: Cannot find MPS file " << fn << "." << std::endl ;
310      continue;
311    }
312    model->solver()->readMps(fn.c_str(), "") ;
313    assert(model->getNumRows() == nRows[m]) ;
314    assert(model->getNumCols() == nCols[m]) ;
315
316    // Careful! We're initialising for the benefit of other code.
317    CoinDrand48(true,1234567);
318    //printf("RAND1 %g %g\n",CoinDrand48(true,1234567),model->randomNumberGenerator()->randomDouble());
319    //printf("RAND1 %g\n",CoinDrand48(true,1234567));
320
321    // Higher limits for the serious problems.
322    int testMaximumNodes = 200000;
323    if (hiSet > 1)
324        testMaximumNodes = 20000000;
325    if (model->getMaximumNodes() > testMaximumNodes) {
326        model->setMaximumNodes(testMaximumNodes);
327    }
328/*
329  Stage 2: Call solver to solve the problem.
330*/
331
332#   ifdef CLP_FACTORIZATION_INSTRUMENT
333    extern double factorization_instrument(int type);
334    double facTime1 = factorization_instrument(0);
335    std::cout
336      << "Factorization - initial solve " << facTime1 << " seconds."
337      << std::endl ;
338      timeTakenFac += facTime1;
339#   endif
340
341    double startTime = CoinCpuTime()+CoinCpuTimeJustChildren();
342
343    // Setup specific to clp
344    OsiClpSolverInterface *siClp =
345        dynamic_cast<OsiClpSolverInterface *>(model->solver()) ;
346    ClpSimplex *modelC = NULL;
347    if (siClp) {
348      modelC = siClp->getModelPtr();
349      ClpMatrixBase * matrix = modelC->clpMatrix();
350      ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
351      if (stuff && stuff[9] && clpMatrix) {
352        // vector matrix!
353        clpMatrix->makeSpecialColumnCopy();
354      }
355
356#     ifdef JJF_ZERO
357      if (clpMatrix) {
358        int numberRows = clpMatrix->getNumRows();
359        int numberColumns = clpMatrix->getNumCols();
360        double * elements = clpMatrix->getMutableElements();
361        const int * row = clpMatrix->getIndices();
362        const CoinBigIndex * columnStart = clpMatrix->getVectorStarts();
363        const int * columnLength = clpMatrix->getVectorLengths();
364        double * smallest = new double [numberRows];
365        double * largest = new double [numberRows];
366        char * flag = new char [numberRows];
367        CoinZeroN(flag, numberRows);
368        for (int i = 0; i < numberRows; i++) {
369            smallest[i] = COIN_DBL_MAX;
370            largest[i] = 0.0;
371        }
372        for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
373            bool isInteger = modelC->isInteger(iColumn);
374            CoinBigIndex j;
375            for (j = columnStart[iColumn];
376                    j < columnStart[iColumn] + columnLength[iColumn]; j++) {
377                int iRow = row[j];
378                double value = fabs(elements[j]);
379                if (!isInteger)
380                    flag[iRow] = 1;
381                smallest[iRow] = CoinMin(smallest[iRow], value);
382                largest[iRow] = CoinMax(largest[iRow], value);
383            }
384        }
385        double * rowLower = modelC->rowLower();
386        double * rowUpper = modelC->rowUpper();
387        bool changed = false;
388        for (int i = 0; i < numberRows; i++) {
389            if (flag[i] && smallest[i] > 10.0 && false) {
390                smallest[i] = 1.0 / smallest[i];
391                if (rowLower[i] > -1.0e20)
392                    rowLower[i] *= smallest[i];
393                if (rowUpper[i] < 1.0e20)
394                    rowUpper[i] *= smallest[i];
395                changed = true;
396            } else {
397                smallest[i] = 0.0;
398            }
399        }
400        if (changed) {
401            printf("SCALED\n");
402            for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
403                CoinBigIndex j;
404                for (j = columnStart[iColumn];
405                        j < columnStart[iColumn] + columnLength[iColumn]; j++) {
406                    int iRow = row[j];
407                    if (smallest[iRow])
408                        elements[j] *= smallest[iRow];
409                }
410            }
411        }
412        delete [] smallest;
413        delete [] largest;
414        delete [] flag;
415      }
416#     endif    // JJF_ZERO
417
418      model->checkModel();
419      modelC->tightenPrimalBounds(0.0, 0, true);
420      model->initialSolve();
421      if (modelC->dualBound() == 1.0e10) {
422        // user did not set - so modify
423        // get largest scaled away from bound
424        ClpSimplex temp = *modelC;
425        temp.dual(0, 7);
426        double largestScaled = 1.0e-12;
427        double largest = 1.0e-12;
428        int numberRows = temp.numberRows();
429        const double * rowPrimal = temp.primalRowSolution();
430        const double * rowLower = temp.rowLower();
431        const double * rowUpper = temp.rowUpper();
432        const double * rowScale = temp.rowScale();
433        int iRow;
434        for (iRow = 0; iRow < numberRows; iRow++) {
435          double value = rowPrimal[iRow];
436          double above = value - rowLower[iRow];
437          double below = rowUpper[iRow] - value;
438          if (above < 1.0e12) {
439              largest = CoinMax(largest, above);
440          }
441          if (below < 1.0e12) {
442              largest = CoinMax(largest, below);
443          }
444          if (rowScale) {
445              double multiplier = rowScale[iRow];
446              above *= multiplier;
447              below *= multiplier;
448          }
449          if (above < 1.0e12) {
450              largestScaled = CoinMax(largestScaled, above);
451          }
452          if (below < 1.0e12) {
453              largestScaled = CoinMax(largestScaled, below);
454          }
455        }
456
457        int numberColumns = temp.numberColumns();
458        const double * columnPrimal = temp.primalColumnSolution();
459        const double * columnLower = temp.columnLower();
460        const double * columnUpper = temp.columnUpper();
461        const double * columnScale = temp.columnScale();
462        int iColumn;
463        for (iColumn = 0; iColumn < numberColumns; iColumn++) {
464          double value = columnPrimal[iColumn];
465          double above = value - columnLower[iColumn];
466          double below = columnUpper[iColumn] - value;
467          if (above < 1.0e12) {
468              largest = CoinMax(largest, above);
469          }
470          if (below < 1.0e12) {
471              largest = CoinMax(largest, below);
472          }
473          if (columnScale) {
474              double multiplier = 1.0 / columnScale[iColumn];
475              above *= multiplier;
476              below *= multiplier;
477          }
478          if (above < 1.0e12) {
479              largestScaled = CoinMax(largestScaled, above);
480          }
481          if (below < 1.0e12) {
482              largestScaled = CoinMax(largestScaled, below);
483          }
484        }
485        std::cout << "Largest (scaled) away from bound " << largestScaled
486                  << " unscaled " << largest << std::endl;
487#       ifdef JJF_ZERO
488        modelC->setDualBound(CoinMax(1.0001e8,
489                                 CoinMin(1000.0*largestScaled,1.00001e10)));
490#       else
491        modelC->setDualBound(CoinMax(1.0001e9,
492                                 CoinMin(1000.0*largestScaled,1.0001e10)));
493#       endif
494      }
495    }    // end clp-specific setup
496/*
497  Cut passes: For small models (n < 500) always do 100 passes, if possible
498  (-100). For larger models, use minimum drop to stop (100, 20).
499*/
500    model->setMinimumDrop(CoinMin(5.0e-2,
501                      fabs(model->getMinimizationObjValue())*1.0e-3+1.0e-4));
502    if (CoinAbs(model->getMaximumCutPassesAtRoot()) <= 100) {
503        if (model->getNumCols() < 500) {
504          model->setMaximumCutPassesAtRoot(-100);
505        } else if (model->getNumCols() < 5000) {
506          model->setMaximumCutPassesAtRoot(100);
507        } else {
508          model->setMaximumCutPassesAtRoot(20);
509        }
510    }
511    // If defaults then increase trust for small models
512    if (model->numberStrong() == 5 && model->numberBeforeTrust() == 10) {
513        int numberColumns = model->getNumCols();
514        if (numberColumns <= 50) {
515            model->setNumberBeforeTrust(1000);
516        } else if (numberColumns <= 100) {
517            model->setNumberBeforeTrust(100);
518        } else if (numberColumns <= 300) {
519            model->setNumberBeforeTrust(50);
520        }
521    }
522    //if (model->getNumCols()>=500) {
523    // switch off Clp stuff
524    //model->setFastNodeDepth(-1);
525    //}
526/*
527  Activate the row cut debugger, if requested.
528*/
529    if (rowCutDebugger[m] == true) {
530      std::string probName ;
531      model->solver()->getStrParam(OsiProbName,probName) ;
532      model->solver()->activateRowCutDebugger(probName.c_str()) ;
533      if (model->solver()->getRowCutDebugger())
534        std::cout << "Row cut debugger activated for " ;
535      else
536        std::cout << "Failed to activate row cut debugger for " ;
537      std::cout << mpsName[m] << "." << std::endl ;
538    }
539    setCutAndHeuristicOptions(*model) ;
540/*
541  More clp-specific setup.
542*/
543    if (siClp) {
544#     ifdef CLP_MULTIPLE_FACTORIZATIONS
545      if (!modelC->factorization()->isDenseOrSmall()) {
546          int denseCode = stuff ? static_cast<int> (stuff[4]) : -1;
547          int smallCode = stuff ? static_cast<int> (stuff[10]) : -1;
548          if (stuff && stuff[8] >= 1) {
549              if (denseCode < 0)
550                  denseCode = 40;
551              if (smallCode < 0)
552                  smallCode = 40;
553          }
554          if (denseCode > 0)
555              modelC->factorization()->setGoDenseThreshold(denseCode);
556          if (smallCode > 0)
557              modelC->factorization()->setGoSmallThreshold(smallCode);
558          if (denseCode >= modelC->numberRows()) {
559              //printf("problem going dense\n");
560              //modelC->factorization()->goDenseOrSmall(modelC->numberRows());
561          }
562      }
563#     endif
564      if (stuff && stuff[8] >= 1) {
565          if (modelC->numberColumns() + modelC->numberRows() <= 10000 &&
566                  model->fastNodeDepth() == -1)
567              model->setFastNodeDepth(-9);
568      }
569    }
570    //OsiObject * obj = new CbcBranchToFixLots(model,0.3,0.0,3,3000003);
571    //model->addObjects(1,&obj);
572    //delete obj;
573/*
574  Finally, the actual call to solve the MIP with branch-and-cut.
575*/
576    model->branchAndBound();
577
578#   ifdef CLP_FACTORIZATION_INSTRUMENT
579    double facTime = factorization_instrument(0);
580    std::cout << "Factorization " << facTime << " seconds." << std::endl ,
581    timeTakenFac += facTime;
582#   endif
583
584/*
585  Stage 3: Do the statistics and check the answer.
586*/
587    double timeOfSolution = CoinCpuTime()+CoinCpuTimeJustChildren()-startTime;
588    std::cout
589      << "Cuts at root node changed objective from "
590      << model->getContinuousObjective() << " to "
591      << model->rootObjectiveAfterCuts() << std::endl ;
592    int numberGenerators = model->numberCutGenerators();
593    for (int iGenerator = 0 ; iGenerator < numberGenerators ; iGenerator++) {
594      CbcCutGenerator *generator = model->cutGenerator(iGenerator);
595#     ifdef CLIQUE_ANALYSIS
596#     ifndef CLP_INVESTIGATE
597      CglImplication *implication =
598            dynamic_cast<CglImplication*>(generator->generator());
599      if (implication) continue;
600#     endif
601#     endif
602      std::cout
603        << generator->cutGeneratorName() << " was tried "
604        << generator->numberTimesEntered() << " times and created "
605        << generator->numberCutsInTotal() << " cuts of which "
606        << generator->numberCutsActive()
607        << " were active after adding rounds of cuts";
608      if (generator->timing())
609        std::cout << " (" << generator->timeInCutGenerator() << " seconds)" ;
610      std::cout << "." << std::endl;
611    }
612    std::cout
613      << model->getNumberHeuristicSolutions()
614      << " solutions found by heuristics." << std::endl ;
615    int numberHeuristics = model->numberHeuristics();
616    for (int iHeuristic = 0 ; iHeuristic < numberHeuristics ; iHeuristic++) {
617      CbcHeuristic *heuristic = model->heuristic(iHeuristic);
618      if (heuristic->numRuns()) {
619        std::cout
620          << heuristic->heuristicName() << " was tried "
621          << heuristic->numRuns() << " times out of "
622          << heuristic->numCouldRun() << " and created "
623          << heuristic->numberSolutionsFound() << " solutions." << std::endl ;
624        }
625    }
626/*
627  Check for the correct answer.
628*/
629    if (!model->status()) {
630
631      double objActual = model->getObjValue() ;
632      double objExpect = objValue[m] ;
633      double tolerance = CoinMin(fabs(objActual),fabs(objExpect)) ;
634      tolerance = CoinMax(1.0e-5,1.0e-5*tolerance) ;
635      //CoinRelFltEq eq(1.0e-3) ;
636
637      std::cout
638        << "cbc_clp (" << mpsName[m] << ") "
639        << std::setprecision(10) << objActual ;
640      if (fabs(objActual-objExpect) < tolerance) {
641        std::cout << std::setprecision(dfltPrecision) << "; okay" ;
642          numProbSolved++;
643      } else  {
644        std::cout
645          << " != " << objExpect << std::setprecision(dfltPrecision)
646          << "; error = " << fabs(objExpect-objActual) ;
647          numberFailures++;
648        //#ifdef COIN_DEVELOP
649        //abort();
650        //#endif
651      }
652    } else {
653      std::cout
654        << "cbc_clp (" << mpsName[m] << ") status not optimal; "
655        << "assuming too many nodes" ;
656    }
657    timeTaken += timeOfSolution;
658    std::cout
659      << " -- (" << model->getNodeCount() << " n / "
660      << model->getIterationCount() << " i / "
661      << timeOfSolution << " s) (subtotal " << timeTaken << " s)"
662      << std::endl;
663    delete model;
664  }
665/*
666  End main loop on MPS problems. Print a summary and calculate the return
667  value.
668*/
669  int returnCode = 0;
670  std::cout
671    << "cbc_clp solved " << numProbSolved << " out of " << numberAttempts;
672  int numberOnNodes = numberAttempts-numProbSolved-numberFailures;
673  if (numberFailures || numberOnNodes) {
674    if (numberOnNodes) {
675      std::cout << " (" << numberOnNodes << " stopped on nodes)";
676      returnCode = numberOnNodes;
677    }
678    if (numberFailures) {
679      std::cout << " (" << numberFailures << " gave bad answer!)";
680      returnCode += 100*numberFailures;
681    }
682  }
683  std::cout
684    << " and took " << timeTaken << " seconds." << std::endl;
685
686  if (testSwitch == -2) {
687    if (numberFailures || numberOnNodes) {
688      std::cout << "****** Unit Test failed." << std::endl ;
689      std::cerr << "****** Unit Test failed." << std::endl ;
690    } else {
691      std::cerr << "****** Unit Test succeeded." << std::endl ;
692    }
693  }
694# ifdef CLP_FACTORIZATION_INSTRUMENT
695  std::cout
696    << "Total factorization time " << timeTakenFac << "seconds." << std::endl ;
697# endif
698  return (returnCode) ;
699}
700
Note: See TracBrowser for help on using the repository browser.