source: stable/0.3/Couenne/src/branch/CouenneVarObject.cpp @ 578

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

added fictitious objective when none is declared in the .nl. Fixing branch on large value with very negative lower and very large upper bounds

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