generateDisjCuts.cpp
Go to the documentation of this file.
1 /* $Id: generateDisjCuts.cpp 945 2013-04-06 20:25:21Z stefan $
2  *
3  * Name: generateDisjCuts.cpp
4  * Author: Pietro Belotti
5  * Purpose: separation method for disjunctive cuts
6  *
7  * (C) Carnegie-Mellon University, 2008-10.
8  * This file is licensed under the Eclipse Public License (EPL)
9  */
10 
11 #include <vector>
12 #include "CoinTime.hpp"
13 
14 #include "CouenneCutGenerator.hpp"
15 #include "CouenneDisjCuts.hpp"
16 #include "CouenneProblem.hpp"
17 #include "CouenneInfeasCut.hpp"
18 
19 using namespace Ipopt;
20 using namespace Couenne;
21 
23 void CouenneDisjCuts::generateCuts (const OsiSolverInterface &si,
24  OsiCuts &cs,
25  const CglTreeInfo info)
26 #if CGL_VERSION_MAJOR == 0 && CGL_VERSION_MINOR <= 57
27  const
28 #endif
29 {
30 
31  // Check if cs contains only one cut and if it is of the form 1 <=
32  // x0 <= -1. That means a previous cut generator has determined that
33  // this node is infeasible and we shouldn't take the pain of running
34  // this CGL.
35 
36  if (isWiped (cs))
37  return;
38 
39  if (jnlst_ -> ProduceOutput (J_DETAILED, J_DISJCUTS))
40  printf ("--- generateDisjCuts: level = %d, pass = %d, intree = %d [%d]\n",
41  info.level, info.pass, info.inTree, depthStopSeparate_);
42 
43  if ((depthStopSeparate_ >= 0) && // if -1 no limit on depth
44  (info.level > depthStopSeparate_)) // check if too deep for adding these cuts
45  return;
46 
47  int nInitCuts = nrootcuts_;
48 
49  if ((info.level <= 0) && !(info.inTree)) {
50 
51  jnlst_ -> Printf (J_ERROR, J_COUENNE,
52  "Disjunctive cuts (root node): ");
53  fflush (stdout);
54  }
55 
56  double time = CoinCpuTime ();
57 
58  // use clone of solver interface
59  OsiSolverInterface *csi = si.clone ();
60 
61  int
62  initRowCuts = cs.sizeRowCuts (),
63  initColCuts = cs.sizeColCuts ();
64 
65  // get disjunctions /////////////////////////////
66  //
67  // start:
68  //
69  // consider problem Ax <= a_0, x in [l,u]
70  //
71  // // A and a_0 may be, depending on depth of BB tree and size of the
72  // // problem, from the linearization at the root node (in this case
73  // // a clone() is needed at the first iteration and branching rules
74  // // should be applied explicitly) or the current LP as passed by
75  // // si.
76  //
77  // Get set of disjunctions (in the form of OsiBranchingObjects) in
78  // number limited by depth, problem size, and sorted according to
79  // branching method (HotInfo if strong branching?)
80  //
81  // preprocessing /////////////////////////////////
82  //
83  // for each disjunction (x_i <= or >= x_i^d) {
84  //
85  // 1) apply left disj., bound tightening returns x in [l_1,u_1]
86  // 2) apply right disj., bound tightening returns x in [l_2,u_2]
87  //
88  // 3) if both infeasible, bail out
89  // 4) if one infeasible only, save column cut, apply tightening to
90  // si and goto start
91  // }
92  //
93  // CGLP //////////////////////////////////////////
94  //
95  // for each disjunction (x_i <= or >= x_i^d) {
96  //
97  // 5) get cuts Bx <= b_0 for left disj.
98  // 6) get cuts Cx <= c_0 for right disj.
99  //
100  // 7) if both infeasible, bail out
101  // 8) if one infeasible only, save column cut, apply tightening to
102  // si and continue
103  //
104  // 9) generate CGLP:
105  // consider subset of Ax <= a_0:
106  // a) active only, or
107  // b) {root linearization cuts} union {currently active}, or
108  // c) all cuts (!)
109  //
110  // 10) solve CGLP
111  //
112  // 11) add corresponding cut, possibly to CGLP as well?
113  // }
114 
115  int maxDisj = (initDisjNumber_ >= 0) ?
116  CoinMin ((int) (csi -> numberObjects () * initDisjPercentage_), initDisjNumber_) :
117  (int) (csi -> numberObjects () * initDisjPercentage_);
118 
119  // number of disjunctions to consider (branching objects)
120  numDisjunctions_ = (depthLevelling_ < 0 || info.level < depthLevelling_) ?
121  (int) (maxDisj) :
122  (int) (maxDisj / (2 + info.level - depthLevelling_));
123 
124  if (numDisjunctions_ < 1) numDisjunctions_ = 1;
125 
126  const int
127  exc_infeasible = 1,
128  exc_normal = 2,
129  max_iterations = 1;
130 
131  couenneCG_ -> Problem () -> domain () -> push (&si, &cs, false);
132 
133  std::vector <std::pair <OsiCuts *, OsiCuts *> > disjunctions;
134 
135  bool infeasNode = false;
136 
137  try {
138 
139  // get disjunctions (rows and cols) ////////////////////////////////////////////////////
140 
141  bool start_over;
142  int iterations = 0;
143 
144  do { // repeat as long as preprocessing or separation of convCuts
145  // shrink the whole problem
146 
147  ++iterations;
148 
149  start_over = false;
150 
151  // preprocess, get column cuts (disjoint bounding boxes)
152  int result = getDisjunctions (disjunctions, *csi, cs, info);
153 
154  if (result == COUENNE_INFEASIBLE) throw exc_infeasible; // fathom node
155  else if (result == COUENNE_TIGHTENED &&
156  iterations < max_iterations) start_over = true; // tightened node
157 
158  if (disjunctions.empty ())
159  throw exc_normal;
160 
161  // generate convexification cuts for each disjunction
162  for (std::vector <std::pair <OsiCuts *, OsiCuts *> >::iterator disjI = disjunctions.begin ();
163  disjI != disjunctions.end (); ++disjI) {
164 
165  if (CoinCpuTime () > couenneCG_ -> Problem () -> getMaxCpuTime ()) {
166  start_over = false;
167  break;
168  }
169 
170  // separate on single disjunction
171 
172  // left
173  result = separateWithDisjunction (disjI -> first, *csi, cs, info);
174  if (result == COUENNE_INFEASIBLE) throw exc_infeasible; // fathom node
175  else if (result == COUENNE_TIGHTENED && iterations < max_iterations) { // tightened node
176  start_over = true;
177  break;
178  }
179 
180  // right
181  result = separateWithDisjunction (disjI -> second, *csi, cs, info);
182  if (result == COUENNE_INFEASIBLE) throw exc_infeasible; // fathom node
183  else if (result == COUENNE_TIGHTENED && iterations < max_iterations) { // tightened node
184  start_over = true;
185  break;
186  }
187  }
188 
189  if (start_over) {
190 
191  for (std::vector <std::pair <OsiCuts *, OsiCuts *> >::iterator disjI = disjunctions.begin ();
192  disjI != disjunctions.end (); ++disjI) {
193  delete disjI -> first;
194  delete disjI -> second;
195  }
196 
197  disjunctions.erase (disjunctions.begin (), disjunctions.end ());
198  }
199 
200  if (!start_over && jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS))
201 
202  // generate convexification cuts for each disjunction
203  for (std::vector <std::pair <OsiCuts *, OsiCuts *> >::iterator disjI = disjunctions.begin ();
204  disjI != disjunctions.end (); ++disjI) {
205 
206  printf ("=========================== CUTS for the LEFT part\n");
207  for (int i=0; i<disjI->first->sizeColCuts (); i++) disjI->first->colCutPtr(i)->print();
208  for (int i=0; i<disjI->first->sizeRowCuts (); i++) disjI->first->rowCutPtr(i)->print();
209  printf ("=========================== CUTS for the RIGHT part\n");
210  for (int i=0; i<disjI->second->sizeColCuts (); i++) disjI->second->colCutPtr(i)->print();
211  for (int i=0; i<disjI->second->sizeRowCuts (); i++) disjI->second->rowCutPtr(i)->print();
212  printf ("===========================\n");
213  }
214 
215  } while (start_over && (iterations < max_iterations));
216 
217  // si contains the tightest bounding box. Use it to update
218  // CouenneProblem's bounds AND add to cs
219 
220  // already done above
221 
222  // CGLP //////////////////////////////////////////////////////////////////////////
223 
224  // maybe one last FBBT before big CGLP?
225 
226  // generate all cuts
227  if (generateDisjCuts (disjunctions, *csi, cs, info) == COUENNE_INFEASIBLE) // node can be fathomed
228  throw exc_infeasible;
229  }
230 
231  catch (int exception) {
232 
233  if (exception == exc_infeasible) { // add infeasible column cut 1 <= x_0 <= -1
234 
235  jnlst_ -> Printf (J_DETAILED, J_DISJCUTS, "--- Disjunctive Cut separator: infeasible node\n");
236  WipeMakeInfeas (cs);
237  infeasNode = true;
238  }
239  }
240 
241  // cleanup
242  for (std::vector <std::pair <OsiCuts *, OsiCuts *> >::iterator disjI = disjunctions.begin ();
243  disjI != disjunctions.end (); ++disjI) {
244 
245  delete disjI -> first;
246  delete disjI -> second;
247  }
248 
249  couenneCG_ -> Problem () -> domain () -> pop ();
250 
251  if (!infeasNode) {
252 
253  // tighten bounds of si based on those tightened of csi
254 
255  CoinPackedVector
256  tighterLower,
257  tighterUpper;
258 
259  const double
260  *oldLo = si. getColLower (), *newLo = csi -> getColLower (),
261  *oldUp = si. getColUpper (), *newUp = csi -> getColUpper ();
262 
263  int ncols = si.getNumCols ();
264 
265  for (int i=0; i<ncols; i++, newLo++, newUp++) {
266 
267  if (*newLo > *oldLo++ + COUENNE_EPS) tighterLower.insert (i, *newLo);
268  if (*newUp < *oldUp++ - COUENNE_EPS) tighterUpper.insert (i, *newUp);
269  }
270 
271  if ((tighterLower.getNumElements () > 0) ||
272  (tighterUpper.getNumElements () > 0)) {
273  OsiColCut tighter;
274  tighter.setLbs (tighterLower);
275  tighter.setUbs (tighterUpper);
276  if (jnlst_ -> ProduceOutput (J_DETAILED, J_DISJCUTS)) {
277  printf ("tightened bounds in disjunctive cuts:");
278  tighter.print ();
279  }
280  cs.insert (tighter);
281  }
282 
283  int deltaNcuts =
284  cs.sizeRowCuts () - initRowCuts +
285  cs.sizeColCuts () - initColCuts;
286 
287  if (info.level <= 0 && !(info.inTree))
288  nrootcuts_ += deltaNcuts;
289  ntotalcuts_ += deltaNcuts;
290 
291  if (jnlst_ -> ProduceOutput (J_DETAILED, J_DISJCUTS)) {
292 
293  if (cs.sizeRowCuts()>initRowCuts) printf ("added %d row cuts\n", cs.sizeRowCuts () - initRowCuts);
294  if (cs.sizeColCuts()>initColCuts) printf ("added %d col cuts\n", cs.sizeColCuts () - initColCuts);
295  }
296 
297  if ((info.level <= 0) && !(info.inTree))
298  jnlst_ -> Printf (J_ERROR, J_COUENNE,
299  "%d cuts\n", CoinMax (0, nrootcuts_ - nInitCuts));
300  else if (deltaNcuts)
301  jnlst_ -> Printf (J_WARNING, J_COUENNE,
302  "In-BB disjunctive cuts: %d row cuts, %d col cuts\n",
303  cs.sizeRowCuts () - initRowCuts,
304  cs.sizeColCuts () - initColCuts);
305  }
306 
307  if (jnlst_ -> ProduceOutput (J_DETAILED, J_DISJCUTS)) {
308 
309  printf ("Disjunctive cuts (%d row + %d col):\n", cs.sizeRowCuts () - initRowCuts, cs.sizeColCuts () - initColCuts);
310  for (int i=initRowCuts; i<cs.sizeRowCuts (); ++i) cs.rowCutPtr (i) -> print ();
311  for (int i=initColCuts; i<cs.sizeColCuts (); ++i) cs.colCutPtr (i) -> print ();
312  }
313 
314  // else {
315  // jnlst_ -> Printf (J_STRONGWARNING, J_COUENNE,
316  // "In-BB disjunctive cuts: %d row cuts, %d col cuts\n",
317  // cs.sizeRowCuts () - initRowCuts,
318  // cs.sizeColCuts () - initColCuts);
319  // if ((info.level <= 0) && !(info.inTree))
320  // jnlst_ -> Printf (J_ERROR, J_COUENNE,
321  // "infeasible node\n");
322  // }
323 
324  delete csi;
325 
326  septime_ += (CoinCpuTime () - time);
327 }
bool isWiped(OsiCuts &cs)
Check whether the previous cut generators have added an infeasible cut.
void fint fint fint real fint real real real real real real real real real fint real fint fint fint real fint fint fint fint * info
const Ipopt::EJournalCategory J_DISJCUTS(Ipopt::J_USER6)
#define COUENNE_EPS
static int
Definition: OSdtoa.cpp:2173
const Ipopt::EJournalCategory J_COUENNE(Ipopt::J_USER8)
void WipeMakeInfeas(OsiCuts &cs)
Add a fictitious cut 1&lt;= x_0 &lt;= -1 as a signal to the node solver that this node is deemed infeasible...