source: trunk/Cbc/examples/repeat.cpp @ 1173

Last change on this file since 1173 was 1173, checked in by forrest, 11 years ago

add $id

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 10.2 KB
Line 
1/* $Id: repeat.cpp 1173 2009-06-04 09:44:10Z forrest $ */
2// Copyright (C) 2005, International Business Machines
3// Corporation and others.  All Rights Reserved.
4#if defined(_MSC_VER)
5// Turn off compiler warning about long names
6#  pragma warning(disable:4786)
7#endif
8
9#include <cassert>
10#include <iomanip>
11
12
13// For Branch and bound
14#include "OsiSolverInterface.hpp"
15#include "CbcModel.hpp"
16#include "CbcStrategy.hpp"
17#include "CbcBranchUser.hpp"
18#include "CbcCompareUser.hpp"
19#include "CbcCutGenerator.hpp"
20#include "CbcHeuristicLocal.hpp"
21#include "OsiClpSolverInterface.hpp"
22
23// Cuts
24
25#include "CglGomory.hpp"
26#include "CglProbing.hpp"
27#include "CglKnapsackCover.hpp"
28#include "CglRedSplit.hpp"
29#include "CglClique.hpp"
30#include "CglFlowCover.hpp"
31#include "CglMixedIntegerRounding2.hpp"
32// Preprocessing
33#include "CglPreProcess.hpp"
34
35// Heuristics
36
37#include "CbcHeuristic.hpp"
38
39#include  "CoinTime.hpp"
40
41//#############################################################################
42
43
44/************************************************************************
45
46This main program reads in an integer model from an mps file.
47
48It then sets up some Cgl cut generators and calls branch and cut.
49
50This is designed to show two things :-
51
521) Use of CbcStrategy to do preprocessing - this should be cleaner than old way
532) Looping round modifying solver and redoing branch and cut.  This uses resetToReferenceCopy
54   which resets stuff like the cutoff and any memory of previous branch and cut.
55
56************************************************************************/
57
58
59int main (int argc, const char *argv[])
60{
61
62  // Define your favorite OsiSolver
63 
64  OsiClpSolverInterface solver1;
65
66  // Read in model using argv[1]
67  // and assert that it is a clean model
68  std::string mpsFileName = "../../Data/Sample/p0033.mps";
69  if (argc>=2) mpsFileName = argv[1];
70  int numMpsReadErrors = solver1.readMps(mpsFileName.c_str(),"");
71  assert(numMpsReadErrors==0);
72  double time1 = CoinCpuTime();
73
74  /* Options are:
75     preprocess to do preprocessing
76     time in minutes
77     if 2 parameters and numeric taken as time
78  */
79  bool preProcess=false;
80  double minutes=-1.0;
81  int nGoodParam=0;
82  for (int iParam=2; iParam<argc;iParam++) {
83    if (!strcmp(argv[iParam],"preprocess")) {
84      preProcess=true;
85      nGoodParam++;
86    } else if (!strcmp(argv[iParam],"time")) {
87      if (iParam+1<argc&&isdigit(argv[iParam+1][0])) {
88        minutes=atof(argv[iParam+1]);
89        if (minutes>=0.0) {
90          nGoodParam+=2;
91          iParam++; // skip time
92        }
93      }
94    }
95  }
96  if (nGoodParam==0&&argc==3&&isdigit(argv[2][0])) {
97    // If time is given then stop after that number of minutes
98    minutes = atof(argv[2]);
99    if (minutes>=0.0) 
100      nGoodParam=1;
101  }
102  if (nGoodParam!=argc-2&&argc>=2) {
103    printf("Usage <file> [preprocess] [time <minutes>] or <file> <minutes>\n");
104    exit(1);
105  }
106  //solver1.getModelPtr()->setLogLevel(0);
107  solver1.messageHandler()->setLogLevel(0);
108  solver1.initialSolve();
109  // Reduce printout
110  solver1.setHintParam(OsiDoReducePrint,true,OsiHintTry);
111  CbcModel model(solver1);
112  model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
113  // Set up some cut generators and defaults
114  // Probing first as gets tight bounds on continuous
115
116  CglProbing generator1;
117  generator1.setUsingObjective(true);
118  generator1.setMaxPass(1);
119  generator1.setMaxPassRoot(5);
120  // Number of unsatisfied variables to look at
121  generator1.setMaxProbe(10);
122  generator1.setMaxProbeRoot(1000);
123  // How far to follow the consequences
124  generator1.setMaxLook(50);
125  generator1.setMaxLookRoot(500);
126  // Only look at rows with fewer than this number of elements
127  generator1.setMaxElements(200);
128  generator1.setRowCuts(3);
129
130  CglGomory generator2;
131  // try larger limit
132  generator2.setLimit(300);
133
134  CglKnapsackCover generator3;
135
136  CglRedSplit generator4;
137  // try larger limit
138  generator4.setLimit(200);
139
140  CglClique generator5;
141  generator5.setStarCliqueReport(false);
142  generator5.setRowCliqueReport(false);
143
144  CglMixedIntegerRounding2 mixedGen;
145  CglFlowCover flowGen;
146 
147  // Add in generators
148  // Experiment with -1 and -99 etc
149  model.addCutGenerator(&generator1,-1,"Probing");
150  model.addCutGenerator(&generator2,-1,"Gomory");
151  model.addCutGenerator(&generator3,-1,"Knapsack");
152  // model.addCutGenerator(&generator4,-1,"RedSplit");
153  model.addCutGenerator(&generator5,-1,"Clique");
154  model.addCutGenerator(&flowGen,-1,"FlowCover");
155  model.addCutGenerator(&mixedGen,-1,"MixedIntegerRounding");
156  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (model.solver());
157  // go faster stripes
158  if (osiclp) {
159    // Turn this off if you get problems
160    // Used to be automatically set
161    osiclp->setSpecialOptions(128);
162    if(osiclp->getNumRows()<300&&osiclp->getNumCols()<500) {
163      //osiclp->setupForRepeatedUse(2,1);
164      osiclp->setupForRepeatedUse(0,1);
165    }
166  } 
167  // Uncommenting this should switch off most CBC messages
168  //model.messagesPointer()->setDetailMessages(10,5,5000);
169  // Allow rounding heuristic
170
171  CbcRounding heuristic1(model);
172  model.addHeuristic(&heuristic1);
173
174  // And local search when new solution found
175
176  CbcHeuristicLocal heuristic2(model);
177  model.addHeuristic(&heuristic2);
178
179  // Redundant definition of default branching (as Default == User)
180  CbcBranchUserDecision branch;
181  model.setBranchingMethod(&branch);
182
183  // Definition of node choice
184  CbcCompareUser compare;
185  model.setNodeComparison(compare);
186
187  // Do initial solve to continuous
188  model.initialSolve();
189
190  // Could tune more
191  double objValue = model.solver()->getObjSense()*model.solver()->getObjValue();
192  double minimumDropA=CoinMin(1.0,fabs(objValue)*1.0e-3+1.0e-4);
193  double minimumDrop= fabs(objValue)*1.0e-4+1.0e-4;
194  printf("min drop %g (A %g)\n",minimumDrop,minimumDropA);
195  model.setMinimumDrop(minimumDrop);
196
197  if (model.getNumCols()<500)
198    model.setMaximumCutPassesAtRoot(-100); // always do 100 if possible
199  else if (model.getNumCols()<5000)
200    model.setMaximumCutPassesAtRoot(100); // use minimum drop
201  else
202    model.setMaximumCutPassesAtRoot(20);
203  model.setMaximumCutPasses(10);
204  //model.setMaximumCutPasses(2);
205
206  // Switch off strong branching if wanted
207  // model.setNumberStrong(0);
208  // Do more strong branching if small
209  if (model.getNumCols()<5000)
210    model.setNumberStrong(10);
211  model.setNumberStrong(20);
212  //model.setNumberStrong(5);
213  model.setNumberBeforeTrust(5);
214  //model.setSizeMiniTree(2);
215
216  model.solver()->setIntParam(OsiMaxNumIterationHotStart,100);
217
218  // If time is given then stop after that number of minutes
219  if (minutes>=0.0) {
220    std::cout<<"Stopping after "<<minutes<<" minutes"<<std::endl;
221    model.setDblParam(CbcModel::CbcMaximumSeconds,60.0*minutes);
222  }
223  // Switch off most output
224  if (model.getNumCols()<3000) {
225    model.messageHandler()->setLogLevel(1);
226    //model.solver()->messageHandler()->setLogLevel(0);
227  } else {
228    model.messageHandler()->setLogLevel(2);
229    model.solver()->messageHandler()->setLogLevel(1);
230  }
231  // Default strategy will leave cut generators as they exist already
232  // so cutsOnlyAtRoot (1) ignored
233  // numberStrong (2) is 5 (default)
234  // numberBeforeTrust (3) is 5 (default is 0)
235  // printLevel (4) defaults (0)
236  CbcStrategyDefault strategy(true,5,5);
237  // Set up pre-processing to find sos if wanted
238  if (preProcess)
239    strategy.setupPreProcessing(2);
240  model.setStrategy(strategy);
241
242  // Go round adding cuts to cutoff last solution
243  // Stop after finding 20 best solutions
244  for (int iPass=0;iPass<20;iPass++) {
245    time1 = CoinCpuTime();
246    // Do complete search
247    model.branchAndBound();
248   
249    std::cout<<mpsFileName<<" took "<<CoinCpuTime()-time1<<" seconds, "
250             <<model.getNodeCount()<<" nodes with objective "
251             <<model.getObjValue()
252             <<(!model.status() ? " Finished" : " Not finished")
253             <<std::endl;
254    // Stop if infeasible
255    if (model.isProvenInfeasible())
256      break;
257    // Print solution if finished - we can't get names from Osi! - so get from OsiClp
258   
259    assert (model.getMinimizationObjValue()<1.0e50);
260    OsiSolverInterface * solver = model.solver();
261    int numberColumns = solver->getNumCols();
262   
263    const double * solution = model.bestSolution();
264    //const double * lower = solver->getColLower();
265    //const double * upper = solver->getColUpper();
266
267    // Get names from solver1 (as OsiSolverInterface may lose)
268    std::vector<std::string> columnNames = *solver1.getModelPtr()->columnNames();
269   
270    int iColumn;
271    std::cout<<std::setiosflags(std::ios::fixed|std::ios::showpoint)<<std::setw(14);
272   
273    std::cout<<"--------------------------------------"<<std::endl;
274    for (iColumn=0;iColumn<numberColumns;iColumn++) {
275      double value=solution[iColumn];
276      if (fabs(value)>1.0e-7&&solver->isInteger(iColumn)) 
277        std::cout<<std::setw(6)<<iColumn<<" "
278                 <<columnNames[iColumn]<<" "
279                 <<value
280          //<<" "<<lower[iColumn]<<" "<<upper[iColumn]
281                 <<std::endl;
282    }
283    std::cout<<"--------------------------------------"<<std::endl;
284 
285    std::cout<<std::resetiosflags(std::ios::fixed|std::ios::showpoint|std::ios::scientific);
286    /* Now add cut to reference copy.
287       resetting to reference copy also gets rid of best solution so we
288       should either save best solution, reset, add cut OR
289       add cut to reference copy then reset - this is doing latter
290    */
291    OsiSolverInterface * refSolver = model.referenceSolver();
292    const double * bestSolution = model.bestSolution();
293    const double * originalLower = refSolver->getColLower();
294    const double * originalUpper = refSolver->getColUpper();
295    CoinPackedVector cut;
296    double rhs = 1.0;
297    for (iColumn=0;iColumn<numberColumns;iColumn++) {
298      double value=bestSolution[iColumn];
299      if (solver->isInteger(iColumn)) {
300        // only works for 0-1 variables
301        assert (originalLower[iColumn]==0.0&&
302                originalUpper[iColumn]==1.0);
303        // double check integer
304        assert (fabs(floor(value+0.5)-value)<1.0e-5);
305        if (value>0.5) {
306          // at 1.0
307          cut.insert(iColumn,-1.0);
308          rhs -= 1.0;
309        } else {
310          // at 0.0
311          cut.insert(iColumn,1.0);
312        }
313      }
314    }
315    // now add cut
316    refSolver->addRow(cut,rhs,COIN_DBL_MAX);
317    model.resetToReferenceSolver();
318  }
319  return 0;
320}   
Note: See TracBrowser for help on using the repository browser.