/home/coin/SVN-release/OS-2.0.1/Couenne/src/disjunctive/getDisjunctions.cpp

Go to the documentation of this file.
00001 /* $Id: getDisjunctions.cpp 141 2009-06-03 04:19:19Z pbelotti $ */
00002 /*
00003  * Name:    getDisjunctions.cpp
00004  * Author:  Pietro Belotti
00005  * Purpose: get nonlinear (and integer?) disjunctions for the problem
00006  *
00007  * (C) Carnegie-Mellon University, 2008. 
00008  * This file is licensed under the Common Public License (CPL)
00009  */
00010 
00011 #include "CouenneCutGenerator.hpp"
00012 #include "CouenneProblem.hpp"
00013 #include "CouenneObject.hpp"
00014 #include "CouenneBranchingObject.hpp"
00015 #include "CouenneDisjCuts.hpp"
00016 
00017 
00019 int CouenneDisjCuts::getDisjunctions (std::vector <std::pair <OsiCuts *, OsiCuts *> > &disjs, 
00020                                       OsiSolverInterface &si, 
00021                                       OsiCuts &cs, 
00022                                       const CglTreeInfo &info) const {
00023 
00024   OsiBranchingInformation brInfo (&si, 
00025                                   true,   // bool normalSolver,
00026                                   false); // bool copySolution=false);
00027 
00028   if (jnlst_ -> ProduceOutput (J_MATRIX, J_DISJCUTS))
00029     for (int i=0; i<si.getNumCols (); i++)
00030       printf ("%3d. %8g  [%8g %8g]\n", i, 
00031               brInfo. solution_ [i], 
00032               brInfo. lower_    [i], 
00033               brInfo. upper_    [i]);
00034 
00035   // get list of best disjunction as though we were branching
00036   branchingMethod_ -> setupList (&brInfo, true); // initialize
00037 
00038   int 
00039     ncols  = si.getNumCols (),
00040     retval = COUENNE_FEASIBLE;
00041 
00042   int nobjs = branchingMethod_ -> numberUnsatisfied ();
00043 
00044   if (nobjs) {
00045 
00046     if (jnlst_ -> ProduceOutput (J_DETAILED, J_DISJCUTS)) 
00047       printf ("---   %d disjunctive objects\n", nobjs);
00048 
00049     const int *candidates = branchingMethod_ -> candidates ();
00050     OsiObject **objs = si. objects ();
00051 
00052     // enumerate first (most infeasible, or most promising) objects of
00053     // the list
00054 
00055     for (int candInd = 0,      num = 0; 
00056          (candInd < nobjs) && (num < numDisjunctions_); 
00057          candInd++) {
00058 
00059       CouenneObject *cObj = dynamic_cast <CouenneObject *> (objs [candidates [candInd]]);
00060 
00061       // TODO: consider also integer objects
00062 
00063       if (!cObj ||                                        // IF    not a nonlinear object
00064           (cObj -> checkInfeasibility (&brInfo) == 0.) || //    or nonlinear, but not violated
00065           (cObj -> isCuttable ()))                        //    or we are on the "good" (convex) side
00066         continue;                                         // THEN skip
00067 
00068       bool     feasLeft = true,  feasRight = true;
00069       OsiCuts *leftCuts = NULL, *rightCuts = NULL;
00070 
00071       const double 
00072         *saveLower = CoinCopyOfArray (si.getColLower (), ncols),
00073         *saveUpper = CoinCopyOfArray (si.getColUpper (), ncols);
00074 
00075       OsiBranchingObject *brObj = cObj -> createBranch (&si, &brInfo, 0); // down!
00076 
00077       if (jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS))
00078         printf ("---   cand [%d] is %d: x_%d <>= %g [%g,%g]\n", 
00079               candInd, 
00080                 candidates [candInd],
00081                 dynamic_cast <CouenneBranchingObject *> (brObj) -> variable () -> Index (), 
00082                 brObj -> value (),
00083                 couenneCG_->Problem()->Lb(dynamic_cast<CouenneBranchingObject*>(brObj)
00084                                           ->variable()->Index ()),
00085                 couenneCG_->Problem()->Ub(dynamic_cast<CouenneBranchingObject*>(brObj)
00086                                           ->variable()->Index ()));
00087 
00088       // left branch ////////////////////////////////////////////////////////////
00089 
00090       if (brObj -> branch (&si) >= COUENNE_INFINITY)
00091         feasLeft = false;                        // left subtree infeasible
00092       else leftCuts = getSingleDisjunction (si); // feasible, store new bounds in OsiCuts
00093 
00094       // restore initial bounds
00095       si.setColLower (saveLower);
00096       si.setColUpper (saveUpper);
00097 
00098       delete brObj;
00099 
00100       // right branch ////////////////////////////////////////////////////////////
00101 
00102       brObj = cObj -> createBranch (&si, &brInfo, 1); // up!
00103 
00104       if (brObj -> branch (&si) >= COUENNE_INFINITY) 
00105         feasRight = false;                        // right subtree infeasible
00106       else rightCuts = getSingleDisjunction (si); // feasible, store new bounds in OsiCuts
00107 
00108       delete brObj;
00109 
00110       // restore initial bounds
00111       si.setColLower (saveLower);
00112       si.setColUpper (saveUpper);
00113 
00114       // done with generating a disjunction. Is any (or both) side infeasible?
00115 
00116       if (feasLeft && feasRight) {
00117 
00118         // feasible, which variables tightened?
00119 
00120         std::pair <OsiCuts *, OsiCuts *> newpair;
00121         newpair.first  = leftCuts;
00122         newpair.second = rightCuts;
00123 
00124         disjs.push_back (newpair);
00125 
00126         num++;
00127 
00128       } else if (feasLeft) {
00129 
00130         if (jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS))
00131           printf ("---   disj: infeasible right\n");
00132         applyColCuts (si, leftCuts); // !!! culprit
00133         retval = COUENNE_TIGHTENED;
00134 
00135       } else if (feasRight) {
00136 
00137         if (jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS))
00138           printf ("---   infeasible left\n");
00139         applyColCuts (si, rightCuts); // !!! culprit
00140         retval = COUENNE_TIGHTENED;
00141 
00142       } else {
00143           if (jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS))
00144             printf ("---   infeasible!!!\n");
00145         retval = COUENNE_INFEASIBLE;
00146         break;
00147       }
00148 
00149       delete [] saveLower;
00150       delete [] saveUpper;
00151     }
00152   }
00153 
00154   // disjunctions (of bounding boxes) generated. Do some extra
00155   // tightening that may result from tightened sides of the
00156   // disjunctions
00157 
00158   if (retval != COUENNE_INFEASIBLE) {
00159 
00160     if (jnlst_ -> ProduceOutput (J_DETAILED, J_DISJCUTS))
00161       printf ("have %d disjunctions\n", disjs.size ());
00162 
00163     // sanity check: for each disjunction, check if left and right
00164     // part are feasible with (possibly improved) bounds of si --
00165     // these may have tightened and now exclude either (or both!)
00166     // sides of a disjunction
00167 
00168     for (std::vector <std::pair <OsiCuts *, OsiCuts *> >::iterator disjI = disjs.begin ();
00169          disjI != disjs.end (); ++disjI) {
00170 
00171       if   (checkDisjSide (si, disjI -> first)  == COUENNE_INFEASIBLE) { // left  side infeasible?
00172         if (checkDisjSide (si, disjI -> second) == COUENNE_INFEASIBLE) { // right side infeasible?
00173 
00174           retval = COUENNE_INFEASIBLE;
00175           break;
00176 
00177         } else {
00178 
00179           if (jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS))
00180             printf ("---   infeasible left [2]!\n");
00181           applyColCuts (si, disjI -> second);
00182           retval = COUENNE_TIGHTENED;
00183         }
00184       } else if (checkDisjSide (si, disjI -> second) == COUENNE_INFEASIBLE) { // right side
00185 
00186         if (jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS))
00187           printf ("---   infeasible right [2]!\n");
00188         applyColCuts (si, disjI -> first);
00189         retval = COUENNE_TIGHTENED;
00190       }
00191     }
00192 
00193     if (retval == COUENNE_INFEASIBLE)
00194       return retval;
00195 
00196     // merge tightened disjunction by intersecting si's bounding box
00197     // with intersections of smallest boxes, one per disjunction,
00198     // containing each both sides of the disjunction
00199 
00200     for (std::vector <std::pair <OsiCuts *, OsiCuts *> >::iterator disjI = disjs.begin ();
00201          disjI != disjs.end (); ++disjI) {
00202 
00203       CoinPackedVector lowerChg, upperChg;
00204 
00205       // find smallest box containing two disjunctions
00206       getBoxUnion (si, disjI -> first, disjI -> second, lowerChg, upperChg);
00207 
00208       if ((lowerChg.getNumElements () > 0) || 
00209           (upperChg.getNumElements () > 0)) {
00210 
00211         OsiColCut cut;
00212         cut.setLbs (lowerChg);
00213         cut.setUbs (upperChg);
00214         applyColCuts (si, &cut);
00215 
00216         cs.insert (cut);
00217       }
00218     }   
00219   }
00220 
00221   return retval;
00222 }
00223 
00224 
00226 void CouenneDisjCuts::applyColCuts (OsiSolverInterface &si, OsiCuts *cuts) const {
00227 
00228   if (jnlst_ -> ProduceOutput (J_MATRIX, J_DISJCUTS)) {
00229     printf ("applying cuts to SI:\n");
00230     // apply column cuts to si
00231     for (int i = cuts -> sizeColCuts (); i--;)
00232       cuts -> colCutPtr (i) -> print ();
00233     printf ("--------------------\n");
00234   }
00235 
00236   // apply column cuts to si
00237   for (int i = cuts -> sizeColCuts (); i--;)
00238     applyColCuts (si, cuts -> colCutPtr (i));
00239 }
00240 
00241 
00243 void CouenneDisjCuts::applyColCuts (OsiSolverInterface &si, OsiColCut *cut) const {
00244 
00245   // apply single column cut to si
00246   if (jnlst_ -> ProduceOutput (J_VECTOR, J_DISJCUTS)) {
00247     printf ("tightening bounds: "); 
00248     cut -> print ();
00249   }
00250 
00251   const CoinPackedVector
00252     &lbs = cut -> lbs (),
00253     &ubs = cut -> ubs ();
00254 
00255   const int    *lind = lbs.getIndices  (), *uind = ubs.getIndices  ();
00256   const double *lval = lbs.getElements (), *uval = ubs.getElements (),
00257     *oldLower = si.getColLower (), 
00258     *oldUpper = si.getColUpper ();
00259 
00260   for (int j=lbs.getNumElements (); j--; lval++, lind++)
00261     if (*lval > oldLower [*lind] + COUENNE_EPS) 
00262       si.setColLower (*lind, *lval);
00263 
00264   for (int j=ubs.getNumElements (); j--; uval++, uind++)
00265     if (*uval < oldUpper [*uind] - COUENNE_EPS) 
00266       si.setColUpper (*uind, *uval);
00267 }

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