source: stable/0.2/Couenne/src/branch/infeasibility.cpp @ 159

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

created new stable branch 0.2 from trunk (rev. 157)

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