/home/coin/SVN-release/OS-2.4.1/Couenne/src/branch/operators/branchExprDiv.cpp

Go to the documentation of this file.
00001 /* $Id: branchExprDiv.cpp 708 2011-06-23 14:04:59Z pbelotti $
00002  *
00003  * Name:    branchExprDiv.cpp
00004  * Author:  Pietro Belotti
00005  * Purpose: select branch for divisions
00006  *
00007  * (C) Carnegie-Mellon University, 2006-11.
00008  * This file is licensed under the Eclipse Public License (EPL)
00009  */
00010 
00011 #include "CouenneExprDiv.hpp"
00012 #include "CouenneExprMul.hpp"
00013 #include "CouenneBranchingObject.hpp"
00014 #include "CouenneObject.hpp"
00015 
00016 using namespace Couenne;
00017 
00020 CouNumber exprDiv::selectBranch (const CouenneObject *obj, 
00021                                  const OsiBranchingInformation *info,
00022                                  expression *&var,
00023                                  double * &brpts, 
00024                                  double * &brDist, // distance of current LP
00025                                                    // point to new convexifications
00026                                  int &way) {
00027 
00028   if (brDist) {free (brDist); brDist = NULL;} // clear it, computeMulBrDist will fill it
00029 
00030   int xi = arglist_ [0]        -> Index (),
00031       yi = arglist_ [1]        -> Index (),
00032       wi = obj -> Reference () -> Index ();
00033 
00034   assert ((xi >= 0) && (yi >= 0) && (wi >= 0));
00035 
00036   brpts  = (double *) realloc (brpts,  sizeof (double));
00037 
00038   // choosing branching variable and point is difficult, use
00039   // proportion in bound intervals
00040 
00041   CouNumber yl = info -> lower_    [yi], 
00042             yu = info -> upper_    [yi],
00043             y0 = info -> solution_ [yi];
00044 
00045   // if [yl,yu] contains 0, create two nodes
00046 
00047   if ((yl < -COUENNE_EPS) && 
00048       (yu >  COUENNE_EPS)) {
00049 
00050     var = arglist_ [1];
00051 
00052     *brpts = 0.;
00053 
00054     way = (y0 > *brpts) ? TWO_RIGHT : TWO_LEFT;
00055 
00056     brDist = computeMulBrDist (info, wi, yi, xi, yi, brpts);
00057 
00058     return CoinMin (brDist [0], brDist [1]);
00059     /*return ((fabs (y0) < COUENNE_EPS) ? 1. :
00060             fabs (info -> solution_ [xi] / y0 - 
00061             info -> solution_ [wi]));*/
00062   }
00063 
00064   // From now on, [yl,yu] may be unlimited in one sense only, and
00065   // interval does not contain 0.
00066   //
00067   // As convexification is still carried out by applying McCormick
00068   // rules to x=w*y (where original operator is w=x/y), try to get
00069   // closer to a situation where both y and w are bounded, if
00070   // necessary by branching on w.
00071   //
00072   // Branch first on y if unbounded, and then on w. As a result of
00073   // bound tightening, if both y and w are bounded, x will be, too.
00074 
00075   if ((yl < -COUENNE_INFINITY) || // and yu < 0
00076       (yu >  COUENNE_INFINITY)) { // and yl > 0
00077 
00078     var = arglist_ [1];
00079 
00080     // if y0 close to bounds, branch away from it
00081     if      (fabs (y0-yl) < COUENNE_NEAR_BOUND) *brpts = y0 + 1. + yl*10.;
00082     else if (fabs (y0-yu) < COUENNE_NEAR_BOUND) *brpts = y0 - 1. + yu*10.;
00083     else                                        *brpts = y0;
00084 
00085     way = (y0 > 0.) ? TWO_LEFT : TWO_RIGHT;
00086 
00087     brDist = computeMulBrDist (info, wi, yi, xi, yi, brpts);
00088 
00089     return CoinMin (brDist [0], brDist [1]);
00090 
00091     //return ((fabs (y0) < COUENNE_EPS) ? 1. : 
00092     //fabs (info -> solution_ [xi] / y0 - info -> solution_ [wi]));
00093   }
00094 
00095   // y is bounded, and y0 nonzero; if w is unbounded, bound w first as
00096   // x will be too.
00097 
00098   CouNumber wl = info -> lower_    [wi], 
00099             wu = info -> upper_    [wi],
00100             w0 = info -> solution_ [wi],
00101             x0 = info -> solution_ [xi];
00102 
00103   // w unbounded in >= one direction ///////////////////////////////////////
00104 
00105   if ((wl < -COUENNE_INFINITY) || 
00106       (wu >  COUENNE_INFINITY)) {
00107 
00108     var = obj -> Reference ();
00109 
00110     if ((wl < -COUENNE_INFINITY) &&
00111         (wu >  COUENNE_INFINITY)) { // unbounded in two directions
00112 
00113       CouNumber 
00114         wreal = x0 / y0, 
00115         wmin  = w0, 
00116         wmax  = wreal; // assume (x0,y0,w0) is below w=x/y
00117 
00118       if (wreal < w0) { // but swap if (x0,y0,w0) is above w=x/y
00119         wmin = wreal;
00120         wmax = w0;
00121       }
00122 
00123       *brpts = wreal;
00124       way = (w0 < wreal) ? TWO_LEFT : TWO_RIGHT;
00125 
00126       brDist = computeMulBrDist (info, wi, yi, xi, wi, brpts);
00127       return CoinMin (brDist [0], brDist [1]);
00128 
00129     } else {
00130 
00131       // unbounded in one direction only, use two way branching
00132 
00133       // if y0 close to bounds, branch away from it
00134       if      (fabs (w0-wl) < COUENNE_NEAR_BOUND) *brpts = w0 + 1 + wl*10;
00135       else if (fabs (w0-wu) < COUENNE_NEAR_BOUND) *brpts = w0 - 1 + wu*10;
00136       else                                        *brpts = w0;
00137 
00138       way = (wl < - COUENNE_INFINITY) ? TWO_RIGHT : TWO_LEFT;
00139 
00140       brDist = computeMulBrDist (info, wi, yi, xi, wi, brpts);
00141       return CoinMin (brDist [0], brDist [1]);
00142     }
00143     //return ((fabs (y0) < COUENNE_EPS) ? 1. : fabs (x0/y0 - w0));
00144   }
00145 
00146   // w and y are bounded (and so is x). Choose between x, y, z
00147   // depending on intervals first and then to vicinity to midpoint
00148   CouNumber 
00149     xl = info -> lower_ [xi],
00150     xu = info -> upper_ [xi],
00151     dx = xu-xl,
00152     dy = yu-yl,
00153     dw = wu-wl;
00154 
00155   // Check largest interval and choose branch variable accordingly.
00156   // Branching point depends on where the current point is, but for
00157   // now just focus on the width of the intervals
00158 
00159   way = TWO_RAND;
00160 
00161   if (dx > dy)
00162     if (dx > dw) {var = arglist_[0];      *brpts = (xl+xu)/2.; /*return fabs (x0-y0*w0);*/}
00163     else         {var = obj->Reference(); *brpts = (wl+wu)/2.; /*return fabs (w0-x0/y0);*/}
00164   else
00165     if (dy > dw) {var = arglist_[1];      *brpts = (yl+yu)/2.; /*return fabs (y0-x0/w0);*/}
00166     else         {var = obj->Reference(); *brpts = (wl+wu)/2.; /*return fabs (w0-x0/y0);*/}
00167 
00168   brDist = computeMulBrDist (info, wi, yi, xi, var -> Index (), brpts);
00169   return CoinMin (brDist [0], brDist [1]);
00170 }

Generated on Thu Nov 10 03:05:43 2011 by  doxygen 1.4.7