source: trunk/Couenne/src/branch/infeasibility.cpp @ 121

Last change on this file since 121 was 121, checked in by pbelotti, 11 years ago

make a NaN infeasibility a nonzero infeasibility. Also, make a fixed variable CouenneVarObject::infeasibility a zero infeasibility (branch is useless and makes [down|up]estimate null).

File size: 5.4 KB
Line 
1/*
2 * Name:    infeasibility.cpp
3 * Authors: Pietro Belotti, Carnegie Mellon University
4 * Purpose: Compute infeasibility of a variable, looking at all expressions it appears in
5 *
6 * (C) Carnegie-Mellon University, 2008.
7 * This file is licensed under the Common Public License (CPL)
8 */
9
10#include "CoinHelperFunctions.hpp"
11
12#include "CouenneProblem.hpp"
13#include "CouenneVarObject.hpp"
14
15/// weights for computing infeasibility
16const CouNumber weiMin = 0.8; // for minimum of infeasibilities of nonlinear objects
17const CouNumber weiMax = 1.3; //     maximum
18const CouNumber weiSum = 0.1; //     sum
19const CouNumber weiAvg = 0.0; //     average
20
21
22/// compute infeasibility of this variable, |w - f(x)| (where w is
23/// the auxiliary variable defined as w = f(x))
24/// TODO: suggest way
25double CouenneVarObject::infeasibility (const OsiBranchingInformation *info, int &way) const {
26
27  assert (reference_);
28  int index = reference_ -> Index ();
29
30  /*printf ("===INFEAS %d = %g [%g,%g]\n",
31          index,
32          info -> solution_ [index],
33          info -> lower_ [index],
34          info -> upper_ [index]);*/
35
36  /*if (info -> lower_ [index] >=
37      info -> upper_ [index] - CoinMin (COUENNE_EPS, feas_tolerance_))
38    return (reference_ -> isInteger ()) ?
39    intInfeasibility (info -> solution_ [index]) : 0.;*/
40
41  problem_ -> domain () -> push
42    (problem_ -> nVars (),
43     info -> solution_, 
44     info -> lower_, 
45     info -> upper_);
46
47  const std::set <int> &dependence = problem_ -> Dependence () [index];
48
49  //////////////////////////////////////////////
50  CouNumber retval = checkInfeasibility (info);
51  //////////////////////////////////////////////
52
53  if (//(retval > CoinMin (COUENNE_EPS, feas_tolerance_)) &&
54      (jnlst_ -> ProduceOutput (J_DETAILED, J_BRANCHING))) {
55    printf ("infeasVar x%d %-10g [", reference_ -> Index (), retval);
56    reference_             -> print (); 
57    if ((dependence.size () == 0) && (reference_ -> Image ())) { // if no list, print image
58      printf (" := ");
59      reference_ -> Image () -> print ();
60    }
61    printf ("]\n");
62  }
63
64  // add term to stop branching on very tiny intervals
65
66  // Compute the up and down estimates here
67  // TODO: We probably only have to do that if pseudo costs option has
68  // been chosen
69
70  const CouenneObject *obj_ignored = NULL; // only care about it in CouenneVarObject.cpp
71
72  CouNumber brkPt = computeBranchingPoint (info, way, obj_ignored);
73
74  if (pseudoMultType_ != PROJECTDIST)
75    setEstimates (info, &retval, &brkPt);
76
77  if (jnlst_ -> ProduceOutput (J_DETAILED, J_BRANCHING)) {
78    printf("index = %d up = %e down = %e bounds [%e,%e] brpt = %e\n", 
79           index, upEstimate_, downEstimate_, 
80           info -> lower_ [index],
81           info -> upper_ [index], brkPt);
82  }
83
84  problem_ -> domain () -> pop (); // domain not used below
85
86  retval = ((retval < CoinMin (COUENNE_EPS, feas_tolerance_)) ? 0. : retval);
87
88  //printf ("infVar x%d ==> returning %g\n", reference_ -> Index (), (reference_ -> isInteger ()) ?
89  //CoinMax (retval, intInfeasibility (info -> solution_ [reference_ -> Index ()])) :
90  //retval);
91
92  return (reference_ -> isInteger ()) ? 
93    CoinMax (retval, intInfeasibility (info -> solution_ [reference_ -> Index ()])) :
94    retval;
95}
96
97
98/// compute infeasibility of this variable, |w - f(x)|, where w is
99/// the auxiliary variable defined as w = f(x)
100double CouenneVarObject::checkInfeasibility (const OsiBranchingInformation * info) const {
101
102  int index = reference_ -> Index ();
103
104  const std::set <int> &dependence = problem_ -> Dependence () [index];
105
106  if (dependence.size () == 0) { // this is a top level auxiliary,
107                                 // nowhere an independent
108
109    // check first if this variable is fixed. If so, it's useless to branch on it
110    if (fabs (info -> upper_ [index] - 
111              info -> lower_ [index]) / 
112        (1. + fabs (info -> solution_ [index])) < COUENNE_EPS)
113      return 0.;
114
115    // otherwise, return a nonzero infeasibility, if necessary. It
116    // might make sense to branch on it
117    const CouenneObject *obj = problem_ -> Objects () [reference_ -> Index ()];
118
119    double retval = (obj -> Reference ()) ? 
120      (1. - 1. / (1. + info -> upper_ [index] - info -> lower_ [index])) *
121      weiSum * obj ->checkInfeasibility (info) : 0.;
122
123    return retval;
124
125    /*return (reference_ -> isInteger ()) ?
126      CoinMax (retval, intInfeasibility (info -> solution_ [reference_ -> Index ()])) :
127      retval;*/
128
129  } else {
130
131    CouNumber
132      infsum = 0.,
133      infmax = 0.,
134      infmin = COIN_DBL_MAX;
135
136    for (std::set <int>::const_iterator i = dependence.begin ();
137         i != dependence.end (); ++i) {
138
139      // *i is the index of an auxiliary that depends on reference_
140
141      const CouenneObject *obj = problem_ -> Objects () [*i];
142      CouNumber infeas = (obj -> Reference ()) ? obj -> checkInfeasibility (info) : 0.;
143
144      if (infeas > infmax) infmax = infeas;
145      if (infeas < infmin) infmin = infeas;
146      infsum += infeas;
147    }
148
149    double retval = 
150      // neglect it if variable has small bound interval (chech
151      // x84=x83/x5 in csched1.nl)
152      (1. - 1. / (1. + info -> upper_ [index] - info -> lower_ [index])) *
153      // to consider maximum, minimum, and sum/avg of the infeasibilities
154      (weiSum * infsum + 
155       weiAvg * (infsum / CoinMax (1., (CouNumber) dependence.size ())) + 
156       weiMin * infmin + 
157       weiMax * infmax);
158
159    return retval;
160
161    /*return (reference_ -> isInteger ()) ?
162      CoinMax (retval, intInfeasibility (info -> solution_ [reference_ -> Index ()])) :
163      retval;*/
164  }
165}
Note: See TracBrowser for help on using the repository browser.