Changeset 553


Ignore:
Timestamp:
May 26, 2007 1:27:17 PM (12 years ago)
Author:
pbelotti
Message:

fixed obj bound management in tightening

Location:
branches/Couenne/Couenne/src
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • branches/Couenne/Couenne/src/convex/CouenneCutGenerator.cpp

    r548 r553  
    5555
    5656  now = CoinCpuTime ();
    57   //  problem_ -> print (std::cout);
    58   //  printf ("======================================\n");
     57  //problem_ -> print (std::cout);
     58  //printf ("======================================\n");
    5959  problem_ -> standardize ();
    6060
     
    6464  septime_ = now;
    6565
    66   //  problem_ -> print (std::cout);
     66  //problem_ -> print (std::cout);
    6767
    6868  //problem_ -> writeMod ("extended-aw.mod", true);
  • branches/Couenne/Couenne/src/convex/boundTightening.cpp

    r548 r553  
    2929  /////////////////////// MIP bound management /////////////////////////////////
    3030
    31   if (objInd < 0) objInd = 0;
     31  if ((objInd >= 0) && babInfo && (babInfo -> babPtr ())) {
    3232
    33   if (babInfo && (babInfo -> babPtr ())) {
    34 
    35     CouNumber primal  = babInfo  -> babPtr () -> model (). getObjValue(),
    36               dual    = babInfo  -> babPtr () -> model (). getBestPossibleObjValue (),
     33    CouNumber UB      = babInfo  -> babPtr () -> model (). getObjValue(),
     34              LB      = babInfo  -> babPtr () -> model (). getBestPossibleObjValue (),
    3735              primal0 = problem_ -> Ub (objInd),
    3836              dual0   = problem_ -> Lb (objInd);
     
    4038    // Bonmin assumes minimization. Hence, primal (dual) is an UPPER
    4139    // (LOWER) bound.
     40   
     41    if ((UB < COUENNE_INFINITY) &&
     42        (UB < primal0 - COUENNE_EPS)) {
     43      // update primal bound (MIP)
    4244
    43     if ((primal <   COUENNE_INFINITY - 1) && (primal < primal0)) { // update primal bound (MIP)
    44       problem_ -> Ub (objInd) = primal;
     45      //      printf ("updating upper %g <-- %g\n", primal0, primal);
     46      problem_ -> Ub (objInd) = UB;
     47      chg_bds [objInd] = 1;
     48    }
     49   
     50    if ((LB   > - COUENNE_INFINITY) &&
     51        (LB   > dual0 + COUENNE_EPS)) { // update dual bound
     52      problem_ -> Lb (objInd) = LB;
    4553      chg_bds [objInd] = 1;
    4654    }
    4755
    48     if ((dual   > - COUENNE_INFINITY + 1) && (dual   > dual0)) { // update dual bound
    49       problem_ -> Lb (objInd) = dual;
    50       chg_bds [objInd] = 1;
    51     }
    52   }
     56    //////////////////////// Reduced cost bound tightening //////////////////////
    5357
    54   //////////////////////// Reduced cost bound tightening //////////////////////
     58    // do it only if a linear convexification is in place already
    5559
    56   // do it only if a linear convexification is in place already
    57 
    58   if ((objInd >= 0) && !firstcall_) {
    59 
     60    if (!firstcall_) {
     61      /*
    6062    CouNumber
    6163      LB = si.getObjValue (),
     
    6365              babInfo  -> babPtr () -> model (). getObjValue() :
    6466              problem_ -> Ub (objInd);
     67      */
    6568
    66     if ((LB > -COUENNE_INFINITY) && (UB < COUENNE_INFINITY)) {
    67       int ncols = si.getNumCols ();
     69      if ((LB > -COUENNE_INFINITY) && (UB < COUENNE_INFINITY)) {
     70        int ncols = si.getNumCols ();
    6871
    69       for (int i=0; i<ncols; i++) {
     72        for (int i=0; i<ncols; i++) {
    7073
    71         CouNumber
    72           x  = si.getColSolution () [i],
    73           rc = si.getReducedCost () [i],
    74           dx = problem_ -> Ub (i) - x;
     74          CouNumber
     75            x  = si.getColSolution () [i],
     76            rc = si.getReducedCost () [i],
     77            dx = problem_ -> Ub (i) - x;
    7578
    76         if ((rc > COUENNE_EPS) && (dx*rc > (UB-LB))) {
    77           // can improve bound on variable w_i
    78           problem_ -> Ub (i) = x + (UB-LB) / rc;
    79           chg_bds [i] = 1;
     79          if ((rc > COUENNE_EPS) && (dx*rc > (UB-LB))) {
     80            // can improve bound on variable w_i
     81            problem_ -> Ub (i) = x + (UB-LB) / rc;
     82            chg_bds [i] = 1;
     83          }
    8084        }
    8185      }
  • branches/Couenne/Couenne/src/convex/createCuts.cpp

    r551 r553  
    1616#include <CouenneProblem.h>
    1717
    18 
    19 //#define MAX_COEFF 1e12
    2018
    2119/// general procedure for inserting a linear cut with up to three
     
    6058    if (i3 >= 0) violation += c3 * x [i3];
    6159
    62     // return NULL if not violated
     60    // return 0 if not violated
    6361
    6462    if (((violation <   COUENNE_EPS) || (sign > 0)) &&
     
    115113    delete [] index;
    116114
    117     // some convexification cuts (as the lower envelopes of convex
    118     // functions) are global, hence here is a tool to make them valid
    119     // throughout the BB tree
    120 
    121     cut -> setGloballyValid (is_global);
     115    cut -> setGloballyValid (is_global); // global?
    122116    cs.insert (cut);
    123117    delete cut;
  • branches/Couenne/Couenne/src/convex/genColCuts.cpp

    r497 r553  
    2121    *indLow = new int [ncols], // indices for OsiColCut
    2222    *indUpp = new int [ncols],
    23     nLow, nUpp = nLow = 0;
     23    nLow, nUpp = nLow = 0,
     24    ind_obj = problem_ -> Obj (0) -> Body () -> Index ();
    2425
    2526  // values fo OsiColCut
     
    3738
    3839    register int index = changed [i];
     40
     41    if (index == ind_obj) continue;
    3942
    4043    CouNumber bd;
  • branches/Couenne/Couenne/src/convex/generateCuts.cpp

    r552 r553  
    2626#define MAX_OBBT_ITER 4
    2727
     28#define LARGE_TOL (LARGE_BOUND/1e6)
     29
     30//
     31void fictitiousBound (OsiCuts &cs,
     32                      CouenneProblem *p,
     33                      bool action) { // true before convexifying, false afterwards
     34
     35  // set trivial dual bound to objective function, if there is none
     36
     37  int ind_obj = p -> Obj (0) -> Body () -> Index ();
     38  if (ind_obj < 0) return;
     39  int sense = (p -> Obj (0) -> Sense () == MINIMIZE) ? -1 : 1;
     40
     41  if (action) {
     42
     43    if (sense < 0) {if (p -> Lb (ind_obj) < - LARGE_BOUND) p -> Lb (ind_obj) = - LARGE_BOUND;}
     44    else            if (p -> Ub (ind_obj) >   LARGE_BOUND) p -> Ub (ind_obj) =   LARGE_BOUND;
     45  }
     46  else
     47    if (sense>0) {if (fabs (p->Ub(ind_obj)-LARGE_BOUND)<LARGE_TOL) p->Ub(ind_obj) = COUENNE_INFINITY;}
     48    else          if (fabs (p->Lb(ind_obj)+LARGE_BOUND)<LARGE_TOL) p->Lb(ind_obj) = -COUENNE_INFINITY;
     49}
     50
    2851
    2952// translate changed bound sparse array into a dense one
     
    4972                                        OsiCuts &cs,
    5073                                        const CglTreeInfo info) const {
     74
     75  //  if (firstcall_)
     76  //    fictitiousBound (si, problem_, true);
     77
    5178  int nInitCuts = cs.sizeRowCuts ();
    5279
     
    135162        const double * beforeUpper = auxinfo -> beforeUpper ();
    136163
    137         if (beforeLower && beforeUpper) {
     164        if (beforeLower || beforeUpper) {
    138165
    139166          // get currentbounds
     
    143170          for (register int i=0; i < ncols; i++)
    144171
    145             if ((   nowLower [i] >= beforeLower [i] + COUENNE_EPS)
    146                 || (nowUpper [i] <= beforeUpper [i] - COUENNE_EPS))
     172            if (beforeLower && (nowLower [i] >= beforeLower [i] + COUENNE_EPS) ||
     173                beforeUpper && (nowUpper [i] <= beforeUpper [i] - COUENNE_EPS))
    147174              chg_bds [i] = 1;
    148175
     
    151178    }
    152179  }
     180
     181  fictitiousBound (cs, problem_, false);
    153182
    154183  int *changed = NULL, nchanged;
     
    168197  //////////////////////
    169198
     199  /*
    170200  if (firstcall_) {
    171201
     
    176206    if (ind_obj >= 0) {
    177207      if (problem_ -> Obj (0) -> Sense () == MINIMIZE) {
    178         if (problem_ -> Lb (ind_obj) < - LARGE_BOUND)
    179           problem_ -> Lb (ind_obj) = - LARGE_BOUND;
     208       if (problem_ -> Lb (ind_obj) < - LARGE_BOUND)
     209         problem_ -> Lb (ind_obj) = - LARGE_BOUND;
    180210      }
    181211      else
    182         if (problem_ -> Ub (ind_obj) > LARGE_BOUND)
    183           problem_ -> Ub (ind_obj) = LARGE_BOUND;
    184     }
    185   }
     212       if (problem_ -> Ub (ind_obj) > LARGE_BOUND)
     213         problem_ -> Ub (ind_obj) = LARGE_BOUND;
     214    }
     215  }
     216  */
    186217
    187218  // change tightened bounds through OsiCuts
     
    232263 end_genCuts:
    233264
     265  if (firstcall_)
     266    fictitiousBound (cs, problem_, true);
     267
    234268  if (firstcall_) {
    235269    firstcall_  = false;
  • branches/Couenne/Couenne/src/convex/obbt.cpp

    r552 r553  
    7575      csi -> setObjective (objcoe);
    7676
    77       // minimize and then maximize x_i on si.
     77      // minimize...
    7878
    79       if (updateBound (csi,  1, problem_ -> Lb (i), isInt)) {
     79      if (updateBound (csi, 1, problem_ -> Lb (i), isInt)) {
    8080        csi -> setColLower (i, problem_ -> Lb (i));
    8181        chg = true;
    8282      }
    8383     
     84      // ...and then maximize x_i on csi.
     85
    8486      if (updateBound (csi, -1, problem_ -> Ub (i), isInt)) {
    85         csi -> setColUpper (i, problem_ -> Ub (i));
     87        csi -> setColUpper  (i, problem_ -> Ub (i));
    8688        chg = true;
    8789      }
  • branches/Couenne/Couenne/src/convex/operators/conv-exprSinCos.cpp

    r551 r553  
    1919
    2020/// convex cuts for sine or cosine
    21 void trigEnvelope (const CouenneCutGenerator *, OsiCuts &,
     21int trigEnvelope (const CouenneCutGenerator *, OsiCuts &,
    2222                   exprAux *, expression *, enum cou_trig);
    2323
    2424
    2525///
    26 void addHexagon (const CouenneCutGenerator *, // cut generator that has called us
     26int addHexagon (const CouenneCutGenerator *, // cut generator that has called us
    2727                 OsiCuts &,      // cut set to be enriched
    2828                 enum cou_trig,  // sine or cosine
     
    3737
    3838#ifdef NEW_TRIG
    39   trigEnvelope (cg, cs, w, w -> Image () -> Argument (), COU_SINE);
     39  if (trigEnvelope (cg, cs, w, w -> Image () -> Argument (), COU_SINE) == 0)
    4040#else
    41   addHexagon (cg, cs, COU_SINE, w, w -> Image () -> Argument());
     41  if (addHexagon (cg, cs, COU_SINE, w, w -> Image () -> Argument()) == 0)
    4242#endif
     43    {
     44
     45    }
    4346}
    4447
     
    5053
    5154#ifdef NEW_TRIG
    52   trigEnvelope (cg, cs, w, w -> Image () -> Argument (), COU_COSINE);
     55  if (trigEnvelope (cg, cs, w, w -> Image () -> Argument (), COU_COSINE) == 0)
    5356#else
    54   addHexagon (cg, cs, COU_COSINE, w, w -> Image () -> Argument());
     57  if (addHexagon (cg, cs, COU_COSINE, w, w -> Image () -> Argument()) == 0)
    5558#endif
     59    {
     60
     61    }
    5662}
    5763
     
    5965/// add lateral edges of the hexagon providing
    6066
    61 void addHexagon (const CouenneCutGenerator *cg, // cut generator that has called us
     67int addHexagon (const CouenneCutGenerator *cg, // cut generator that has called us
    6268                 OsiCuts &cs,       // cut set to be enriched
    6369                 enum cou_trig tt,  // sine or cosine
     
    6773  unary_function fn = (tt == COU_SINE) ? sin : cos;
    6874
     75  // retrieve argument bounds
    6976  expression *lbe, *ube;
    7077  arg -> getBounds (lbe, ube);
     
    7380            ub = (*ube) ();
    7481
    75   int x_ind = arg -> Index ();
    76   int w_ind = aux -> Index ();
     82  delete lbe;
     83  delete ube;
     84
     85  int ncuts = 0,
     86    x_ind = arg -> Index (),
     87    w_ind = aux -> Index ();
    7788
    7889  if (fabs (ub - lb) < COUENNE_EPS) {
     
    8394    else                {f = cos (x0); fp = -sin (x0);}
    8495
    85     cg -> createCut (cs, f - fp*x0, 0, w_ind, 1., x_ind, -fp);
    86     return;
     96    return cg -> createCut (cs, f - fp*x0, 0, w_ind, 1., x_ind, -fp);
    8797  }
    8898
     
    92102  // left
    93103  if (lb > -COUENNE_INFINITY) { // if not unbounded
    94     cg -> createCut (cs, fn (lb) - lb, -1, w_ind, 1., x_ind, -1.); // up:  w - x <= f lb - lb
    95     cg -> createCut (cs, fn (lb) + lb, +1, w_ind, 1., x_ind,  1.); // dn:  w + x >= f lb + lb
     104    ncuts += cg -> createCut (cs, fn (lb) - lb, -1, w_ind, 1., x_ind, -1.); // up:  w-x <= f lb - lb
     105    ncuts += cg -> createCut (cs, fn (lb) + lb, +1, w_ind, 1., x_ind,  1.); // dn:  w+x >= f lb + lb
    96106  }
    97107
    98108  // right
    99109  if (ub <  COUENNE_INFINITY) { // if not unbounded
    100     cg -> createCut (cs, fn (ub) - ub, +1, w_ind, 1., x_ind, -1.); // dn: w - x >= f ub - ub
    101     cg -> createCut (cs, fn (ub) + ub, -1, w_ind, 1., x_ind,  1.); // up: w + x <= f ub + ub
    102   }
    103 
    104   delete lbe;
    105   delete ube;
     110    ncuts += cg -> createCut (cs, fn (ub) - ub, +1, w_ind, 1., x_ind, -1.); // dn: w - x >= f ub - ub
     111    ncuts += cg -> createCut (cs, fn (ub) + ub, -1, w_ind, 1., x_ind,  1.); // up: w + x <= f ub + ub
     112  }
     113
     114  return ncuts;
    106115}
    107116
     
    119128/// real linearization of sine/cosine
    120129
    121 void trigEnvelope (const CouenneCutGenerator *cg, // cut generator that has called us
     130int trigEnvelope (const CouenneCutGenerator *cg, // cut generator that has called us
    122131                   OsiCuts &cs,                   // cut set to be enriched
    123132                   exprAux *w,
     
    132141            ub = (*ube) (),
    133142            // if cosine, scale variables to pretend this is a sine problem
    134             displacement = (which_trig == COU_COSINE) ? M_PI_2 : 0.;
     143            displ = (which_trig == COU_COSINE) ? M_PI_2 : 0.;
    135144
    136145  delete lbe;
    137146  delete ube;
    138147
    139   int xi = arg -> Index (),
    140       wi = w   -> Index ();
     148  int ncuts = 0,
     149    xi = arg -> Index (),
     150    wi = w   -> Index ();
    141151
    142152  if (fabs (ub - lb) < COUENNE_EPS) {
     
    147157    else                        {f = cos (x0); fp = -sin (x0);}
    148158
    149     cg -> createCut (cs, f - fp*x0, 0, wi, 1., xi, -fp);
    150     return;
     159    return cg -> createCut (cs, f - fp*x0, 0, wi, 1., xi, -fp);
    151160  }
    152161
     
    156165       skip_dn = false;
    157166
    158   if (lb > -COUENNE_INFINITY) bayEnvelope (cg, cs, wi, xi, lb, ub, displacement, skip_up, skip_dn);
    159   if (ub <  COUENNE_INFINITY) bayEnvelope (cg, cs, wi, xi, ub, lb, displacement, skip_up, skip_dn);
     167  if (lb > -COUENNE_INFINITY) ncuts += bayEnvelope (cg, cs, wi, xi, lb, ub, displ, skip_up, skip_dn);
     168  if (ub <  COUENNE_INFINITY) ncuts += bayEnvelope (cg, cs, wi, xi, ub, lb, displ, skip_up, skip_dn);
     169
     170  return ncuts;
    160171}
    161172
     
    181192    sinrx0 = sin (rx0), zero;
    182193
    183   int
    184     //    up   = (modulo (rx0, 2*M_PI) < M_PI) ? +1 : -1,
     194  int ncuts = 0,
    185195    up   = (rx0 < M_PI) ? +1 : -1,
    186196    left = (x0  < x1)   ? +1 : -1;
     
    201211    // out of the "belly": tangent. If on upper bay we consider the
    202212    // lower half-plane, and viceversa --> use -up
    203     cg -> addTangent (cs, wi, xi, x0, sin (rx0), cos (rx0), -up);
     213    ncuts += cg -> addTangent (cs, wi, xi, x0, sin (rx0), cos (rx0), -up);
    204214
    205215    // leftmost extreme to search for tangent point
     
    211221                        (rx0, extr0, extr0 + M_PI_2))) <= 0)) { // before closest leaning point
    212222      if (!*s1) // -> chord, if not already added in previous call
    213         *s1 = (cg -> addSegment (cs, wi, xi, x0, sin (rx0), x1,         sin (rx1), up) > 0);
    214     } else     cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), up);
     223        *s1 = ((ncuts += cg -> addSegment (cs, wi, xi, x0, sin (rx0), x1,  sin (rx1), up)) > 0);
     224    } else      ncuts += cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), up);
    215225  }
    216226  else {
     
    224234      if (up * (sin (rx1) - sinrx0 - cosrx0 * (rx1-rx0)) < 0)
    225235        // (b,sinb) below tangent --> tangent
    226         cg -> addTangent (cs, wi, xi, x0, sinrx0, cosrx0, -up);
     236        ncuts += cg -> addTangent (cs, wi, xi, x0, sinrx0, cosrx0, -up);
    227237      else {    // up: either chord or leaning plane
    228238        CouNumber searchpt = M_PI_2 * (2 + 3*left - up);
     
    230240        if (left * (rx1 - tpt) < 0) {
    231241          if (!*s0)
    232             *s0 = cg -> addSegment (cs, wi, xi, x0, sin (rx0), x1,         sin (rx1), -up) > 0;
    233         } else    cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), -up);
     242            *s0 = ((ncuts += cg -> addSegment (cs, wi, xi, x0, sin (rx0), x1, sin (rx1), -up)) > 0);
     243        } else      ncuts += cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), -up);
    234244      }
    235245    } else {
    236246      CouNumber searchpt = M_PI_2 * (2 + 3*left - up);
    237247      tpt = trigNewton (rx0, searchpt, searchpt + left * M_PI_2);
    238       cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), -up);
     248      ncuts += cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), -up);
    239249    }
    240250
    241251    // down: other chord or leaning plane
    242252    if ((left * (rx1 - (zero + M_PI)) < 0) ||
    243         (left * (rx1 - (tpt = trigNewton (rx0,
    244                                           (2 +   left - up) * M_PI_2,
    245                                           (2 + 2*left - up) * M_PI_2))) < 0)) {
     253        (left * (rx1 - (tpt = trigNewton (rx0, (2 +   left - up) * M_PI_2,
     254                                               (2 + 2*left - up) * M_PI_2))) < 0)) {
    246255      if (!*s1)
    247         *s1 = (cg -> addSegment (cs, wi, xi, x0, sin (rx0), x1,         sin (rx1), up) > 0);
    248     } else     cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), up);
    249   }
    250 
    251   return 0;
    252 }
     256        *s1 = ((ncuts += cg -> addSegment (cs, wi, xi, x0, sin (rx0), x1, sin (rx1), up)) > 0);
     257    } else      ncuts += cg -> addSegment (cs, wi, xi, x0, sin (rx0), base + tpt, sin (tpt), up);
     258  }
     259
     260  return ncuts;
     261}
  • branches/Couenne/Couenne/src/problem/impliedBounds.cpp

    r539 r553  
    1717      nvar = nVars ();
    1818
    19   //printf ("implied========================\n");
     19  /*printf ("=====================implied\n");
     20  for (int i=0; i < nVars () + nAuxs (); i++)
     21    if ((i < nVars ()) || (Aux (i-nVars ()) -> Multiplicity () > 0))
     22    printf ("x%d: [%g,%g]\n", i, lb_ [i], ub_ [i]);*/
    2023
    21   /*  for (int i=0; i < nVars () + nAuxs (); i++)
    22       printf ("x%d: [%g,%g]\n", i, lb_ [i], ub_ [i]);*/
     24  for (int i=nAuxs (); i--;)
    2325
    24   for (int i=nAuxs (); i--;) {
     26    //    if (i+nVars () != 79) // !!!
     27{
    2528
    2629    //    for (int j=0; j<nAuxs () + nVars (); j++ )
     
    2831
    2932    if (lb_ [nvar+i] > ub_ [nvar+i] + COUENNE_EPS) {
    30       //      printf ("w_%d has infeasible bounds [%g,%g]\n",
    31       //        i+nvar, lb_ [nvar+i], ub_ [nvar+i]);
     33      /*      printf ("w_%d has infeasible bounds [%g,%g]\n",
     34              i+nvar, lb_ [nvar+i], ub_ [nvar+i]);*/
    3235      return -1;
    3336    }
     
    4043
    4144   if (auxiliaries_ [i] -> Image () -> impliedBound (nvar+i, lb_, ub_, chg_bds) > COUENNE_EPS) {
    42     //      printf ("implied: [%g,%g] -> [%g,%g]", l0, u0, lb_ [nvar+i], ub_ [nvar+i]);
     45     //      printf ("impli %2d [%g,%g] -> [%g,%g] ", nvar+i, l0, u0, lb_ [nvar+i], ub_ [nvar+i]);
     46     /*printf ("impli %2d ", nvar+i);
    4347
    44 
    45     //      auxiliaries_ [i] -> print (std::cout); printf (" := ");
    46     //      auxiliaries_ [i] -> Image () -> print (std::cout); printf ("\n");
     48      auxiliaries_ [i] -> print (std::cout); printf (" := ");
     49      auxiliaries_ [i] -> Image () -> print (std::cout); printf ("\n");*/
    4750      nchg++;
    4851    }
    4952  }
    50 
    51   /*  for (int i=0; i<nvar; i++)
     53  /*
     54  for (int i=0; i<nvar; i++)
    5255    printf (" [%g, %g]\n",
    5356            expression::Lbound (i),
     
    6467    auxiliaries_ [i] -> Image () -> print (std::cout); fflush (stdout);
    6568    printf ("\n");
    66     }*/
    67 
     69    }
     70  */
    6871  return nchg;
    6972}
  • branches/Couenne/Couenne/src/problem/problemIO.cpp

    r547 r553  
    4646    out << " ] " << std::endl;
    4747  }
     48
     49  printf ("bounds:\n");
     50  for (std::vector <exprVar *>::iterator i = variables_.begin ();
     51       i != variables_.end (); i++) {
     52    (*i) -> print (out);
     53    out << " in ["
     54        << lb_ [(*i) -> Index ()] << ','
     55        << ub_ [(*i) -> Index ()] << "]\n";
     56  }
     57
    4858  printf ("end\n");
    4959}
  • branches/Couenne/Couenne/src/problem/tightenBounds.cpp

    r539 r553  
    3030  /*printf ("tighten========================\n");
    3131
    32   for (int i=0; i<nVars(); i++)
    33     printf ("x_%d [%g, %g]\n", i,
    34             expression::Lbound (i),
    35             expression::Ubound (i));*/
     32  for (int i=0; i < nVars () + nAuxs (); i++)
     33    if ((i < nVars ()) || (Aux (i-nVars ()) -> Multiplicity () > 0))
     34      printf ("x_%d [%g, %g]\n", i,
     35              expression::Lbound (i),
     36              expression::Ubound (i));*/
    3637
    3738  for (register int i = nVars (), j=0;
    38        j < naux; j++) {
     39       j < naux; j++)
    3940
    40     //    for (int j=0; j<nAuxs () + nVars (); j++ )
    41     //      printf ("+++ %d: [%.4f %.4f]\n", j, lb_ [j], ub_ [j]);
     41    if (Aux (j) -> Multiplicity () > 0) {
    4242
    43     /*    printf ("w_%d [%g, %g] ", i+j,
     43      //    for (int j=0; j<nAuxs () + nVars (); j++ )
     44      //      printf ("+++ %d: [%.4f %.4f]\n", j, lb_ [j], ub_ [j]);
     45
     46      /*    printf ("w_%d [%g, %g] ", i+j,
    4447            expression::Lbound (i+j),
    4548            expression::Ubound (i+j));*/
    4649
    47     CouNumber ll = (*(Aux (j) -> Lb ())) (),
    48               uu = (*(Aux (j) -> Ub ())) ();
     50      CouNumber ll = (*(Aux (j) -> Lb ())) (),
     51        uu = (*(Aux (j) -> Ub ())) ();
    4952
    50     //    printf (" ---> [%g, %g]\n ", ll, uu);
     53      //    printf (" ---> [%g, %g]\n ", ll, uu);
    5154
    52     /*auxiliaries_ [j] -> print (std::cout);
    53     printf (" := ");
    54     auxiliaries_ [j] -> Image () -> print (std::cout); fflush (stdout);
    55     printf ("\n");*/
     55      /*auxiliaries_ [j] -> print (std::cout);
     56        printf (" := ");
     57        auxiliaries_ [j] -> Image () -> print (std::cout); fflush (stdout);
     58        printf ("\n");*/
    5659
    57     if (ll > uu + COUENNE_EPS) {
     60      if (ll > uu + COUENNE_EPS) {
     61        /*
     62        printf ("w_%d has infeasible bounds [%g,%g]: ", i+j, ll, uu);
     63        Aux (j) -> Lb () -> print (std::cout); printf (" --- ");
     64        Aux (j) -> Ub () -> print (std::cout); printf ("\n");
     65        */
     66        return -1; // declare this node infeasible
     67      }
    5868
    59       //      printf ("w_%d has infeasible bounds [%g,%g]: ", i+j, ll, uu);
    60       //      Aux (j) -> Lb () -> print (std::cout); printf (" - ");
    61       //      Aux (j) -> Ub () -> print (std::cout); printf ("\n");
     69      bool chg = false;
    6270
    63       return -1; // declare this node infeasible
     71      // check if lower bound got higher   
     72      if ((ll > - COUENNE_INFINITY) && (ll >= lb_ [i+j] + COUENNE_EPS)) {
     73
     74        /*printf ("update lbound %d: %.10e >= %.10e + %.12e\n",
     75          i+j, ll, lb_ [i+j], COUENNE_EPS);*/
     76
     77        /*printf ("update lbound %d: %g >= %g\n",
     78          i+j, ll, lb_ [i+j]);*/
     79
     80        /*printf ("propa %2d [%g,(%g)] -> [%g,(%g)] (%g)",
     81                i+j, lb_ [i+j], ub_ [i+j], ll, uu, lb_ [i+j] - ll);
     82        Aux (j)             -> print (std::cout); printf (" := ");
     83        Aux (j) -> Image () -> print (std::cout); printf ("\n");*/
     84
     85        lb_ [i+j] = ll;
     86        chg = true;
     87      }
     88
     89      // check if upper bound got lower
     90      if ((uu < COUENNE_INFINITY) && (uu <= ub_ [i+j] - COUENNE_EPS)) {
     91
     92        /*printf ("update ubound %d: %.10e <= %.10e - %.12e (%.12e)\n",
     93          i+j, uu, ub_ [i+j], COUENNE_EPS, uu - ub_ [i+j]);*/
     94        /*printf ("update ubound %d: %g >= %g\n",
     95          i+j, uu, ub_ [i+j]);*/
     96
     97        /*printf ("propa %2d [(%g),%g] -> [(%g),%g] (%g)",
     98                i+j, lb_ [i+j], ub_ [i+j], ll, uu, ub_ [i+j] - uu);
     99        Aux (j)             -> print (std::cout); printf (" := ");
     100        Aux (j) -> Image () -> print (std::cout); printf ("\n");*/
     101
     102        ub_ [i+j] = uu;
     103        chg = true;
     104      }
     105
     106      if (chg && chg_bds && !(chg_bds [i+j])) {
     107        nchg++;
     108        chg_bds [i+j] = 1;
     109      }
     110
     111      // useless if assume expression::[lu]b_ etc already point to
     112      // problem::[lu]b_
     113      expression::update (NULL, lb_, ub_);
    64114    }
    65 
    66     bool chg = false;
    67 
    68     // check if lower bound got higher   
    69     if ((ll > - COUENNE_INFINITY + 1) && (ll >= lb_ [i+j] + COUENNE_EPS)) {
    70 
    71       /*printf ("update lbound %d: %.10e >= %.10e + %.12e\n",
    72         i+j, ll, lb_ [i+j], COUENNE_EPS);*/
    73 
    74       /*printf ("update lbound %d: %g >= %g\n",
    75         i+j, ll, lb_ [i+j]);*/
    76 
    77       //      printf (" propaga: [%g,(%g)] -> [%g,(%g)] ", lb_ [i+j], ub_ [i+j], ll, uu);
    78       //      Aux (j)             -> print (std::cout); printf (" := ");
    79       //      Aux (j) -> Image () -> print (std::cout); printf ("\n");
    80 
    81       lb_ [i+j] = ll;
    82       chg = true;
    83     }
    84 
    85     // check if upper bound got lower
    86     if ((uu < COUENNE_INFINITY - 1) && (uu <= ub_ [i+j] - COUENNE_EPS)) {
    87 
    88       /*printf ("update ubound %d: %.10e <= %.10e - %.12e (%.12e)\n",
    89         i+j, uu, ub_ [i+j], COUENNE_EPS, uu - ub_ [i+j]);*/
    90       /*printf ("update ubound %d: %g >= %g\n",
    91         i+j, uu, ub_ [i+j]);*/
    92 
    93       //      printf (" propaga: [(%g),%g] -> [(%g),%g] ", lb_ [i+j], ub_ [i+j], ll, uu);
    94       //      Aux (j)             -> print (std::cout); printf (" := ");
    95       //      Aux (j) -> Image () -> print (std::cout); printf ("\n");
    96 
    97       ub_ [i+j] = uu;
    98       chg = true;
    99     }
    100 
    101     if (chg && chg_bds && !(chg_bds [i+j])) {
    102       nchg++;
    103       chg_bds [i+j] = 1;
    104     }
    105 
    106     // useless if assume expression::[lu]b_ etc already point to
    107     // problem::[lu]b_
    108     expression::update (NULL, lb_, ub_);
    109   }
    110115
    111116  return nchg;
Note: See TracChangeset for help on using the changeset viewer.