00001 /* $Id: computeMulBrDist.cpp 141 2009-06-03 04:19:19Z pbelotti $ */ 00002 /* 00003 * Name: computeMulBrDist.cpp 00004 * Author: Pietro Belotti 00005 * Purpose: compute distance to new convexifications generated by branching on product 00006 * 00007 * (C) Carnegie-Mellon University, 2008. 00008 * This file is licensed under the Common Public License (CPL) 00009 */ 00010 00011 #include "CouennePrecisions.hpp" 00012 #include "CouenneTypes.hpp" 00013 #include "CouenneObject.hpp" 00014 00015 #include "exprMul.hpp" 00016 #include "funtriplets.hpp" 00017 #include "projections.hpp" 00018 00019 // compute distance from future convexifications in set \f$\{(x,y,w): 00020 // w = xy\}\f$ with x,y,w bounded 00021 double *computeMulBrDist (const OsiBranchingInformation *info, 00022 int xi, int yi, int wi, 00023 int brind, double *brpt, int nPts) { 00024 00025 // use rule of thumb to compute distance: fix two of the three 00026 // variables and compute distance between current LP point and curve 00027 // w=x/y (z is the branching variable, x or y) 00028 // 00029 // 1) fix x,y: distances are w - x/y and ||(w-x/y, z-zb)||_2 00030 // 2) fix x,w: y - x/w and ||(y-x/w, z-zb)||_2 00031 // 3) fix y,w: x - y*w and ||(x-y*w, z-zb)||_2 00032 00033 CouNumber 00034 x0 = info -> solution_ [xi], //xl = info -> lower_ [xi], xu = info -> upper_ [xi], 00035 y0 = info -> solution_ [yi], //yl = info -> lower_ [yi], yu = info -> upper_ [yi], 00036 w0 = info -> solution_ [wi]; //wl = info -> lower_ [wi], wu = info -> upper_ [wi]; 00037 00038 double *dist = (double *) malloc (2 * sizeof (double)); 00039 00040 // two main cases: 00041 // 00042 // 1) wi is the branching index: the bounding box is divided in two 00043 // by the rule w <= wb and w >= wb. Finding the distances from the 00044 // current point (x0,y0,w0) to the two semi-surfaces depends on 00045 // which side w0 stands. 00046 // 00047 // 2) xi or yi are the branching index: reduce to the problem of 00048 // finding the distance from a point (x0,y0,w0) to the same surface 00049 // within the (two) element of a partition of the bounding box. 00050 00051 00052 // case 1 //////////////////////////////////////////////////////////////////// 00053 // 00054 // Depending on whether 00055 // 00056 // a) w0 <=/>= wb 00057 // b) w0 <=/>= x0*y0 00058 // c) wb <=/>= 0 00059 // 00060 // there are eight (!) cases, each with a similar computation for 00061 // the two distances. We unify it below. 00062 00063 // [for now give a simplified version with rough distance computation] 00064 00065 if (brind == wi) { 00066 00067 // easy implementation: for own side, w0 - x0 y0; for other side, 00068 // horiz/vert distance to curve x0 w0 = wb 00069 00070 bool side = ((x0*y0 > *brpt) && (*brpt > 0) || 00071 (x0*y0 < *brpt) && (*brpt < 0)); 00072 00073 dist [side ? 1 : 0] = CoinMax 00074 (fabs (w0 - x0*y0), CoinMin 00075 ((fabs (y0) > COUENNE_EPS) ? fabs (x0 - *brpt / y0) : 0., 00076 (fabs (x0) > COUENNE_EPS) ? fabs (y0 - *brpt / x0) : 0.)); 00077 00078 dist [side ? 0 : 1] = fabs (w0 - x0*y0); 00079 } 00080 00081 // case 2 //////////////////////////////////////////////////////////////////// 00082 00083 else { 00084 00085 CouNumber diff = info -> solution_ [brind] - *brpt; 00086 bool side = (diff > 0.); 00087 00088 dist [side ? 0 : 1] = CoinMax (fabs (w0 - x0 * y0), fabs (diff)); 00089 dist [side ? 1 : 0] = fabs (w0 - x0 *y0); 00090 } 00091 00092 //dist [0] = dist [1] = 1; 00093 return dist; 00094 }