source: trunk/Couenne/src/branch/CouenneVarObject.cpp @ 591

Last change on this file since 591 was 591, checked in by pbelotti, 9 years ago

probing now applied even in absence of NLP solution

  • Property svn:keywords set to Author Date Id Revision
File size: 10.9 KB
Line 
1/* $Id: CouenneVarObject.cpp 591 2011-05-31 17:03:03Z pbelotti $
2 *
3 * Name:    CouenneVarObject.cpp
4 * Authors: Pietro Belotti, Carnegie Mellon University
5 * Purpose: Base object for variables (to be used in branching)
6 *
7 * (C) Carnegie-Mellon University, 2008-11.
8 * This file is licensed under the Eclipse Public License (EPL)
9 */
10
11#include "CoinHelperFunctions.hpp"
12
13#include "BonBabSetupBase.hpp"
14
15#include "CouenneProblem.hpp"
16#include "CouenneVarObject.hpp"
17#include "CouenneBranchingObject.hpp"
18#include "CouenneComplObject.hpp"
19
20using namespace Couenne;
21
22/// Constructor with information for branching point selection strategy
23CouenneVarObject::CouenneVarObject (CouenneCutGenerator *c,
24                                    CouenneProblem *p,
25                                    exprVar *ref, 
26                                    Bonmin::BabSetupBase *base, 
27                                    JnlstPtr jnlst,
28                                    int varSelection):
29
30  // Do not set variable (expression).
31  // This way, no expression-dependent strategy is chosen
32  CouenneObject (c, p, ref, base, jnlst),
33  varSelection_ (varSelection) {
34
35  if (jnlst_ -> ProduceOutput (J_SUMMARY, J_BRANCHING)) {
36    printf ("created Variable Object: "); 
37    reference_ -> print (); 
38    printf (" with %s strategy [clamp=%g, alpha=%g]\n", 
39            (strategy_ == LP_CLAMPED)   ? "lp-clamped" : 
40            (strategy_ == LP_CENTRAL)   ? "lp-central" : 
41            (strategy_ == BALANCED)     ? "balanced"   : 
42            (strategy_ == MIN_AREA)     ? "min-area"   : 
43            (strategy_ == MID_INTERVAL) ? "mid-point"  : 
44            (strategy_ == NO_BRANCH)    ? "no-branching (null infeasibility)" : 
45                                          "no strategy",
46            lp_clamp_, alpha_);
47  }
48}
49
50
51/// apply the branching rule
52OsiBranchingObject *CouenneVarObject::createBranch (OsiSolverInterface *si, 
53                                                    const OsiBranchingInformation *info, 
54                                                    int way) const {
55
56  // Before anything, use the LP point if so instructed
57
58  problem_ -> domain () -> push
59    (problem_ -> nVars (),
60     info -> solution_,
61     info -> lower_,
62     info -> upper_, false); // don't have to alloc+copy
63
64  // For some obscure reason, this only seems to work if
65  // strong/reliability branching is not used. I suppose it has to do
66  // with omitted pseudocost multiplier update, or some other hidden
67  // part of the code.
68
69  if ((varSelection_ == Bonmin::BabSetupBase::OSI_SIMPLE) && 
70      ((strategy_ == CouenneObject::LP_CLAMPED) ||
71       (strategy_ == CouenneObject::LP_CENTRAL) ||
72       (strategy_ == CouenneObject::MID_INTERVAL))) {
73
74    int indVar = reference_ -> Index ();
75
76    CouNumber
77      brpt  = info -> solution_ [indVar],
78      l     = info -> lower_    [indVar],
79      u     = info -> upper_    [indVar];
80     
81    // these vanilla assignments would drive Couenne crazy. If any of
82    // the bounds is (in absolute value) above 1e20, branching values
83    // can be detrimental for bound tightening and convexification
84    // cuts, for instance.
85
86    // Actually, even smaller values (such as 1e10) can trigger bad BB
87    // tree development (see wall in test files, not globallib/wall)
88
89#define LARGE_VALUE 1e8
90
91    if ((l < - LARGE_VALUE) && 
92        (u >   LARGE_VALUE) &&
93        (fabs (brpt) > LARGE_VALUE / 10))
94      brpt = 0.;
95
96    if (l < - COUENNE_INFINITY) l = -1. - 2 * fabs (brpt);
97    if (u >   COUENNE_INFINITY) u = +1. + 2 * fabs (brpt);
98
99    CouNumber width = lp_clamp_ * (u-l);
100
101    switch (strategy_) {
102
103    case CouenneObject::LP_CLAMPED:   brpt = CoinMax (l + width, CoinMin (brpt, u - width));        break;
104    case CouenneObject::LP_CENTRAL:   if ((brpt < l + width) || 
105                                          (brpt > u - width)) brpt = .5 * (l+u);                    break;
106    case CouenneObject::MID_INTERVAL: brpt = midInterval (brpt, 
107                                                          info -> lower_ [indVar], 
108                                                          info -> upper_ [indVar]);                 break;
109    default: assert (false); // this will never be used
110    }
111
112    OsiBranchingObject *obj = new CouenneBranchingObject (si, this, jnlst_, cutGen_, problem_, reference_, 
113                                                          TWO_LEFT, brpt, doFBBT_, doConvCuts_);
114
115    problem_ -> domain () -> pop ();
116
117    return obj;
118  }
119
120  // now deal with the more complicated branching selections
121
122  // The infeasibility on an (independent) variable x_i is given by
123  // something more elaborate than |w-f(x)|, that is, a function of
124  // all infeasibilities of all expressions which depend on x_i.
125
126  // problem_ -> domain () -> push
127  //   (problem_ -> nVars (),
128  //    info -> solution_,
129  //    info -> lower_,
130  //    info -> upper_, false); // don't have to alloc+copy
131
132  int bestWay;
133  const CouenneObject *criticalObject = NULL; // should create the branchingObject
134
135  CouNumber bestPt = computeBranchingPoint (info, bestWay, criticalObject);
136
137  ///////////////////////////////////////////
138
139  int indVar = reference_ -> Index ();
140
141  jnlst_ -> Printf (J_ITERSUMMARY, J_BRANCHING, ":::: creating branching on x_%d @%g [%g,%g]\n", 
142                    indVar, 
143                    info -> solution_ [indVar],
144                    info -> lower_    [indVar],
145                    info -> upper_    [indVar]);
146
147  OsiBranchingObject *brObj = criticalObject ? 
148    criticalObject -> createBranch (si, info, way) :
149    new CouenneBranchingObject (si, this, jnlst_, cutGen_, problem_, reference_, 
150                                way, bestPt, doFBBT_, doConvCuts_);
151
152  problem_ -> domain () -> pop ();
153
154  return brObj;
155}
156
157
158/// compute branching point (used in createBranch ())
159CouNumber CouenneVarObject::computeBranchingPoint(const OsiBranchingInformation *info,
160                                                  int& bestWay,
161                                                  const CouenneObject *&criticalObject) const {
162  criticalObject = NULL;
163
164  if (jnlst_ -> ProduceOutput (J_DETAILED, J_BRANCHING)) {
165    printf ( "---------- computeBRPT for "); 
166    reference_ -> print (); 
167    printf (" [%g,%g]\n", 
168            info -> lower_ [reference_ -> Index ()],
169            info -> upper_ [reference_ -> Index ()]);
170  }
171
172  expression *brVar = NULL; // branching variable
173
174  CouNumber
175    brdistDn = 0.,
176    brdistUp = 0.,
177    bestPt   = 0.,
178   *brPts    = NULL, // branching point(s)
179   *brDist   = NULL, // distances from LP point to each new convexification
180    maxdist = - COIN_DBL_MAX;
181
182  bool chosen = false;
183
184  bestWay = TWO_LEFT;
185  int whichWay = TWO_LEFT,
186    index = reference_ -> Index ();
187
188  std::set <int> deplist = problem_ -> Dependence () [index];
189
190  for (std::set <int>::iterator i = deplist.begin (); i != deplist.end (); ++i) {
191
192    const CouenneObject *obj = problem_ -> Objects () [*i];
193
194    CouNumber improv = 0.;
195
196    assert (obj -> Reference ());
197
198    if (jnlst_ -> ProduceOutput (J_MATRIX, J_BRANCHING)) {
199      printf ("  ** "); 
200      obj -> Reference () -> print (); 
201      if (reference_ -> Image ()) {printf (" := "); obj -> Reference () -> Image () -> print ();}
202      printf ("\n");
203    }
204
205    if (obj -> Reference ()) {
206      if (obj -> Reference () -> Image ())
207        improv = obj -> Reference () -> Image ()
208          -> selectBranch (obj, info,                      // input parameters
209                           brVar, brPts, brDist, whichWay); // result: who, where, distances, direction
210      else {
211        brVar = obj -> Reference ();
212        brPts  = (double *) realloc (brPts, sizeof (double)); 
213        brDist = (double *) realloc (brDist, 2 * sizeof (double)); 
214
215        double point = info -> solution_ [obj -> Reference () -> Index ()];
216
217        *brPts = point;
218        improv = 0.;
219
220        if (point > floor (point)) {improv =                  brDist [0] = point - floor (point);}
221        if (point < ceil  (point)) {improv = CoinMin (improv, brDist [1] = ceil (point) - point);}
222
223        point -= floor (point);
224        whichWay = (point < 0.45) ? TWO_LEFT : (point > 0.55) ? TWO_RIGHT : TWO_RAND;
225      }
226    }
227
228    if (jnlst_ -> ProduceOutput (J_MATRIX, J_BRANCHING)) {
229      printf ("  --> "); 
230      if (brVar) {
231        brVar -> print (); 
232        if (brPts) 
233          printf (" at %g, improv %g <%g>, indices = %d,%d\n", 
234                  *brPts, improv, maxdist, index, brVar -> Index ());
235      }
236    }
237
238    if (brVar &&
239        (brVar -> Index () == index) &&    // it's us!
240        (fabs (improv) > maxdist) &&       // this branching seems to induce a higher improvement
241        (fabs (*brPts) < COU_MAX_COEFF)) { // and branching pt is limited
242
243      criticalObject = (problem_ -> Objects () [*i]); // set this object as the branch creator
244
245      brdistDn = brDist [0];
246      brdistUp = brDist [1];
247      chosen = true;
248      bestPt = *brPts;
249      maxdist = improv;
250      bestWay = whichWay;
251    }
252  }
253
254  // no hits on this VarObject's variable, that is, this variable was
255  // never chosen
256
257  if (!chosen) {
258
259    bestPt = info -> solution_ [index];
260    brVar  = reference_;
261
262    CouNumber
263      l     = info -> lower_ [index], 
264      u     = info -> upper_ [index],
265      width = lp_clamp_ * (u-l);
266
267    switch (strategy_) {
268
269    case CouenneObject::LP_CLAMPED:
270      bestPt = CoinMax (l + width, CoinMin (bestPt, u - width));
271      break;
272    case CouenneObject::LP_CENTRAL: 
273      if ((bestPt < l + width) || (bestPt > u - width))
274        bestPt = (l+u)/2;
275      break;
276    case CouenneObject::MID_INTERVAL: 
277    default: 
278      // all other cases (balanced, min-area)
279      bestPt = midInterval (bestPt, l, u);
280
281      if (jnlst_ -> ProduceOutput (J_MATRIX, J_BRANCHING)) {
282        if (CoinMin (fabs (bestPt - l), fabs (bestPt - u)) < 1e-3) {
283          printf ("computed failsafe %g [%g,%g] for ", bestPt, l,u);
284          reference_ -> print (); printf ("\n");
285        }
286      }
287      break;
288    }
289
290#define LARGE_VALUE 1e8
291    if ((l < - LARGE_VALUE) && 
292        (u >   LARGE_VALUE) &&
293        (fabs (bestPt) > LARGE_VALUE / 10))
294      bestPt = 0.;
295
296    brPts  = (double *) realloc (brPts, sizeof (double));
297    *brPts = bestPt;
298
299    if (jnlst_ -> ProduceOutput (J_MATRIX, J_BRANCHING)) {
300      printf ("  ::: failsafe:  %g [%g,%g] for ", 
301              bestPt, info -> lower_ [index], info -> upper_ [index]); 
302      reference_ -> print ();
303      printf ("\n");
304    }
305
306  } else {
307
308      if (jnlst_ -> ProduceOutput (J_MATRIX, J_BRANCHING)) {
309        if (CoinMin (fabs (bestPt - info -> lower_ [index]), 
310                     fabs (bestPt - info -> upper_ [index])) < 1e-3) {
311          printf ("  computed %g [%g,%g] for ", 
312                  bestPt, info -> lower_ [index], info -> upper_ [index]); 
313          reference_ -> print ();
314          printf ("\n");
315        }
316      }
317  }
318
319  if (pseudoMultType_ == PROJECTDIST) {
320
321    if (chosen) {
322      downEstimate_ = brdistDn;
323      upEstimate_   = brdistUp;
324    } 
325    else downEstimate_ = upEstimate_ = 1.;
326  }
327
328  if (brPts)  free (brPts);
329  if (brDist) free (brDist);
330
331  return bestPt;
332}
333
334
335#define TOL 0.
336
337/// fix nonlinear coordinates of current integer-nonlinear feasible solution
338double CouenneVarObject::feasibleRegion (OsiSolverInterface *solver, 
339                                         const OsiBranchingInformation *info) const {
340  int index = reference_ -> Index ();
341
342  assert (index >= 0);
343
344  double val = info -> solution_ [index];
345
346  // fix that variable to its current value
347  solver -> setColLower (index, val-TOL);
348  solver -> setColUpper (index, val+TOL);
349
350  return 0.;
351}
352
353
354/// are we on the bad or good side of the expression?
355bool CouenneVarObject::isCuttable () const {
356
357  const std::set <int>                &deplist = problem_ -> Dependence () [reference_ -> Index ()];
358  const std::vector <CouenneObject *> &objects = problem_ -> Objects ();
359
360  for (std::set <int>::const_iterator depvar = deplist. begin ();
361       depvar != deplist. end (); ++depvar)
362    if (!(objects [*depvar] -> isCuttable ()))
363      return false;
364
365  return (!(reference_ -> isInteger ()));
366}
367
Note: See TracBrowser for help on using the repository browser.