/home/coin/SVN-release/OS-2.0.1/Couenne/src/convex/createCuts.cpp

Go to the documentation of this file.
00001 /* $Id: createCuts.cpp 141 2009-06-03 04:19:19Z pbelotti $ */
00002 /*
00003  * Name:    createCuts.cpp
00004  * Author:  Pietro Belotti
00005  * Purpose: a standard cut creator for use with convexification
00006  *
00007  * (C) Carnegie-Mellon University, 2006-08. 
00008  * This file is licensed under the Common Public License (CPL)
00009  */
00010 
00011 #include "OsiRowCut.hpp"
00012 
00013 #include "CouenneTypes.hpp"
00014 #include "CouennePrecisions.hpp"
00015 #include "CouenneCutGenerator.hpp"
00016 #include "CouenneProblem.hpp"
00017 
00018 
00020 bool badCoeff (CouNumber coe) {
00021 
00022   coe = fabs (coe);
00023   return ((coe > COU_MAX_COEFF) || ((coe < COU_MIN_COEFF) && (coe > 0.)));
00024 }
00025 
00026 
00029 
00030 int CouenneCutGenerator::createCut (OsiCuts &cs,
00031                                     CouNumber lb, CouNumber ub, 
00032                                     int i1, CouNumber c1,
00033                                     int i2, CouNumber c2,
00034                                     int i3, CouNumber c3,
00035                                     bool is_global) const {
00036   bool numerics = false;
00037 
00038   // a maximum of three terms allowed here. If index -1 (default) the
00039   // term is not considered
00040 
00041   int nterms = 0;
00042 
00043   // CAUTION: this can make the problem infeasible...
00044   if (fabs (c3) <= 1.e-21) {                                    i3 = -1;} // shift coeff/index to
00045   if (fabs (c2) <= 1.e-21) {                  c2 = c3; i2 = i3; i3 = -1;} // keep consistency
00046   if (fabs (c1) <= 1.e-21) {c1 = c2; i1 = i2; c2 = c3; i2 = i3; i3 = -1;}
00047   // why 1.0e-21? Look at CoinPackedMatrix.cpp:2188
00048 
00049 #if 0
00050   if (i1 >= 0) {if (fabs (c1) > COU_MAX_COEFF) numerics = true; nterms++;} else c1 = 0;
00051   if (i2 >= 0) {if (fabs (c2) > COU_MAX_COEFF) numerics = true; nterms++;} else c2 = 0;
00052   if (i3 >= 0) {if (fabs (c3) > COU_MAX_COEFF) numerics = true; nterms++;} else c3 = 0;
00053 #else
00054   if (i1 >= 0){if (badCoeff (c1)) numerics = true; nterms++;} else c1 = 0;
00055   if (i2 >= 0){if (badCoeff (c2)) numerics = true; nterms++;} else c2 = 0;
00056   if (i3 >= 0){if (badCoeff (c3)) numerics = true; nterms++;} else c3 = 0;
00057 #endif
00058 
00059   if (!nterms) // nonsense cut
00060     return 0;
00061 
00062   // cut has large coefficients/rhs, bail out
00063   if (numerics
00064       //|| ((fabs (lb) < COU_MIN_COEFF) || 
00065       //(fabs (ub) < COU_MIN_COEFF))
00066       || ((fabs (lb) > COU_MAX_COEFF) &&
00067           (fabs (ub) > COU_MAX_COEFF))) {
00068 
00069     jnlst_->Printf(J_STRONGWARNING, J_CONVEXIFYING,
00070                    "### Discarding cut, large coeff/rhs: %g (%d), %g (%d), %g (%d); [%g,%g]\n", 
00071                    c1, i1, c2, i2, c3, i3, lb, ub);
00072     return 0;
00073   }
00074 
00075   if (!firstcall_ && addviolated_) { // need to check violation 
00076 
00077     const CouNumber *x = problem_ -> X ();
00078 
00079     // compute violation
00080     CouNumber violation = 0;
00081 
00082     if (i1 >= 0) violation += c1 * x [i1];
00083     if (i2 >= 0) violation += c2 * x [i2];
00084     if (i3 >= 0) violation += c3 * x [i3];
00085 
00086     // return 0 if not violated
00087 
00088     if ((violation < ub + 0 * COUENNE_EPS) &&
00089         (violation > lb - 0 * COUENNE_EPS))
00090       return 0;
00091   }
00092 
00093   // check if cut violates optimal solution (irrespective of the
00094   // branching rules applied, so handle with care)
00095 
00096   CouNumber *best = problem_ -> bestSol ();
00097 
00098   if (best &&
00099       ((i1 < 0) || ((best [i1] >= problem_ -> Lb (i1)) && (best [i1] <= problem_ -> Ub (i1)))) &&
00100       ((i2 < 0) || ((best [i2] >= problem_ -> Lb (i2)) && (best [i2] <= problem_ -> Ub (i2)))) &&
00101       ((i3 < 0) || ((best [i3] >= problem_ -> Lb (i3)) && (best [i3] <= problem_ -> Ub (i3))))) {
00102 
00103     CouNumber lhs = 0.;
00104 
00105     if (i1 >= 0) lhs += c1 * best [i1];
00106     if (i2 >= 0) lhs += c2 * best [i2];
00107     if (i3 >= 0) lhs += c3 * best [i3];
00108 
00109     if (lhs > ub + COUENNE_EPS)
00110       jnlst_->Printf(J_STRONGWARNING, J_CONVEXIFYING,
00111                      "### cut (%d,%d,%d) (%g,%g,%g) violates optimum: %g >= %g [%g]\n", 
00112                      i1,i2,i3, c1,c2,c3, lhs, ub, lhs - ub); 
00113 
00114     if (lhs < lb - COUENNE_EPS)
00115       jnlst_->Printf(J_STRONGWARNING, J_CONVEXIFYING,
00116                      "### cut (%d,%d,%d) (%g,%g,%g) violates optimum: %g <= %g [%g]\n", 
00117                      i1,i2,i3, c1,c2,c3, lhs, lb, lb - lhs);
00118   }
00119 
00120   // You are here if:
00121   //
00122   // 1) this is the first call to CouenneCutGenerator::generateCuts(), or
00123   // 2) you also want unviolated cuts, or
00124   // 3) the cut is violated
00125 
00126   // two cases: cut is of the form w1 [<|>]= alpha, hence a column
00127   // cut, or it is of the form (a w1 + b w2 + c w3 [<|>]= alpha), a
00128   // row cut
00129 
00130   if ((i2 < 0) && (i3 < 0)) { // column cut /////////////////////////////////////////
00131 
00132     if (   (fabs (c1) < COUENNE_EPS)
00133         && (fabs (lb) > COU_MAX_COEFF * COUENNE_EPS)
00134         && (fabs (ub) > COU_MAX_COEFF * COUENNE_EPS)) {
00135 
00136       jnlst_->Printf(J_STRONGWARNING, J_CONVEXIFYING,
00137                      "#### nonsense column cut: %e <= %e w_%d <= %e\n", 
00138                      lb, c1, i1, ub);
00139       return 0;
00140     }
00141 
00142     OsiColCut *cut = new OsiColCut;
00143 
00144     CouNumber 
00145       ll = lb / c1, 
00146       uu = ub / c1;
00147 
00148     if (c1 < 0) {
00149       CouNumber tmp = ll;
00150       ll = uu;
00151       uu = tmp;
00152     }
00153 
00154     CouNumber &curL = problem_ -> Lb (i1),
00155               &curU = problem_ -> Ub (i1);
00156 
00157     if (uu < COUENNE_INFINITY) {
00158       cut -> setUbs (1, &i1, &uu); 
00159       if (uu < curU - COUENNE_EPS) {
00160         curU = uu; // TODO: chg_bds
00161       }
00162     }
00163 
00164     if (ll > -COUENNE_INFINITY) {
00165       cut -> setLbs (1, &i1, &ll); 
00166       if (ll > curL + COUENNE_EPS) {
00167         curL = ll; // idem
00168       }
00169     }
00170 
00171     cut -> setGloballyValid (is_global); // global?
00172 
00173     cs.insert (cut);
00174     delete cut;
00175 
00176   } else {
00177 
00178     // row cut //////////////////////////////////////////////////////////////////////
00179 
00180     CouNumber *coeff = new CouNumber [nterms]; 
00181     int       *index = new int       [nterms];
00182     OsiRowCut *cut   = new OsiRowCut;
00183 
00184     int nt = 0;
00185 
00186     if (i1 >= 0) {coeff [nt] = c1; index [nt++] = i1;}
00187     if (i2 >= 0) {coeff [nt] = c2; index [nt++] = i2;}
00188     if (i3 >= 0) {coeff [nt] = c3; index [nt++] = i3;}
00189 
00190     if (lb > -COUENNE_INFINITY) cut -> setLb (lb);
00191     if (ub <  COUENNE_INFINITY) cut -> setUb (ub);
00192 
00193     cut -> setRow (nterms, index, coeff);
00194 
00195     delete [] coeff;
00196     delete [] index;
00197 
00198     cut -> setGloballyValid (is_global); // global?
00199 
00200     cs.insert (cut);
00201     delete cut;
00202   }
00203 
00204   return 1;
00205 }
00206 
00207 
00210 
00211 int CouenneCutGenerator::createCut (OsiCuts &cs,
00212                                     CouNumber rhs, int sign, 
00213                                     int i1, CouNumber c1,
00214                                     int i2, CouNumber c2,
00215                                     int i3, CouNumber c3,
00216                                     bool is_global)       const {
00217 
00218   return createCut (cs, (sign >= 0) ? rhs : - COIN_DBL_MAX,
00219                         (sign <= 0) ? rhs :   COIN_DBL_MAX,
00220                     i1, c1, i2, c2, i3, c3, is_global);
00221 }

Generated on Thu Oct 8 03:02:56 2009 by  doxygen 1.4.7