00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include "CouenneCutGenerator.hpp"
00012 #include "CouenneProblem.hpp"
00013 #include "CouenneExprVar.hpp"
00014 #include "CouenneProblemElem.hpp"
00015 #include "BonBabInfos.hpp"
00016 #include "BonCbc.hpp"
00017
00018
00019 #define THRES_IMPROVED 0
00020
00021 using namespace Couenne;
00022
00023
00024
00025 bool CouenneProblem::btCore (t_chg_bounds *chg_bds) const {
00026
00027 if (!chg_bds) {
00028
00029 chg_bds = new t_chg_bounds [nVars ()];
00030
00031 for (int i=0; i < nVars (); i++)
00032
00033 if (Var (i) -> Multiplicity () > 0) {
00034
00035 chg_bds [i].setLower (t_chg_bounds::CHANGED);
00036 chg_bds [i].setUpper (t_chg_bounds::CHANGED);
00037 }
00038 }
00039
00040
00041
00042 int ntightened = 0,
00043 nbwtightened = 0,
00044 niter = 0;
00045
00046 bool first = true;
00047
00048 installCutOff ();
00049
00050
00051 bool contains_optimum = false;
00052
00053 if (optimum_ != NULL) {
00054 contains_optimum = true;
00055 for (int i=0; i < nVars (); i++)
00056 if ((optimum_ [i] < Lb (i) * (1 - COUENNE_EPS) - COUENNE_EPS) ||
00057 (optimum_ [i] > Ub (i) * (1 + COUENNE_EPS) + COUENNE_EPS)) {
00058
00059
00060
00061
00062 contains_optimum = false;
00063 break;
00064 }
00065 }
00066
00067 if (max_fbbt_iter_) do {
00068
00069 if (CoinCpuTime () > maxCpuTime_)
00070 break;
00071
00072
00073
00074
00075 ntightened = ((nbwtightened > 0) || first) ?
00076 tightenBounds (chg_bds) : 0;
00077
00078
00079
00080
00081
00082 nbwtightened = ((ntightened > 0) || ((ntightened==0) && first)) ? impliedBounds (chg_bds) : 0;
00083
00084 if (first)
00085 first = false;
00086
00087 if ((ntightened < 0) || (nbwtightened < 0)) {
00088 Jnlst () -> Printf (Ipopt::J_ITERSUMMARY, J_BOUNDTIGHTENING, "infeasible BT\n");
00089 return false;
00090 }
00091
00092
00093
00094
00095 if (contains_optimum) {
00096 for (int i=0; i<nVars (); i++)
00097 if ((optimum_[i] < Lb(i) - COUENNE_EPS * (1. + CoinMin (fabs(optimum_ [i]), fabs (Lb(i))))) ||
00098 (optimum_[i] > Ub(i) + COUENNE_EPS * (1. + CoinMin (fabs(optimum_ [i]), fabs (Ub(i)))))) {
00099 printf ("bound tightening CUTS optimum: x%d [%e,%e] val = %e, violation = %e\n",
00100 i, Lb (i), Ub (i), optimum_ [i],
00101 CoinMax (- optimum_ [i] + Lb (i),
00102 optimum_ [i] - Ub (i)));
00103 contains_optimum = false;
00104 }
00105 }
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124 } while (((ntightened > 0) || (nbwtightened > 0)) &&
00125 (ntightened + nbwtightened > THRES_IMPROVED) &&
00126 ((max_fbbt_iter_ < 0) || (niter++ < max_fbbt_iter_)));
00127
00128 fbbtReachedIterLimit_ = ((max_fbbt_iter_ > 0) && (niter >= max_fbbt_iter_));
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141 for (int i = 0, j = nVars (); j--; i++)
00142 if (Var (i) -> Multiplicity () > 0) {
00143
00144
00145 if ((Lb (i) > Ub (i) + COUENNE_BOUND_PREC * (1 + CoinMin (fabs (Lb (i)), fabs (Ub (i))))) ||
00146 (Ub (i) < - MAX_BOUND) ||
00147 (Lb (i) > MAX_BOUND)) {
00148
00149 Jnlst () -> Printf (Ipopt::J_ITERSUMMARY, J_BOUNDTIGHTENING, "final test: infeasible BT\n");
00150 return false;
00151 }
00152
00153
00154 if (Lb (i) > Ub (i)) {
00155 CouNumber swap = Lb (i);
00156 Lb (i) = Ub (i);
00157 Ub (i) = swap;
00158 }
00159 }
00160
00161 return true;
00162 }
00163
00164
00167
00168 bool CouenneProblem::boundTightening (t_chg_bounds *chg_bds,
00169 Bonmin::BabInfo * babInfo) const {
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179 Jnlst () -> Printf (Ipopt::J_ITERSUMMARY, J_BOUNDTIGHTENING,
00180 "Feasibility-based Bound Tightening\n");
00181
00182 int objInd = Obj (0) -> Body () -> Index ();
00183
00185
00186 if ((objInd >= 0) && babInfo && (babInfo -> babPtr ())) {
00187
00188 CouNumber UB = babInfo -> babPtr () -> model (). getObjValue(),
00189 LB = babInfo -> babPtr () -> model (). getBestPossibleObjValue (),
00190 primal0 = Ub (objInd),
00191 dual0 = Lb (objInd);
00192
00193
00194
00195
00196 if ((UB < COUENNE_INFINITY) &&
00197 (UB < primal0 - COUENNE_EPS)) {
00198
00199 Ub (objInd) = UB;
00200 chg_bds [objInd].setUpper(t_chg_bounds::CHANGED);
00201 }
00202
00203 if ((LB > - COUENNE_INFINITY) &&
00204 (LB > dual0 + COUENNE_EPS)) {
00205 Lb (objInd) = LB;
00206 chg_bds [objInd].setLower(t_chg_bounds::CHANGED);
00207 }
00208 }
00209
00210 return btCore (chg_bds);
00211
00212
00213
00214
00215 }
00216
00217
00219 int CouenneProblem::redCostBT (const OsiSolverInterface *psi,
00220 t_chg_bounds *chg_bds) const {
00221
00222 int
00223 nchanges = 0,
00224 objind = Obj (0) -> Body () -> Index ();
00225
00226 if (objind < 0)
00227 return 0;
00228
00229 CouNumber
00230 UB = getCutOff (),
00231 LB = Lb (objind);
00232
00234
00235 if ((LB > -COUENNE_INFINITY) &&
00236 (UB < COUENNE_INFINITY)) {
00237
00238 const double
00239 *X = psi -> getColSolution (),
00240 *L = psi -> getColLower (),
00241 *U = psi -> getColUpper (),
00242 *RC = psi -> getReducedCost ();
00243
00244 if (jnlst_ -> ProduceOutput (Ipopt::J_MATRIX, J_BOUNDTIGHTENING)) {
00245 printf ("REDUCED COST BT (LB=%g, UB=%g):\n", LB, UB);
00246 for (int i=0, j=0; i < nVars (); i++) {
00247 if (Var (i) -> Multiplicity () <= 0)
00248 continue;
00249 printf ("%3d %7e [%7e %7e] c %7e ", i, X [i], L [i], U [i], RC [i]);
00250 if (!(++j % 3))
00251 printf ("\n");
00252 }
00253 printf ("-----------\n");
00254 }
00255
00256 int ncols = psi -> getNumCols ();
00257
00258 for (int i=0; i<ncols; i++)
00259 if ((i != objind) &&
00260 (Var (i) -> Multiplicity () > 0)) {
00261
00262 CouNumber
00263 x = X [i],
00264 l = L [i],
00265 u = U [i],
00266 rc = RC [i];
00267
00268 #define CLOSE_TO_BOUNDS 1.e-15
00269
00270 if ((fabs (rc) < CLOSE_TO_BOUNDS) ||
00271 (fabs (l-u) < CLOSE_TO_BOUNDS))
00272 continue;
00273
00274 bool isInt = Var (i) -> isInteger ();
00275
00276 if ((rc >= 0.) &&
00277 (fabs (x-l) <= CLOSE_TO_BOUNDS)) {
00278
00279 if (LB + (u-l)*rc > UB) {
00280
00281 CouNumber newUb = l + (UB-LB) / rc;
00282 newUb = !isInt ? newUb : floor (newUb + COUENNE_EPS);
00283
00284 if (newUb < Ub (i)) {
00285
00286 Ub (i) = newUb;
00287
00288 nchanges++;
00289 chg_bds [i].setLower (t_chg_bounds::CHANGED);
00290 }
00291 }
00292
00293 } else if ((rc <= 0.) &&
00294 (fabs (x-u) <= CLOSE_TO_BOUNDS)) {
00295
00296 if (LB - (u-l) * rc > UB) {
00297
00298 CouNumber newLb = u + (UB-LB) / rc;
00299
00300 newLb = !isInt ? newLb : ceil (newLb - COUENNE_EPS);
00301
00302 if (newLb > Lb (i)) {
00303
00304 Lb (i) = newLb;
00305
00306 nchanges++;
00307 chg_bds [i].setUpper (t_chg_bounds::CHANGED);
00308 }
00309 }
00310 }
00311 }
00312
00313 if (jnlst_ -> ProduceOutput (Ipopt::J_MATRIX, J_BOUNDTIGHTENING)) {
00314 printf ("AFTER reduced cost bt:\n");
00315 for (int i=0, j=0; i < nVars (); ++i) {
00316 if (Var (i) -> Multiplicity () <= 0)
00317 continue;
00318 printf ("%3d [%7e %7e] ", i, Lb (i), Ub (i));
00319 if (!(++j % 4))
00320 printf ("\n");
00321 }
00322 printf ("-----------\n");
00323 }
00324 }
00325
00326 return nchanges;
00327 }