00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include "CoinHelperFunctions.hpp"
00012 #include "OsiClpSolverInterface.hpp"
00013 #include "OsiCuts.hpp"
00014 #include "CoinTime.hpp"
00015
00016 #include "CouenneFixPoint.hpp"
00017
00018 #include "CouenneProblem.hpp"
00019 #include "CouennePrecisions.hpp"
00020 #include "CouenneExprVar.hpp"
00021 #include "CouenneInfeasCut.hpp"
00022
00023 using namespace Ipopt;
00024 using namespace Couenne;
00025
00027
00028 void CouenneFixPoint::generateCuts (const OsiSolverInterface &si,
00029 OsiCuts &cs,
00030 const CglTreeInfo treeInfo)
00031 #if CGL_VERSION_MAJOR == 0 && CGL_VERSION_MINOR <= 57
00032 const
00033 #endif
00034 {
00035
00044
00045 if (firstCall_)
00046 firstCall_ = false;
00047 else
00048 if (!(problem_ -> fbbtReachedIterLimit ()))
00049 return;
00050
00051 if (isWiped (cs))
00052 return;
00053
00054 if (treeInfo.inTree &&
00055 treeInfo.level > 0 &&
00056 treeInfo.pass > 1)
00057 return;
00058
00059 int nInitTightened = nTightened_;
00060
00061 if (treeInfo.inTree &&
00062 treeInfo.level <= 0) {
00063
00064 problem_ -> Jnlst () -> Printf (J_ERROR, J_COUENNE, "Fixed Point FBBT: ");
00065 fflush (stdout);
00066 }
00067
00068 ++nRuns_;
00069
00070 double now = CoinCpuTime ();
00071
00072 problem_ -> domain () -> push (&si, &cs);
00073
00092
00093 OsiSolverInterface *fplp = NULL;
00094
00095 if (true) {
00096
00097
00098 fplp = new OsiClpSolverInterface;
00099 }
00100
00117
00119
00120 const CoinPackedMatrix *A = si. getMatrixByRow ();
00121
00122 const int
00123 n = si. getNumCols (),
00124 m = si. getNumRows (),
00125 nCuts = cs.sizeRowCuts (),
00126 *ind = A -> getIndices ();
00127
00128 const double
00129 *lb = si. getColLower (),
00130 *ub = si. getColUpper (),
00131 *rlb = si. getRowLower (),
00132 *rub = si. getRowUpper (),
00133 *coe = A -> getElements ();
00134
00135 if (problem_ -> Jnlst () -> ProduceOutput (J_ERROR, J_BOUNDTIGHTENING))
00136 for (int i=0; i<n; i++)
00137 printf ("----------- x_%d in [%g,%g]\n", i, lb [i], ub [i]);
00138
00139
00140 fplp -> messageHandler () -> setLogLevel (0);
00141
00142
00143 for (int i=0; i<n; i++) {
00144 bool isActive = problem_ -> Var (i) -> Multiplicity () > 0;
00145 fplp -> addCol (0, NULL, NULL, lb [i], ub [i], isActive ? -1. : 0.);
00146 }
00147
00148 for (int i=0; i<n; i++) {
00149 bool isActive = problem_ -> Var (i) -> Multiplicity () > 0;
00150 fplp -> addCol (0, NULL, NULL, lb [i], ub [i], isActive ? +1. : 0.);
00151 }
00152
00153 if (extendedModel_) {
00154
00155 for (int j=0; j<m; j++) fplp -> addCol (0, NULL, NULL, rlb [j], COIN_DBL_MAX, 0.);
00156 for (int j=0; j<m; j++) fplp -> addCol (0, NULL, NULL, -COIN_DBL_MAX, rub [j], 0.);
00157 }
00158
00159
00160
00161 for (int j=0; j<m; j++) {
00162
00163 int nEl = A -> getVectorSize (j);
00164
00165 if (!nEl)
00166 continue;
00167
00168 if (problem_ -> Jnlst () -> ProduceOutput (J_ERROR, J_BOUNDTIGHTENING)) {
00169
00170 printf ("row %4d, %4d elements: ", j, nEl);
00171
00172 for (int i=0; i<nEl; i++) {
00173 printf ("%+g x%d ", coe [i], ind [i]);
00174 fflush (stdout);
00175 }
00176
00177 printf ("\n");
00178 }
00179
00180
00181
00182 if (extendedModel_ || rlb [j] > -COUENNE_INFINITY)
00183 for (int i=0; i<nEl; i++)
00184 createRow (-1, ind [i], n, fplp, ind, coe, rlb [j], nEl, extendedModel_, j, m + nCuts);
00185
00186 if (extendedModel_ || rub [j] < COUENNE_INFINITY)
00187 for (int i=0; i<nEl; i++)
00188 createRow (+1, ind [i], n, fplp, ind, coe, rub [j], nEl, extendedModel_, j, m + nCuts);
00189
00190
00191
00192 if (extendedModel_) {
00193 createRow (-1, 2*n + j, n, fplp, ind, coe, rlb [j], nEl, extendedModel_, j, m + nCuts);
00194 createRow (+1, 2*n + m + j, n, fplp, ind, coe, rub [j], nEl, extendedModel_, j, m + nCuts);
00195 }
00196
00197 ind += nEl;
00198 coe += nEl;
00199 }
00200
00201
00202
00203 for (int j = 0, jj = nCuts; jj--; j++) {
00204
00205
00206
00207 OsiRowCut *cut = cs.rowCutPtr (j);
00208
00209 const CoinPackedVector &row = cut -> row ();
00210
00211 const int
00212 nEl = row.getNumElements (),
00213 *ind = row.getIndices ();
00214
00215 const double *coe = row.getElements ();
00216
00217 if (extendedModel_) {
00218 fplp -> addCol (0, NULL, NULL, cut -> lb (), COIN_DBL_MAX, 0.);
00219 fplp -> addCol (0, NULL, NULL, -COIN_DBL_MAX, cut -> ub (), 0.);
00220 }
00221
00222 if (extendedModel_ || cut -> lb () > -COUENNE_INFINITY)
00223 for (int i=0; i<nEl; i++)
00224 createRow (-1, ind [i], n, fplp, ind, coe, cut -> lb (), nEl, extendedModel_, m + j, m + nCuts);
00225
00226 if (extendedModel_ || cut -> ub () < COUENNE_INFINITY)
00227 for (int i=0; i<nEl; i++)
00228 createRow (+1, ind [i], n, fplp, ind, coe, cut -> ub (), nEl, extendedModel_, m + j, m + nCuts);
00229
00230
00231
00232 if (extendedModel_) {
00233 createRow (-1, 2*n + j, n, fplp, ind, coe, cut -> lb (), nEl, extendedModel_, m + j, m + nCuts);
00234 createRow (+1, 2*n + m + nCuts + j, n, fplp, ind, coe, cut -> ub (), nEl, extendedModel_, m + j, m + nCuts);
00235 }
00236
00237 ind += nEl;
00238 coe += nEl;
00239 }
00240
00241
00242
00243 if (extendedModel_)
00244
00245 for (int j=0; j<m; j++) {
00246
00247 int ind [2] = {2*n + j, 2*n + m + j};
00248 double coe [2] = {1., -1.};
00249
00250 CoinPackedVector row (2, ind, coe);
00251 fplp -> addRow (row, -COIN_DBL_MAX, 0.);
00252 }
00253
00256
00257 fplp -> setObjSense (-1.);
00258
00259
00260
00261
00262 fplp -> initialSolve ();
00263
00264 if (fplp -> isProvenOptimal ()) {
00265
00266
00267
00268 const double
00269 *newLB = fplp -> getColSolution (),
00270 *newUB = newLB + n,
00271 *oldLB = si. getColLower (),
00272 *oldUB = si. getColUpper ();
00273
00274
00275
00276 int
00277 *indLB = new int [n],
00278 *indUB = new int [n],
00279 ntightenedL = 0,
00280 ntightenedU = 0;
00281
00282 double
00283 *valLB = new double [n],
00284 *valUB = new double [n];
00285
00286 for (int i=0; i<n; i++) {
00287
00288 if (problem_ -> Jnlst () -> ProduceOutput (J_ERROR, J_BOUNDTIGHTENING))
00289 printf ("x%d: [%g,%g] --> [%g,%g]\n", i,
00290 oldLB [i], oldUB [i],
00291 newLB [i], newUB [i]);
00292
00293 if (newLB [i] > oldLB [i] + COUENNE_EPS) {
00294
00295 indLB [ntightenedL] = i;
00296 valLB [ntightenedL++] = newLB [i];
00297
00298 ++nTightened_;
00299 }
00300
00301 if (newUB [i] < oldUB [i] - COUENNE_EPS) {
00302
00303 indUB [ntightenedU] = i;
00304 valUB [ntightenedU++] = newUB [i];
00305
00306 ++nTightened_;
00307 }
00308 }
00309
00310 if (ntightenedL || ntightenedU) {
00311
00312 OsiColCut newBound;
00313
00314 newBound.setLbs (ntightenedL, indLB, valLB);
00315 newBound.setUbs (ntightenedU, indUB, valUB);
00316
00317 cs.insert (newBound);
00318 }
00319
00320 delete [] indLB;
00321 delete [] indUB;
00322 delete [] valLB;
00323 delete [] valUB;
00324
00325 CPUtime_ += CoinCpuTime () - now;
00326
00327 if (treeInfo.inTree &&
00328 treeInfo.level <= 0) {
00329
00330 problem_ -> Jnlst () -> Printf (J_ERROR, J_COUENNE, "%d bounds tightened (%g seconds)\n",
00331 nTightened_ - nInitTightened, CoinCpuTime () - now);
00332 }
00333
00334 } else
00335 if (treeInfo.inTree &&
00336 treeInfo.level <= 0)
00337 problem_ -> Jnlst () -> Printf (J_ERROR, J_COUENNE, " FPLP infeasible or unbounded.\n");
00338
00339 delete fplp;
00340
00341 problem_ -> domain () -> pop ();
00342 }
00343
00344
00345
00346
00347
00348
00349
00350
00351
00352
00353
00354
00355
00356
00357
00358
00359
00360
00361 void CouenneFixPoint::createRow (int sign,
00362 int indexVar,
00363 int nVars,
00364 OsiSolverInterface *p,
00365 const int *indices,
00366 const double *coe,
00367 const double rhs,
00368 const int nEl,
00369 bool extMod,
00370 int indCon,
00371 int nCon) const {
00372
00474
00475
00476 if (problem_ -> Jnlst () -> ProduceOutput (J_ERROR, J_BOUNDTIGHTENING)) {
00477
00478 printf ("creating constraint from: ");
00479
00480 for (int i=0; i<nEl; i++)
00481 printf ("%+g x%d ", coe [i], indices [i]);
00482
00483 printf ("%c= %g for variable index %d: ", sign > 0 ? '<' : '>', rhs, indexVar);
00484 }
00485
00486 int nTerms = nEl;
00487
00488 if (extMod)
00489 nTerms++;
00490
00491 int *iInd = new int [nTerms];
00492 double *elem = new double [nTerms];
00493
00494
00495
00496 CoinCopyN (coe, nEl, elem);
00497
00498 if (extMod) {
00499 elem [nEl] = -1.;
00500 iInd [nEl] = 2*nVars + indCon + ((sign > 0) ? nCon : 0);
00501 }
00502
00503
00504
00505 for (int k=0; k<nEl; k++) {
00506
00507 int curInd = indices [k];
00508
00509 iInd [k] = curInd;
00510
00511
00512
00513 if (curInd == indexVar) {
00514 if (((sign > 0) && (coe [k] > 0.)) ||
00515 ((sign < 0) && (coe [k] < 0.)))
00516
00517 iInd [k] += nVars;
00518
00519 } else if (((coe [k] > 0.) && (sign < 0)) ||
00520 ((coe [k] < 0.) && (sign > 0)))
00521
00522 iInd [k] += nVars;
00523 }
00524
00525 CoinPackedVector vec (nTerms, iInd, elem);
00526
00527 double
00528 lb = sign > 0 ? -COIN_DBL_MAX : extMod ? 0. : rhs,
00529 ub = sign < 0 ? +COIN_DBL_MAX : extMod ? 0. : rhs;
00530
00531 p -> addRow (vec, lb, ub);
00532
00533
00534
00535 if (problem_ -> Jnlst () -> ProduceOutput (J_ERROR, J_BOUNDTIGHTENING)) {
00536
00537 for (int i=0; i<nEl; i++)
00538 printf ("%+g x%d ", elem [i], iInd [i]);
00539
00540 printf ("in [%g,%g]\n", lb, ub);
00541 }
00542
00543
00544
00545
00546
00547 delete [] iInd;
00548 delete [] elem;
00549 }