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

Go to the documentation of this file.
00001 /* $Id: problem.cpp 154 2009-06-16 18:52:53Z pbelotti $ */
00002 /*
00003  * Name:    problem.cpp
00004  * Author:  Pietro Belotti
00005  * Purpose: methods of the class CouenneProblem
00006  *
00007  * (C) Carnegie-Mellon University, 2006-08. 
00008  * This file is licensed under the Common Public License (CPL)
00009  */
00010 
00011 #include <vector>
00012 
00013 #include "CoinHelperFunctions.hpp"
00014 #include "CoinTime.hpp"
00015 
00016 #include "CouenneTypes.hpp"
00017 
00018 #include "expression.hpp"
00019 #include "exprConst.hpp"
00020 #include "exprGroup.hpp"
00021 #include "exprClone.hpp"
00022 #include "exprAux.hpp"
00023 #include "lqelems.hpp"
00024 #include "CouenneProblem.hpp"
00025 #include "CouenneProblemElem.hpp"
00026 #include "lqelems.hpp"
00027 
00028 // tricky... smaller values cut the optimum in OS unitTest
00029 const CouNumber SafeCutoff = 1e-4; 
00030 
00033 
00034 void CouenneProblem::initAuxs () const {
00035 
00036   domain_.current () -> resize (nVars ());
00037 
00038   // initially, auxiliary variables are unbounded, their bounds only
00039   // depending on their function
00040 
00041   int nvars = nVars ();
00042 
00043   for (int i=0; i < nvars; i++) {
00044 
00045     int indvar = variables_ [i] -> Index ();
00046 
00047     if ((variables_ [i] -> Type () == AUX) &&                   // this is an auxiliary
00048         (indvar >= nOrigVars_) || // and not an original, originally
00049         (variables_ [i] -> Multiplicity () == 0))               // or a useless one
00050       //int index = variables_ [i] -> Index ();
00051       Lb (indvar) = - (Ub (indvar) = COIN_DBL_MAX);
00052   }
00053 
00054   // first initialize with values from constraints
00055 
00056   //Jnlst()->Printf(Ipopt::J_VECTOR, J_PROBLEM, "Initial bounds for aux (initAuxs):\n");
00057 
00058   for (std::vector <CouenneConstraint *>::const_iterator con = constraints_.begin ();
00059        con != constraints_.end (); ++con) {
00060 
00061     CouNumber
00062       lb = (*((*con) -> Lb ())) (),
00063       ub = (*((*con) -> Ub ())) ();
00064 
00065     int index = (*con) -> Body () -> Index ();
00066 
00067     assert (index >= 0);
00068 
00069     if ((Lb (index) = CoinMax (Lb (index), lb)) <= -COUENNE_INFINITY) Lb (index) = -COIN_DBL_MAX;
00070     if ((Ub (index) = CoinMin (Ub (index), ub)) >=  COUENNE_INFINITY) Ub (index) =  COIN_DBL_MAX;
00071   }
00072 
00073   // only one loop is sufficient here, since auxiliary variable are
00074   // defined in such a way that w_i does NOT depend on w_j if i<j.
00075 
00076   Jnlst () -> Printf (Ipopt::J_MOREMATRIX, J_PROBLEM, "InitAux -- assigning bounds\n");
00077 
00078   for (int j=0, i=nVars (); i--; j++) {
00079 
00080     int ord = numbering_ [j];
00081 
00082     // ignore these variables!
00083     if (variables_ [ord] -> Multiplicity () == 0) {
00084       Lb (ord) = - (Ub (ord) = COIN_DBL_MAX);
00085       X (ord) = 0.;
00086       continue;
00087     }
00088 
00089     // and handle only those with nonzero multiplicity
00090     if (variables_ [ord] -> Type () == AUX) {
00091 
00092       Jnlst () -> Printf (Ipopt::J_MOREMATRIX, J_PROBLEM, 
00093                           "w_%04d [%10g,%10g] ", ord, Lb (ord), Ub (ord));
00094 
00095       CouNumber l, u;
00096 
00097       variables_ [ord] -> Image () -> getBounds (l, u);
00098 
00099       /*printf ("printing bounds: [%g %g]\n", Lb (ord), Ub (ord));
00100       variables_ [ord] -> Lb () -> print (); printf ("\n");
00101       variables_ [ord] -> Ub () -> print (); printf ("\n");*/
00102 
00103       Jnlst () -> Printf (Ipopt::J_MOREMATRIX, J_PROBLEM, 
00104                           " ( --> w_%04d [%10g,%10g] ) vs [%10g %10g]", 
00105                           ord, l, u, Lb (ord), Ub (ord));
00106 
00107       // set bounds 
00108       if ((Lb (ord) = CoinMax (Lb (ord), l)) <= -COUENNE_INFINITY) Lb (ord) = -COIN_DBL_MAX;
00109       if ((Ub (ord) = CoinMin (Ub (ord), u)) >=  COUENNE_INFINITY) Ub (ord) =  COIN_DBL_MAX;
00110       //if ((lb_ [ord] = (*(aux -> Lb ())) ()) <= -COUENNE_INFINITY) lb_ [ord] = -DBL_MAX;
00111       //if ((ub_ [ord] = (*(aux -> Ub ())) ()) >=  COUENNE_INFINITY) ub_ [ord] =  DBL_MAX;
00112 
00113       Jnlst () -> Printf (Ipopt::J_MOREMATRIX, J_PROBLEM, 
00114                           " --> [%10g,%10g]\n", Lb (ord), Ub (ord));
00115 
00116       bool integer = variables_ [ord] -> isInteger ();
00117 
00118       if (integer) {
00119         Lb (ord) = ceil  (Lb (ord) - COUENNE_EPS);
00120         Ub (ord) = floor (Ub (ord) + COUENNE_EPS);
00121       }
00122 
00123       X (ord) = CoinMax (Lb (ord), CoinMin (Ub (ord), (*(variables_ [ord] -> Image ())) ()));
00124     }
00125   }
00126 
00127   restoreUnusedOriginals (); // no argument as the local x vector will be used
00128 }
00129 
00130 
00133 void CouenneProblem::getAuxs (CouNumber * x) const {
00134 
00135   // set point at x, don't copy
00136   domain_.push (nVars (), x, domain_.lb (), domain_.ub (), false);
00137 
00138   // set auxiliary w to f(x). This procedure is exact even though the
00139   // auxiliary variables have an incomplete image, i.e. they have been
00140   // decomposed previously, since they are updated with increasing
00141   // index.
00142 
00143   for (int j=0, i=nVars (); i--; j++) {
00144 
00145     int index = numbering_ [j];
00146     exprVar *var = variables_ [index];
00147 
00148     if (var -> Multiplicity () > 0) {
00149 
00150       CouNumber l, u;
00151 
00152       if (var -> Type () == AUX)
00153         var -> Image () -> getBounds (l,u);
00154       else {
00155         l = Lb (index);
00156         u = Ub (index);
00157       }
00158 
00159       if (var -> Type () == AUX) {
00160         X (index) =  // addresses of x[] and X() are equal
00161           CoinMax (l, CoinMin (u, (*(var -> Image ())) ())); 
00162       }
00163     } else X (index) = 0.;
00164   }
00165 
00166   restoreUnusedOriginals ();
00167 
00168   domain_.pop ();
00169 }
00170 
00171 
00174 
00175 void CouenneProblem::fillObjCoeff (double *&obj) {
00176 
00177   // linearized objective can be an exprAux, an exprSub, an exprGroup,
00178   // or an exprSum. In the last two cases, the components are
00179   // variables or constants
00180 
00181   expression *body = objectives_ [0] -> Body ();
00182   //int sense = objectives_ [0] -> Sense ();
00183 
00184   switch (body -> code ()) {
00185 
00186   case COU_EXPRVAR:   //
00187     obj [body -> Index ()] = 1; //(sense == MINIMIZE) ? 1 : -1;
00188     break;
00189 
00190   case COU_EXPRSUB: { // 
00191 
00192     expression **arglist = body -> ArgList ();
00193 
00194     obj [arglist [0] -> Index ()] =  1; //(sense == MINIMIZE) ?  1 : -1;
00195     obj [arglist [1] -> Index ()] = -1; //(sense == MINIMIZE) ? -1 :  1;
00196 
00197   } break;
00198 
00199   case COU_EXPRGROUP: { // 
00200 
00201     exprGroup *eg    = dynamic_cast <exprGroup *> (body -> isaCopy () ? 
00202                                                    body -> Copy () :
00203                                                    body);
00204 
00205     const exprGroup::lincoeff &lcoe = eg -> lcoeff ();
00206 
00207     //    if (sense == MINIMIZE) while (*index >= 0) obj [*index++] =  *coeff++;
00208     //    else                   while (*index >= 0) obj [*index++] = -*coeff++;      
00209 
00210     for (int n = lcoe.size (), i=0; n--; i++)
00211       //exprGroup::lincoeff::iterator el = lcoe.begin (); el != lcoe.end (); ++el)
00212       obj [lcoe [i]. first -> Index ()] = lcoe [i]. second;
00213     //(sense == MINIMIZE) ? 
00214     //(lcoe [i]. second) : 
00215     //-(lcoe [i]. second);
00216 
00217   } // no break, as exprGroup is derived from exprSum
00218 
00219   case COU_EXPRSUM: { // 
00220 
00221     expression **arglist = body -> ArgList ();
00222 
00223     for (int i = body -> nArgs (); i--;)
00224       switch ((arglist [i]) -> code ()) {
00225 
00226       case COU_EXPRCONST: 
00227         break;
00228 
00229       case COU_EXPRVAR: 
00230         obj [arglist [i] -> Index ()] = 1; //(sense == MINIMIZE) ? 1 : -1;
00231         break;
00232 
00233       case COU_EXPRMUL: {
00234 
00235         expression **mulArgList = arglist [i] -> ArgList ();
00236         int index = mulArgList [0] -> Index ();
00237 
00238         if (index >= 0) obj [index]                      = mulArgList [1] -> Value ();
00239         else            obj [mulArgList [1] -> Index ()] = mulArgList [0] -> Value ();
00240       } break;
00241 
00242       default: 
00243         Jnlst()->Printf(Ipopt::J_ERROR, J_PROBLEM,
00244                         "Couenne: invalid element of sum\nAborting\n");
00245         exit (-1);
00246       }
00247   } break;
00248 
00249   case COU_EXPRCONST: break; // a constant objective
00250 
00251   default:
00252     Jnlst()->Printf(Ipopt::J_WARNING, J_PROBLEM,
00253                     "Couenne: warning, objective function not recognized\n");
00254     break;
00255   }
00256 }
00257 
00258 
00260 void CouenneProblem::setCutOff (CouNumber cutoff) const {
00261 
00262   int indobj = objectives_ [0] -> Body () -> Index ();
00263 
00264   // AW: Should we use the value of the objective variable computed by 
00265   //     Couenne here?
00266   if ((indobj >= 0) && (cutoff < pcutoff_ -> getCutOff () - COUENNE_EPS)) {
00267 
00268     Jnlst () -> Printf (Ipopt::J_DETAILED, J_PROBLEM,
00269                         "Setting new cutoff %.10e for optimization variable index %d val = %.10e\n",
00270                         cutoff, indobj,
00271                         pcutoff_ -> getCutOff ());
00272 
00273     if (Var (indobj) -> isInteger ())
00274       pcutoff_    -> setCutOff (floor (cutoff + COUENNE_EPS));
00275     else pcutoff_ -> setCutOff (cutoff + SafeCutoff * (1. + fabs(cutoff)));
00276   }
00277 } // tolerance needed to retain feasibility
00278 
00279 
00282 void CouenneProblem::installCutOff () const {
00283 
00284   int indobj = objectives_ [0] -> Body () -> Index ();
00285 
00286   if (indobj >= 0) {
00287 
00288     // all problem are assumed to be minimization
00289     double cutoff = pcutoff_ -> getCutOff();
00290 
00291     if (cutoff < Ub (indobj))
00292       Ub (indobj) = cutoff;
00293 
00294   } else jnlst_ -> Printf (J_SUMMARY, J_PROBLEM, 
00295                            "Warning, could not install cutoff - negative objective index\n");
00296 }
00297 
00298 
00299 // clear all spurious variables pointers not referring to the variables_ vector
00300 void CouenneProblem::realign () {
00301 
00302   // link variables to problem's domain
00303   for (std::vector <exprVar *>::iterator i = variables_.begin ();
00304        i != variables_.end (); ++i) {
00305 
00306     (*i) -> linkDomain (&domain_);
00307     (*i) -> realign (this);
00308     if ((*i) -> Type () == AUX)
00309       (*i) -> Image () -> realign (this);
00310   }
00311 
00312   // link variables to problem's domain
00313   for (std::vector <CouenneObjective *>::iterator i = objectives_.begin ();
00314        i != objectives_.end (); ++i) 
00315     (*i) -> Body () -> realign (this);
00316 
00317 
00318   // link variables to problem's domain
00319   for (std::vector <CouenneConstraint *>::iterator i = constraints_.begin ();
00320        i != constraints_.end (); ++i)
00321     (*i) -> Body () -> realign (this);
00322 }
00323 
00324 
00326 void CouenneProblem::registerOptions (Ipopt::SmartPtr <Bonmin::RegisteredOptions> roptions) {
00327 
00328   roptions -> SetRegisteringCategory ("Couenne options", Bonmin::RegisteredOptions::CouenneCategory);
00329 
00330   roptions -> AddNumberOption
00331     ("art_cutoff",
00332      "Artificial cutoff",
00333      COIN_DBL_MAX,
00334      "Default value is infinity.");
00335 
00336   roptions -> AddNumberOption
00337     ("opt_window",
00338      "Window around known optimum",
00339      COIN_DBL_MAX,
00340      "Default value is infinity.");
00341 
00342   roptions -> AddNumberOption
00343     ("feas_tolerance",
00344      "Tolerance for constraints/auxiliary variables",
00345      feas_tolerance_default,
00346      "Default value is zero.");
00347 
00348   roptions -> AddStringOption2 
00349     ("feasibility_bt",
00350      "Feasibility-based (cheap) bound tightening",
00351      "yes",
00352      "no","",
00353      "yes","");
00354 
00355   roptions -> AddStringOption2 
00356     ("redcost_bt",
00357      "Reduced cost bound tightening",
00358      "yes",
00359      "no","",
00360      "yes","");
00361 
00362   roptions -> AddStringOption2 
00363     ("use_quadratic",
00364      "Use quadratic expressions and related exprQuad class",
00365      "no",
00366      "no","Use an auxiliary for each bilinear term",
00367      "yes","Create one only auxiliary for a quadrati expression");
00368 
00369   roptions -> AddStringOption2 
00370     ("optimality_bt",
00371      "Optimality-based (expensive) bound tightening",
00372      "yes",
00373      "no","",
00374      "yes","");
00375 
00376   roptions -> AddLowerBoundedIntegerOption
00377     ("log_num_obbt_per_level",
00378      "Specify the frequency (in terms of nodes) for optimality-based bound tightening.",
00379      -1,1,
00380      "\
00381 If -1, apply at every node (expensive!). \
00382 If 0, apply at root node only. \
00383 If k>=0, apply with probability 2^(k - level), level being the current depth of the B&B tree.");
00384 
00385   roptions -> AddStringOption2 
00386     ("aggressive_fbbt",
00387      "Aggressive feasibility-based bound tightening (to use with NLP points)",
00388      "yes",
00389      "no","",
00390      "yes","");
00391 
00392   roptions -> AddLowerBoundedIntegerOption
00393     ("log_num_abt_per_level",
00394      "Specify the frequency (in terms of nodes) for aggressive bound tightening.",
00395      -1,2,
00396      "\
00397 If -1, apply at every node (expensive!). \
00398 If 0, apply at root node only. \
00399 If k>=0, apply with probability 2^(k - level), level being the current depth of the B&B tree.");
00400 
00401   roptions -> AddNumberOption
00402     ("art_lower",
00403      "Artificial lower bound",
00404      -COIN_DBL_MAX,
00405      "Default value is -COIN_DBL_MAX.");
00406 
00407   roptions -> AddStringOption3
00408     ("branching_object",
00409      "type of branching object for variable selection",
00410      "var_obj",
00411      "vt_obj",   "use Violation Transfer from Tawarmalani and Sahinidis",
00412      "var_obj",  "use one object for each variable",
00413      "expr_obj", "use one object for each nonlinear expression");
00414 
00415   roptions -> AddStringOption2 
00416     ("delete_redundant",
00417      "Eliminate redundant variables, which appear in the problem as x_k = x_h",
00418      "yes",
00419      "no","Keep redundant variables, making the problem a bit larger",
00420      "yes","Eliminate redundant (the problem will be equivalent, only smaller)");
00421 }

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