1 // $Id$
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 "OsiClpSolverInterface.hpp"
19
20 // Cuts
21
22 #include "CglGomory.hpp"
23 #include "CglProbing.hpp"
24 #include "CglKnapsackCover.hpp"
25 #include "CglRedSplit.hpp"
26 #include "CglClique.hpp"
27 #include "CglFlowCover.hpp"
28 #include "CglMixedIntegerRounding2.hpp"
29 // Preprocessing
30 #include "CglPreProcess.hpp"
31
32 // Heuristics
33
34 #include "CbcHeuristic.hpp"
35
36 #include "CoinTime.hpp"
37
38 //#############################################################################
39
40 /************************************************************************
41
42 This main program reads in an integer model from an mps file.
43
44 It then sets up some Cgl cut generators and calls branch and cut.
45
46 Then it uses solution as hot start
47
48 ************************************************************************/
49
main(int argc,const char * argv[])50 int main(int argc, const char *argv[])
51 {
52
53 // Define your favorite OsiSolver
54
55 OsiClpSolverInterface solver1;
56
57 // Read in model using argv[1]
58 // and assert that it is a clean model
59 std::string mpsFileName;
60 #if defined(SAMPLEDIR)
61 mpsFileName = SAMPLEDIR "/p0033.mps";
62 #else
63 if (argc < 2) {
64 fprintf(stderr, "Do not know where to find sample MPS files.\n");
65 exit(1);
66 }
67 #endif
68 if (argc >= 2)
69 mpsFileName = argv[1];
70 int numMpsReadErrors = solver1.readMps(mpsFileName.c_str(), "");
71 if (numMpsReadErrors != 0) {
72 printf("%d errors reading MPS file\n", numMpsReadErrors);
73 return numMpsReadErrors;
74 }
75 double time1 = CoinCpuTime();
76
77 /* Options are:
78 preprocess to do preprocessing
79 time in minutes
80 if 2 parameters and numeric taken as time
81 */
82 bool preProcess = false;
83 double minutes = -1.0;
84 int nGoodParam = 0;
85 for (int iParam = 2; iParam < argc; iParam++) {
86 if (!strcmp(argv[iParam], "preprocess")) {
87 preProcess = true;
88 nGoodParam++;
89 } else if (!strcmp(argv[iParam], "time")) {
90 if (iParam + 1 < argc && isdigit(argv[iParam + 1][0])) {
91 minutes = atof(argv[iParam + 1]);
92 if (minutes >= 0.0) {
93 nGoodParam += 2;
94 iParam++; // skip time
95 }
96 }
97 }
98 }
99 if (nGoodParam == 0 && argc == 3 && isdigit(argv[2][0])) {
100 // If time is given then stop after that number of minutes
101 minutes = atof(argv[2]);
102 if (minutes >= 0.0)
103 nGoodParam = 1;
104 }
105 if (nGoodParam != argc - 2 && argc >= 2) {
106 printf("Usage <file> [preprocess] [time <minutes>] or <file> <minutes>\n");
107 exit(1);
108 }
109 //solver1.getModelPtr()->setLogLevel(0);
110 solver1.messageHandler()->setLogLevel(0);
111 solver1.initialSolve();
112 // Reduce printout
113 solver1.setHintParam(OsiDoReducePrint, true, OsiHintTry);
114 CbcModel model(solver1);
115 model.solver()->setHintParam(OsiDoReducePrint, true, OsiHintTry);
116 // Set up some cut generators and defaults
117 // Probing first as gets tight bounds on continuous
118
119 CglProbing generator1;
120 generator1.setUsingObjective(true);
121 generator1.setMaxPass(1);
122 generator1.setMaxPassRoot(5);
123 // Number of unsatisfied variables to look at
124 generator1.setMaxProbe(10);
125 generator1.setMaxProbeRoot(1000);
126 // How far to follow the consequences
127 generator1.setMaxLook(50);
128 generator1.setMaxLookRoot(500);
129 // Only look at rows with fewer than this number of elements
130 generator1.setMaxElements(200);
131 generator1.setRowCuts(3);
132
133 CglGomory generator2;
134 // try larger limit
135 generator2.setLimit(300);
136
137 CglKnapsackCover generator3;
138
139 CglRedSplit generator4;
140 // try larger limit
141 generator4.setLimit(200);
142
143 CglClique generator5;
144 generator5.setStarCliqueReport(false);
145 generator5.setRowCliqueReport(false);
146
147 CglMixedIntegerRounding2 mixedGen;
148 CglFlowCover flowGen;
149
150 // Add in generators
151 // Experiment with -1 and -99 etc
152 model.addCutGenerator(&generator1, -1, "Probing");
153 model.addCutGenerator(&generator2, -1, "Gomory");
154 model.addCutGenerator(&generator3, -1, "Knapsack");
155 // model.addCutGenerator(&generator4,-1,"RedSplit");
156 model.addCutGenerator(&generator5, -1, "Clique");
157 model.addCutGenerator(&flowGen, -1, "FlowCover");
158 model.addCutGenerator(&mixedGen, -1, "MixedIntegerRounding");
159 OsiClpSolverInterface *osiclp = dynamic_cast< OsiClpSolverInterface * >(model.solver());
160 // go faster stripes
161 if (osiclp) {
162 // Turn this off if you get problems
163 // Used to be automatically set
164 osiclp->setSpecialOptions(128);
165 if (osiclp->getNumRows() < 300 && osiclp->getNumCols() < 500) {
166 //osiclp->setupForRepeatedUse(2,1);
167 osiclp->setupForRepeatedUse(0, 1);
168 }
169 }
170 // Uncommenting this should switch off most CBC messages
171 //model.messagesPointer()->setDetailMessages(10,5,5000);
172 // Allow rounding heuristic
173
174 CbcRounding heuristic1(model);
175 model.addHeuristic(&heuristic1);
176
177 // Redundant definition of default branching (as Default == User)
178 CbcBranchUserDecision branch;
179 model.setBranchingMethod(&branch);
180
181 // Definition of node choice
182 CbcCompareUser compare;
183 model.setNodeComparison(compare);
184
185 // Do initial solve to continuous
186 model.initialSolve();
187
188 // Could tune more
189 double objValue = model.solver()->getObjSense() * model.solver()->getObjValue();
190 double minimumDropA = CoinMin(1.0, fabs(objValue) * 1.0e-3 + 1.0e-4);
191 double minimumDrop = fabs(objValue) * 1.0e-4 + 1.0e-4;
192 printf("min drop %g (A %g)\n", minimumDrop, minimumDropA);
193 model.setMinimumDrop(minimumDrop);
194
195 if (model.getNumCols() < 500)
196 model.setMaximumCutPassesAtRoot(-100); // always do 100 if possible
197 else if (model.getNumCols() < 5000)
198 model.setMaximumCutPassesAtRoot(100); // use minimum drop
199 else
200 model.setMaximumCutPassesAtRoot(20);
201 model.setMaximumCutPasses(10);
202 //model.setMaximumCutPasses(2);
203
204 // Switch off strong branching if wanted
205 // model.setNumberStrong(0);
206 // Do more strong branching if small
207 if (model.getNumCols() < 5000)
208 model.setNumberStrong(10);
209 model.setNumberStrong(20);
210 //model.setNumberStrong(5);
211 model.setNumberBeforeTrust(5);
212 //model.setSizeMiniTree(2);
213
214 model.solver()->setIntParam(OsiMaxNumIterationHotStart, 100);
215
216 // Switch off most output
217 if (model.getNumCols() < 3000) {
218 model.messageHandler()->setLogLevel(1);
219 //model.solver()->messageHandler()->setLogLevel(0);
220 } else {
221 model.messageHandler()->setLogLevel(2);
222 model.solver()->messageHandler()->setLogLevel(1);
223 }
224 model.messageHandler()->setLogLevel(6);
225 model.solver()->messageHandler()->setLogLevel(1);
226 // If time is given then stop after that number of minutes
227 if (minutes >= 0.0) {
228 std::cout << "Stopping after " << minutes << " minutes" << std::endl;
229 model.setDblParam(CbcModel::CbcMaximumSeconds, 60.0 * minutes);
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 int numberColumns = model.solver()->getNumCols();
242 double *bestSolution = NULL;
243 int *priority = new int[numberColumns];
244 // Do two passes
245 for (int iPass = 0; iPass < 2; iPass++) {
246 time1 = CoinCpuTime();
247 // Do hot start on second pass
248 if (bestSolution) {
249 model.setHotstartSolution(bestSolution, priority);
250 delete[] bestSolution;
251 bestSolution = NULL;
252 delete[] priority;
253 model.setMaximumNodes(40000);
254 } else {
255 model.setMaximumNodes(40000);
256 }
257 // Do complete search
258 model.branchAndBound();
259
260 std::cout << mpsFileName << " took " << CoinCpuTime() - time1 << " seconds, "
261 << model.getNodeCount() << " nodes with objective "
262 << model.getObjValue()
263 << (!model.status() ? " Finished" : " Not finished")
264 << std::endl;
265 // Print solution if finished - we can't get names from Osi! - so get from OsiClp
266
267 assert(model.getMinimizationObjValue() < 1.0e50);
268 OsiSolverInterface *solver = model.solver();
269 int numberColumns = solver->getNumCols();
270
271 const double *solution = solver->getColSolution();
272 // save solution
273 if (!iPass) {
274 bestSolution = CoinCopyOfArray(solution, numberColumns);
275 for (int i = 0; i < numberColumns; i++) {
276 if (solution[i] > 0.5)
277 priority[i] = -1;
278 else
279 priority[i] = 2;
280 }
281 }
282 //const double * lower = solver->getColLower();
283 //const double * upper = solver->getColUpper();
284
285 // Get names from solver1 (as OsiSolverInterface may lose)
286 std::vector< std::string > columnNames = *solver1.getModelPtr()->columnNames();
287
288 int iColumn;
289 std::cout << std::setiosflags(std::ios::fixed | std::ios::showpoint) << std::setw(14);
290
291 std::cout << "--------------------------------------" << std::endl;
292 for (iColumn = 0; iColumn < numberColumns; iColumn++) {
293 double value = solution[iColumn];
294 if (fabs(value) > 1.0e-7 && solver->isInteger(iColumn))
295 std::cout << std::setw(6) << iColumn << " "
296 << columnNames[iColumn] << " "
297 << value
298 //<<" "<<lower[iColumn]<<" "<<upper[iColumn]
299 << std::endl;
300 }
301 std::cout << "--------------------------------------" << std::endl;
302
303 std::cout << std::resetiosflags(std::ios::fixed | std::ios::showpoint | std::ios::scientific);
304 model.resetToReferenceSolver();
305 }
306 return 0;
307 }
308