source: trunk/Cbc/examples/sample2.cpp

Last change on this file was 2469, checked in by unxusr, 5 months ago

formatting

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 11.3 KB
Line 
1// $Id: sample2.cpp 2469 2019-01-06 23:17:46Z forrest $
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#include <cassert>
7#include <iomanip>
8
9#include "CoinPragma.hpp"
10
11// For Branch and bound
12#include "OsiSolverInterface.hpp"
13#include "CbcModel.hpp"
14#include "CbcBranchUser.hpp"
15#include "CbcCompareUser.hpp"
16#include "CbcCutGenerator.hpp"
17#include "CbcStrategy.hpp"
18#include "CbcHeuristicLocal.hpp"
19#include "OsiClpSolverInterface.hpp"
20
21// Cuts
22
23#include "CglGomory.hpp"
24#include "CglProbing.hpp"
25#include "CglKnapsackCover.hpp"
26#include "CglRedSplit.hpp"
27#include "CglClique.hpp"
28#include "CglFlowCover.hpp"
29#include "CglMixedIntegerRounding2.hpp"
30// Preprocessing
31#include "CglPreProcess.hpp"
32
33// Heuristics
34
35#include "CbcHeuristic.hpp"
36
37#include "CoinTime.hpp"
38
39//#############################################################################
40
41/************************************************************************
42
43This main program reads in an integer model from an mps file.
44
45It then sets up some Cgl cut generators and calls branch and cut.
46
47Branching is simple binary branching on integer variables.
48
49Node selection is depth first until first solution is found and then
50based on objective and number of unsatisfied integer variables.
51In this example the functionality is the same as default but it is
52a user comparison function.
53
54Variable branching selection is on maximum minimum-of-up-down change
55after strong branching on 5 variables closest to 0.5.
56
57A simple rounding heuristic is used.
58
59Preprocessing may be selected in two ways - the second is preferred now
60
61************************************************************************/
62#define PREPROCESS 2
63
64int main(int argc, const char *argv[])
65{
66
67  // Define your favorite OsiSolver
68
69  OsiClpSolverInterface solver1;
70
71  // Read in model using argv[1]
72  // and assert that it is a clean model
73  std::string dirsep(1, CoinFindDirSeparator());
74  std::string mpsFileName;
75#if defined(SAMPLEDIR)
76  mpsFileName = SAMPLEDIR;
77  mpsFileName += dirsep + "p0033.mps";
78#else
79  if (argc < 2) {
80    fprintf(stderr, "Do not know where to find sample MPS files.\n");
81    exit(1);
82  }
83#endif
84  if (argc >= 2)
85    mpsFileName = argv[1];
86  int numMpsReadErrors = solver1.readMps(mpsFileName.c_str(), "");
87  if (numMpsReadErrors != 0) {
88    printf("%d errors reading MPS file\n", numMpsReadErrors);
89    return numMpsReadErrors;
90  }
91  double time1 = CoinCpuTime();
92
93  /* Options are:
94     preprocess to do preprocessing
95     time in minutes
96     if 2 parameters and numeric taken as time
97  */
98  bool preProcess = false;
99  double minutes = -1.0;
100  int nGoodParam = 0;
101  for (int iParam = 2; iParam < argc; iParam++) {
102    if (!strcmp(argv[iParam], "preprocess")) {
103      preProcess = true;
104      nGoodParam++;
105    } else if (!strcmp(argv[iParam], "time")) {
106      if (iParam + 1 < argc && isdigit(argv[iParam + 1][0])) {
107        minutes = atof(argv[iParam + 1]);
108        if (minutes >= 0.0) {
109          nGoodParam += 2;
110          iParam++; // skip time
111        }
112      }
113    }
114  }
115  if (nGoodParam == 0 && argc == 3 && isdigit(argv[2][0])) {
116    // If time is given then stop after that number of minutes
117    minutes = atof(argv[2]);
118    if (minutes >= 0.0)
119      nGoodParam = 1;
120  }
121  if (nGoodParam != argc - 2 && argc >= 2) {
122    printf("Usage <file> [preprocess] [time <minutes>] or <file> <minutes>\n");
123    exit(1);
124  }
125  solver1.initialSolve();
126  // Reduce printout
127  solver1.setHintParam(OsiDoReducePrint, true, OsiHintTry);
128  // See if we want preprocessing
129  OsiSolverInterface *solver2 = &solver1;
130#if PREPROCESS == 1
131  CglPreProcess process;
132  if (preProcess) {
133    /* Do not try and produce equality cliques and
134       do up to 5 passes */
135    solver2 = process.preProcess(solver1, false, 5);
136    if (!solver2) {
137      printf("Pre-processing says infeasible\n");
138      exit(2);
139    }
140    solver2->resolve();
141  }
142#endif
143  CbcModel model(*solver2);
144  model.solver()->setHintParam(OsiDoReducePrint, true, OsiHintTry);
145  // Set up some cut generators and defaults
146  // Probing first as gets tight bounds on continuous
147
148  CglProbing generator1;
149  generator1.setUsingObjective(true);
150  generator1.setMaxPass(1);
151  generator1.setMaxPassRoot(5);
152  // Number of unsatisfied variables to look at
153  generator1.setMaxProbe(10);
154  generator1.setMaxProbeRoot(1000);
155  // How far to follow the consequences
156  generator1.setMaxLook(50);
157  generator1.setMaxLookRoot(500);
158  // Only look at rows with fewer than this number of elements
159  generator1.setMaxElements(200);
160  generator1.setRowCuts(3);
161
162  CglGomory generator2;
163  // try larger limit
164  generator2.setLimit(300);
165
166  CglKnapsackCover generator3;
167
168  CglRedSplit generator4;
169  // try larger limit
170  generator4.setLimit(200);
171
172  CglClique generator5;
173  generator5.setStarCliqueReport(false);
174  generator5.setRowCliqueReport(false);
175
176  CglMixedIntegerRounding2 mixedGen;
177  CglFlowCover flowGen;
178
179  // Add in generators
180  // Experiment with -1 and -99 etc
181  model.addCutGenerator(&generator1, -1, "Probing");
182  model.addCutGenerator(&generator2, -1, "Gomory");
183  model.addCutGenerator(&generator3, -1, "Knapsack");
184  // model.addCutGenerator(&generator4,-1,"RedSplit");
185  model.addCutGenerator(&generator5, -1, "Clique");
186  model.addCutGenerator(&flowGen, -1, "FlowCover");
187  model.addCutGenerator(&mixedGen, -1, "MixedIntegerRounding");
188  // Say we want timings
189  int numberGenerators = model.numberCutGenerators();
190  int iGenerator;
191  for (iGenerator = 0; iGenerator < numberGenerators; iGenerator++) {
192    CbcCutGenerator *generator = model.cutGenerator(iGenerator);
193    generator->setTiming(true);
194  }
195  OsiClpSolverInterface *osiclp = dynamic_cast< OsiClpSolverInterface * >(model.solver());
196  // go faster stripes
197  if (osiclp) {
198    // Turn this off if you get problems
199    // Used to be automatically set
200    osiclp->setSpecialOptions(128);
201    if (osiclp->getNumRows() < 300 && osiclp->getNumCols() < 500) {
202      //osiclp->setupForRepeatedUse(2,0);
203      osiclp->setupForRepeatedUse(0, 0);
204    }
205  }
206  // Uncommenting this should switch off all CBC messages
207  // model.messagesPointer()->setDetailMessages(10,10000,NULL);
208  // Allow rounding heuristic
209
210  CbcRounding heuristic1(model);
211  model.addHeuristic(&heuristic1);
212
213  // And local search when new solution found
214
215  CbcHeuristicLocal heuristic2(model);
216  model.addHeuristic(&heuristic2);
217
218  // Redundant definition of default branching (as Default == User)
219  CbcBranchUserDecision branch;
220  model.setBranchingMethod(&branch);
221
222  // Definition of node choice
223  CbcCompareUser compare;
224  model.setNodeComparison(compare);
225
226  // Do initial solve to continuous
227  model.initialSolve();
228
229  // Could tune more
230  double objValue = model.solver()->getObjSense() * model.solver()->getObjValue();
231  double minimumDropA = CoinMin(1.0, fabs(objValue) * 1.0e-3 + 1.0e-4);
232  double minimumDrop = fabs(objValue) * 1.0e-4 + 1.0e-4;
233  printf("min drop %g (A %g)\n", minimumDrop, minimumDropA);
234  model.setMinimumDrop(minimumDrop);
235
236  if (model.getNumCols() < 500)
237    model.setMaximumCutPassesAtRoot(-100); // always do 100 if possible
238  else if (model.getNumCols() < 5000)
239    model.setMaximumCutPassesAtRoot(100); // use minimum drop
240  else
241    model.setMaximumCutPassesAtRoot(20);
242  model.setMaximumCutPasses(10);
243  //model.setMaximumCutPasses(2);
244
245  // Switch off strong branching if wanted
246  // model.setNumberStrong(0);
247  // Do more strong branching if small
248  if (model.getNumCols() < 5000)
249    model.setNumberStrong(10);
250  model.setNumberStrong(20);
251  //model.setNumberStrong(5);
252  model.setNumberBeforeTrust(5);
253
254  model.solver()->setIntParam(OsiMaxNumIterationHotStart, 100);
255
256  // If time is given then stop after that number of minutes
257  if (minutes >= 0.0) {
258    std::cout << "Stopping after " << minutes << " minutes" << std::endl;
259    model.setDblParam(CbcModel::CbcMaximumSeconds, 60.0 * minutes);
260  }
261  // Switch off most output
262  if (model.getNumCols() < 3000) {
263    model.messageHandler()->setLogLevel(1);
264    //model.solver()->messageHandler()->setLogLevel(0);
265  } else {
266    model.messageHandler()->setLogLevel(2);
267    model.solver()->messageHandler()->setLogLevel(1);
268  }
269  //model.messageHandler()->setLogLevel(2);
270  //model.solver()->messageHandler()->setLogLevel(2);
271  //model.setPrintFrequency(50);
272  //#define DEBUG_CUTS
273#ifdef DEBUG_CUTS
274  // Set up debugger by name (only if no preprocesing)
275  if (!preProcess) {
276    std::string problemName;
277    model.solver()->getStrParam(OsiProbName, problemName);
278    model.solver()->activateRowCutDebugger(problemName.c_str());
279  }
280#endif
281#if PREPROCESS == 2
282  // Default strategy will leave cut generators as they exist already
283  // so cutsOnlyAtRoot (1) ignored
284  // numberStrong (2) is 5 (default)
285  // numberBeforeTrust (3) is 5 (default is 0)
286  // printLevel (4) defaults (0)
287  CbcStrategyDefault strategy(true, 5, 5);
288  // Set up pre-processing to find sos if wanted
289  if (preProcess)
290    strategy.setupPreProcessing(2);
291  model.setStrategy(strategy);
292#endif
293  // Do complete search
294
295  model.branchAndBound();
296
297  std::cout << mpsFileName << " took " << CoinCpuTime() - time1 << " seconds, "
298            << model.getNodeCount() << " nodes with objective "
299            << model.getObjValue()
300            << (!model.status() ? " Finished" : " Not finished")
301            << std::endl;
302
303  // Print more statistics
304  std::cout << "Cuts at root node changed objective from " << model.getContinuousObjective()
305            << " to " << model.rootObjectiveAfterCuts() << std::endl;
306
307  for (iGenerator = 0; iGenerator < numberGenerators; iGenerator++) {
308    CbcCutGenerator *generator = model.cutGenerator(iGenerator);
309    std::cout << generator->cutGeneratorName() << " was tried "
310              << generator->numberTimesEntered() << " times and created "
311              << generator->numberCutsInTotal() << " cuts of which "
312              << generator->numberCutsActive() << " were active after adding rounds of cuts";
313    if (generator->timing())
314      std::cout << " ( " << generator->timeInCutGenerator() << " seconds)" << std::endl;
315    else
316      std::cout << std::endl;
317  }
318  // Print solution if finished - we can't get names from Osi! - so get from OsiClp
319
320  if (model.getMinimizationObjValue() < 1.0e50) {
321#if PREPROCESS == 1
322    // post process
323    OsiSolverInterface *solver;
324    if (preProcess) {
325      process.postProcess(*model.solver());
326      // Solution now back in solver1
327      solver = &solver1;
328    } else {
329      solver = model.solver();
330    }
331#else
332    OsiSolverInterface *solver = model.solver();
333#endif
334    int numberColumns = solver->getNumCols();
335
336    const double *solution = solver->getColSolution();
337
338    // Get names from solver1 (as OsiSolverInterface may lose)
339    std::vector< std::string > columnNames = *solver1.getModelPtr()->columnNames();
340
341    int iColumn;
342    std::cout << std::setiosflags(std::ios::fixed | std::ios::showpoint) << std::setw(14);
343
344    std::cout << "--------------------------------------" << std::endl;
345    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
346      double value = solution[iColumn];
347      if (fabs(value) > 1.0e-7 && solver->isInteger(iColumn))
348        std::cout << std::setw(6) << iColumn << " "
349                  << columnNames[iColumn] << " "
350                  << value << std::endl;
351    }
352    std::cout << "--------------------------------------" << std::endl;
353
354    std::cout << std::resetiosflags(std::ios::fixed | std::ios::showpoint | std::ios::scientific);
355  }
356  return 0;
357}
Note: See TracBrowser for help on using the repository browser.