source: stable/2.8/Cbc/src/unitTestClp.cpp @ 1883

Last change on this file since 1883 was 1883, checked in by stefan, 5 years ago

sync with trunk rev 1882

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