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

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

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