1 //-----------------------------------------------------
2 // Simple example usage of the cut generation library.
3 //
4 // This sample program iteratively tightens a
5 // given formulation by adding conic cuts, then calls
6 // branch-and-bound to solve the tightened
7 // formulation.
8 //
9 // usage:
10 //   conic_cgl mpsFileName objectiveSense
11 // where:
12 //   mpsFileName: Name of an mps file (without the
13 //                file extension)
14 //   objectiveSense: min for minimization,
15 //                   max for maximization.
16 // example:
17 //   conic_cgl ../../Data/Sample/p0033 min
18 //-----------------------------------------------------
19 
20 // TODO:
21 //-> write OsiConicCuts
22 //-> write CglConicMIR
23 
24 // STDLIB headers
25 #include <cassert>
26 #include <iostream>
27 #include <string>
28 #include <cassert>
29 // CoinUtils headers
30 #include "CoinError.hpp"
31 #include "CoinWarmStartBasis.hpp"
32 // OSI headers
33 #include "OsiCuts.hpp"
34 #include "OsiClpSolverInterface.hpp"
35 // OSICONIC header
36 #include <OsiConicSolverInterface.hpp>
37 #include <OsiConicCuts.hpp>
38 // COLA headers
39 #include <ColaModel.hpp>
40 // CGL headers
41 #include "CglKnapsackCover.hpp"
42 #include "CglSimpleRounding.hpp"
43 #include "CglGMI.hpp"
44 #include "CglGomory.hpp"
45 #include "CglMixedIntegerRounding.hpp"
46 #include "CglMixedIntegerRounding2.hpp"
47 // Conic CGL headers
48 #include "CglConicMIR.hpp"
49 
50 using std::cerr;
51 using std::cout;
52 using std::endl;
53 using std::string;
54 
main(int argc,const char * argv[])55 int main(int argc, const char *argv[])
56 {
57   // If no parms specified then use these
58   string mpsFileName;
59 #if defined(SAMPLEDIR)
60   mpsFileName = SAMPLEDIR "/p0033.mps";
61 #else
62   if (argc == 1) {
63     fprintf(stderr, "Do not know where to find sample MPS files.\n");
64     exit(1);
65   }
66 #endif
67   string objSense = "min";
68 
69   // Make sure a file name and objective sense or nothing
70   // were specified
71   if ( argc!=1 && argc!=3 ) {
72     cerr <<"Incorrect number of command line parameters." <<endl;
73     cerr <<"  usage:" <<endl;
74     cerr <<"    "<<argv[0] <<" mpsFileName objectiveSense" <<endl;
75     cerr <<"  where:" <<endl;
76     cerr <<"    mpsFileName: Name of an mps file" <<endl;
77     cerr <<"                 without \".mps\" file extension" <<endl;
78     cerr <<"    objectiveSense: min for minimization," <<endl;
79     cerr <<"                    max for maximization." <<endl;
80     return 1;
81   }
82 
83   // Make sure valid objective sense was specified
84   if (argc==3) {
85     mpsFileName = argv[1];
86     objSense = argv[2];
87     if( objSense!="min" && objSense!="max" ){
88       cerr <<"Unrecognized objective sense specifed" <<endl;
89       cerr <<"  specified value: \"" <<argv[2] <<"\"" <<endl;
90 
91       cerr <<"  valid values: \"min\" for minimization," <<endl;
92       cerr <<"                \"max\" for maximization." <<endl;
93       return 1;
94     }
95   }
96 
97   try {
98     // Instantiate a specific solver interface
99     OsiConicSolverInterface * si = new ColaModel();
100     // Read file describing problem
101     si->readMps(mpsFileName.c_str(),"mps");
102     // Set objective min to max
103     // First make sure min or max were specified
104     if (objSense=="min")
105       si->setObjSense(1.0);
106     else
107       si->setObjSense(-1.0);
108     // Solve continuous problem
109     si->initialSolve();
110     // Original number of rows (so we can take off inactive cuts)
111     int numberRows = si->getNumRows();
112     // Save the orig lp/soco relaxation value for
113     // comparisons later
114     double origLpObj = si->getObjValue();
115     // Track the total number of cuts applied
116     int totalNumberApplied = 0;
117     // Instantiate cut generators
118     CglKnapsackCover cg1;
119     CglSimpleRounding cg2;
120     CglGMI cg3;
121     CglGomory cg4;
122     CglMixedIntegerRounding cg5;
123     CglMixedIntegerRounding2 cg6;
124     //CglConicMIR cg7;
125     //---------------------------------------------------
126     // Keep applying cuts until
127     //   1. no more cuts are generated
128     // or
129     //   2. the objective function value doesn't change
130     //---------------------------------------------------
131     bool equalObj;
132     CoinRelFltEq eq(0.0001);
133     OsiSolverInterface::ApplyCutsReturnCode acRc;
134     double obj;
135     do {
136       // Get current solution value
137       obj = si->getObjValue();
138       // Generate and apply cuts
139       OsiCuts milp_cuts;
140       OsiSolverInterface * ssi = dynamic_cast<OsiSolverInterface*>(si);
141       cg1.generateCuts(*ssi,milp_cuts);
142       cg2.generateCuts(*ssi,milp_cuts);
143       cg3.generateCuts(*ssi,milp_cuts);
144       cg4.generateCuts(*ssi,milp_cuts);
145       cg5.generateCuts(*ssi,milp_cuts);
146       cg6.generateCuts(*ssi,milp_cuts);
147       //cg7.generateCuts(*si,cuts);
148       acRc = si->applyCuts(milp_cuts,0.0);
149       // Print applyCuts return code
150       cout <<endl <<endl;
151       cout <<milp_cuts.sizeCuts() <<" cuts were generated" <<endl;
152       cout <<"  " <<acRc.getNumInconsistent() <<" were inconsistent" <<endl;
153       cout <<"  " <<acRc.getNumInconsistentWrtIntegerModel()
154            <<" were inconsistent for this problem" <<endl;
155       cout <<"  " <<acRc.getNumInfeasible() <<" were infeasible" <<endl;
156       cout <<"  " <<acRc.getNumIneffective() <<" were ineffective" <<endl;
157       cout <<"  " <<acRc.getNumApplied() <<" were applied" <<endl;
158       cout <<endl <<endl;
159       // Increment the counter of total cuts applied
160       totalNumberApplied += acRc.getNumApplied();
161       // If no cuts were applied, then done
162       if (acRc.getNumApplied()==0)
163 	break;
164       // Resolve
165       si->resolve();
166       cout <<endl;
167       cout <<"After applying cuts, objective value changed from "
168 	   << obj << " to " << si->getObjValue() << endl << endl;
169       // Take off slack cuts
170       int numberRowsNow = si->getNumRows();
171       int * del = new int [numberRowsNow-numberRows];
172       const CoinWarmStartBasis* basis = dynamic_cast<const CoinWarmStartBasis*>(si->getWarmStart());
173       assert (basis);
174       int nDelete=0;
175       for (int i=numberRows;i<numberRowsNow;i++) {
176 	CoinWarmStartBasis::Status status = basis->getArtifStatus(i);
177 	if (status == CoinWarmStartBasis::basic)
178 	  del[nDelete++] = i;
179       }
180       delete basis;
181       if (nDelete) {
182 	si->deleteRows(nDelete,del);
183 	// should take zero iterations
184 	si->resolve();
185 	cout << nDelete << " rows deleted as basic - resolve took "
186 	     << si->getIterationCount() <<" iterations"
187 	     <<endl;
188       }
189       delete [] del;
190       // -----------------------------------------------
191       // Set Boolean flag to true if new objective is
192       // almost equal to prior value.
193       //
194       // The test is:
195       // abs(oldObj-newObj) <= 0.0001*(CoinMax(abs(oldObj),abs(newObj))+1.);
196       // see CoinRelFloatEqual.h
197       // -----------------------------------------------
198       equalObj = eq(si->getObjValue(), obj);
199     } while (!equalObj);
200     // Print total number of cuts applied,
201     // and total improvement in the lp objective value
202     cout <<endl <<endl;
203     cout << "----------------------------------------------------------"
204          <<endl;
205     cout << "Cut generation phase completed:" <<endl;
206     cout << "   " << totalNumberApplied << " cuts were applied in total,"
207          <<endl;
208     cout << "   changing the lp objective value from " << origLpObj
209          << " to " << si->getObjValue() <<endl;
210     cout << "----------------------------------------------------------"
211          <<endl;
212     cout <<endl <<endl;
213     // If you want to solve problem, change "#if 0" to "#if 1"
214 #if 0
215     // Solve MIP Problem
216     si->branchAndBound();
217     // Print the solution
218     cout << "The objective function value: " << si->getObjValue() <<endl;
219     const double * soln = si->getColSolution();
220     int i;
221     for ( i=0; i<si->getNumCols(); i ++ ) {
222       cout << " x[" << i << "] = " << soln[i] << endl;
223     }
224 #endif
225   }
226   catch (CoinError e) {
227     cout << e.className() << "::" << e.methodName() << " - " << e.message() << endl;
228   }
229   return 0;
230 }
231