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 "CbcModel.hpp"
13 #include "CbcBranchActual.hpp"
14 #include "OsiClpSolverInterface.hpp"
15 
16 // Time
17 #include "CoinTime.hpp"
18 
19 /************************************************************************
20 
21 This main program reads in an integer model from an mps file.
22 
23 It then tries to find SOS structure
24 
25 *************************************************************************/
main(int argc,const char * argv[])26 int main(int argc, const char *argv[])
27 {
28 
29   // Define your favorite OsiSolver
30 
31   OsiClpSolverInterface solver1;
32 
33   // Read in model using argv[1]
34   // and assert that it is a clean model
35   std::string mpsFileName;
36 #if defined(MIPLIB3DIR)
37   mpsFileName = MIPLIB3DIR "/10teams";
38 #else
39   if (argc < 2) {
40     fprintf(stderr, "Do not know where to find miplib3 MPS files.\n");
41     exit(1);
42   }
43 #endif
44   if (argc >= 2)
45     mpsFileName = argv[1];
46   int numMpsReadErrors = solver1.readMps(mpsFileName.c_str(), "");
47   if (numMpsReadErrors != 0) {
48     printf("%d errors reading MPS file\n", numMpsReadErrors);
49     return numMpsReadErrors;
50   }
51 
52   int iRow, iColumn;
53   int numberColumns = solver1.getNumCols();
54   int numberRows = solver1.getNumRows();
55   // get row copy
56   const CoinPackedMatrix *matrix = solver1.getMatrixByRow();
57   const double *element = matrix->getElements();
58   const int *column = matrix->getIndices();
59   const CoinBigIndex *rowStart = matrix->getVectorStarts();
60   const int *rowLength = matrix->getVectorLengths();
61   const double *rowLower = solver1.getRowLower();
62   const double *rowUpper = solver1.getRowUpper();
63   const double *columnLower = solver1.getColLower();
64 
65   // Look for possible SOS
66   int numberSOS = 0;
67   int *mark = new int[numberColumns];
68   CoinFillN(mark, numberColumns, -1);
69   for (iRow = 0; iRow < numberRows; iRow++) {
70     if (rowLower[iRow] == 1.0 && rowUpper[iRow] == 1.0) {
71       bool goodRow = true;
72       for (int j = rowStart[iRow]; j < rowStart[iRow] + rowLength[iRow]; j++) {
73         int iColumn = column[j];
74         if (element[j] != 1.0 || !solver1.isInteger(iColumn) || mark[iColumn] >= 0 || columnLower[iColumn]) {
75           goodRow = false;
76           break;
77         }
78       }
79       if (goodRow) {
80         // mark all
81         for (int j = rowStart[iRow]; j < rowStart[iRow] + rowLength[iRow]; j++) {
82           int iColumn = column[j];
83           mark[iColumn] = numberSOS;
84         }
85         numberSOS++;
86       }
87     }
88   }
89   std::cout << numberSOS << " SOS" << std::endl;
90   if (!numberSOS)
91     return 0;
92     /*  This example does not look to find the correct order.  SOS are much more
93       powerful if there is a genuine order e.g. size or time.
94 
95       There are two pieces of code here.
96       1) Leave integrality conditions and add SOS as extra.  They should have
97       a higher priority.
98       2) Take off integrality conditions and do as SOS of type 2.  This is artificial
99       in this case.
100   */
101     //#define SOS2
102 #ifndef SOS2
103   CbcModel model(solver1);
104   // Do sets and priorities
105   CbcObject **objects = new CbcObject *[numberSOS];
106   int numberIntegers = model.numberIntegers();
107   /* model may not have created objects
108      If none then create
109   */
110   if (!numberIntegers || !model.numberObjects()) {
111     model.findIntegers(true);
112     numberIntegers = model.numberIntegers();
113   }
114   int *priority = new int[numberSOS];
115   // Set SOS priorities high
116   CoinFillN(priority, numberSOS, 1);
117   // Set up SOS
118   int *which = new int[numberColumns];
119   for (int iSOS = 0; iSOS < numberSOS; iSOS++) {
120     int n = 0;
121     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
122       if (mark[iColumn] == iSOS)
123         which[n++] = iColumn;
124     }
125     // NULL uses 0,1,2 .. as weights
126     objects[iSOS] = new CbcSOS(&model, n, which, NULL, iSOS, 1);
127   }
128 #else
129   // take off integers
130   for (iColumn = 0; iColumn < numberColumns; iColumn++)
131     solver1.setContinuous(iColumn);
132   CbcModel model(solver1);
133   // Do sets and priorities
134   CbcObject **objects = new CbcObject *[numberSOS];
135   // Set up SOS
136   int *which = new int[numberColumns];
137   for (int iSOS = 0; iSOS < numberSOS; iSOS++) {
138     int n = 0;
139     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
140       if (mark[iColumn] == iSOS)
141         which[n++] = iColumn;
142     }
143     // NULL uses 0,1,2 .. as weights
144     objects[iSOS] = new CbcSOS(&model, n, which, NULL, iSOS, 2);
145   }
146 #endif
147   delete[] mark;
148   delete[] which;
149   model.addObjects(numberSOS, objects);
150   for (iColumn = 0; iColumn < numberSOS; iColumn++)
151     delete objects[iColumn];
152   delete[] objects;
153 #ifndef SOS2
154   model.passInPriorities(priority, true);
155   delete[] priority;
156 #endif
157 
158   // If time is given then stop after that number of minutes
159   if (argc > 2) {
160     int minutes = atoi(argv[2]);
161     std::cout << "Stopping after " << minutes << " minutes" << std::endl;
162     assert(minutes >= 0);
163     model.setDblParam(CbcModel::CbcMaximumSeconds, 60.0 * minutes);
164   }
165   // Switch off most output
166   model.solver()->setHintParam(OsiDoReducePrint, true, OsiHintTry);
167   if (model.getNumCols() < 3000) {
168     model.messageHandler()->setLogLevel(1);
169     //model.solver()->messageHandler()->setLogLevel(0);
170   } else {
171     model.messageHandler()->setLogLevel(2);
172     model.solver()->messageHandler()->setLogLevel(1);
173   }
174   model.messageHandler()->setLogLevel(1);
175 
176   double time1 = CoinCpuTime();
177 
178   // Do complete search
179 
180   model.branchAndBound();
181 
182   std::cout << mpsFileName << " took " << CoinCpuTime() - time1 << " seconds, "
183             << model.getNodeCount() << " nodes with objective "
184             << model.getObjValue()
185             << (!model.status() ? " Finished" : " Not finished")
186             << std::endl;
187 
188   // Print solution - we can't get names from Osi!
189 
190   if (model.getMinimizationObjValue() < 1.0e50) {
191     int numberColumns = model.solver()->getNumCols();
192 
193     const double *solution = model.solver()->getColSolution();
194 
195     int iColumn;
196     std::cout << std::setiosflags(std::ios::fixed | std::ios::showpoint) << std::setw(14);
197 
198     std::cout << "--------------------------------------" << std::endl;
199     for (iColumn = 0; iColumn < numberColumns; iColumn++) {
200       double value = solution[iColumn];
201 #ifndef SOS2
202       if (fabs(value) > 1.0e-7 && model.solver()->isInteger(iColumn))
203         std::cout << std::setw(6) << iColumn << " " << value << std::endl;
204 #else
205       if (fabs(value) > 1.0e-7)
206         std::cout << std::setw(6) << iColumn << " " << value << std::endl;
207 #endif
208     }
209     std::cout << "--------------------------------------" << std::endl;
210 
211     std::cout << std::resetiosflags(std::ios::fixed | std::ios::showpoint | std::ios::scientific);
212   }
213   return 0;
214 }
215