CouenneAggrProbing.cpp
Go to the documentation of this file.
1 /* $Id: CouenneAggrProbing.cpp 945 2013-04-06 20:25:21Z stefan $
2  *
3  * Name: CouenneAggrProbing.cpp
4  * Author: Giacomo Nannicini
5  * Purpose: Aggressive probing
6  *
7  * (C) Giacomo Nannicini, 2010.
8  * This file is licensed under the Eclipse Public License (EPL)
9  */
10 
11 #include "CouenneAggrProbing.hpp"
12 #include "CouenneProblemElem.hpp"
13 #include "CouenneExprVar.hpp"
14 #include "CouenneExprOpp.hpp"
15 //#include "BonCbc.hpp"
16 #include "CouenneBab.hpp"
17 #include "CouenneCutGenerator.hpp"
18 #include <string>
19 
20 #define COUENNE_AGGR_PROBING_FINITE_BOUND 1.0e+10
21 #define COUENNE_AGGR_PROBING_MIN_INTERVAL 1.0e-2
22 #define COUENNE_AGGR_PROBING_BND_RELAX COUENNE_EPS
23 
24 using namespace Couenne;
25 
28 {
29  couenne_ = setup;
31  maxTime_ = COIN_DBL_MAX;
32  maxFailedSteps_ = 10;
33  maxNodes_ = 1000;
35  restoreCutoff_ = false;
36 
37 }
38 
40  couenne_ = new CouenneSetup(*rhs.couenne_);
41  numCols_ = rhs.numCols_;
42  maxTime_ = rhs.maxTime_;
44  maxNodes_ = rhs.maxNodes_;
47 }
48 
50 }
51 
53  // Nothing for the moment, but will be added later as needed
54 }
55 
57  return maxTime_;
58 }
59 
60 void CouenneAggrProbing::setMaxTime(double value){
61  maxTime_ = value;
62 }
63 
65  return maxFailedSteps_;
66 }
67 
69  maxFailedSteps_ = value;
70 }
71 
73  return maxNodes_;
74 }
75 
77  maxNodes_ = value;
78 }
79 
81  return restoreCutoff_;
82 }
83 
85  restoreCutoff_ = value;
86 }
87 
88 double CouenneAggrProbing::probeVariable(int index, bool probeLower){
89 
90  // Useful objects for easy access
91  OsiSolverInterface* nlp = couenne_->nonlinearSolver();
92  OsiSolverInterface* lp = couenne_->continuousSolver();
93  CouenneProblem* problem = couenne_->couennePtr()->Problem();
94 
95  // Save initial bounds
96  double initUpper = lp->getColUpper()[index];
97  double initLower = lp->getColLower()[index];
98 
99  double* initLowerLp = new double[numCols_];
100  double* initUpperLp = new double[numCols_];
101 
102  memcpy(initLowerLp, lp->getColLower(), numCols_*sizeof(double));
103  memcpy(initUpperLp, lp->getColUpper(), numCols_*sizeof(double));
104 
105  if (initUpper < initLower + COUENNE_EPS){
106  // Variable is fixed, so we can't tighten
107  return ((probeLower) ? initLower : initUpper);
108  }
109 
110  // Index of the aux variable representing the objective function
111  int indobj = problem->Obj(0)->Body()->Index();
112 
113  // Initial cutoff value
114  double initCutoff = problem->getCutOff (); //Ub()[indobj];
115 
116  double* initCutoffSol = NULL;
117 
118  if (restoreCutoff_ && problem->getCutOff() < COUENNE_INFINITY){
119  initCutoffSol = new double[numCols_];
120  memcpy(initCutoffSol, problem->getCutOffSol(), numCols_*sizeof(double));
121  }
122 
123  // Save parameters
124  Bonmin::BabSetupBase::NodeComparison initNodeComparison =
130  //couenne_->nodeComparisonMethod() = Bonmin::BabSetupBase::bestBound;
133  problem->setCheckAuxBounds(true);
134 
137  couenne_->heuristics().clear();
138 
139  double currentBound = (probeLower) ? initLower : initUpper;
140  double startTime = CoinCpuTime();
141  int failedSteps = 0;
142  double intervalSize = 0.0;
143  double tryBound = 0.0;
144 
145  int iter = 0;
146 
147  if (probeLower)
148  std::cout << "Probing lower on var " << index << std::endl;
149  else
150  std::cout << "Probing upper on var " << index << std::endl;
151 
152  if ((fabs(currentBound) > COUENNE_AGGR_PROBING_FINITE_BOUND) &&
153  ((probeLower && initUpper > -COUENNE_AGGR_PROBING_FINITE_BOUND) ||
154  (!probeLower && initLower < COUENNE_AGGR_PROBING_FINITE_BOUND))){
155  // The bound is too large to apply the standard probing method;
156  // try to reduce it to a finite value. We only do this if we want
157  // to probe a variable on an infinite (or close to) bound, and the
158  // other bound of the variable is sufficiently far away
159  if (probeLower){
161  lp->setColLower(index, currentBound);
162  problem->Lb()[index] = currentBound;
163  lp->setColUpper(index, tryBound);
164  problem->Ub()[index] = tryBound;
165  if (index < problem->nOrigVars()){
166  nlp->setColLower(index, currentBound);
167  nlp->setColUpper(index, tryBound);
168  }
169  }
170  else{
172  lp->setColLower(index, tryBound);
173  problem->Lb()[index] = tryBound;
174  lp->setColUpper(index, currentBound);
175  problem->Ub()[index] = currentBound;
176  if (index < problem->nOrigVars()){
177  nlp->setColLower(index, tryBound);
178  nlp->setColUpper(index, currentBound);
179  }
180  }
181 
184  CoinMin(maxTime_-(CoinCpuTime()-startTime),
185  maxTime_*0.5));
186 
187  if (restoreCutoff_){
188  problem->resetCutOff(initCutoff);
189  if (indobj >= 0)
190  problem->Ub(indobj) = initCutoff;
191  problem->installCutOff();
192  }
193 
194  std::cout << "Iteration " << iter << ", current bound " << currentBound
195  << ", try bound " << tryBound << std::endl;
196 
198  // Bonmin::Bab bb;
199  // bb.setUsingCouenne(true);
200 
201  CouenneBab bb;
202 
203  bb(couenne_);
204  if (bb.model().isProvenInfeasible()){
206  currentBound = tryBound;
207  std::cout << "Probing succeeded; brought to finite" << std::endl;
208  }
209  else{
211  std::cout << "Probing failed; still infinity, exit" << std::endl;
212  }
213  iter++;
214  }
215 
216  // Now that we have a finite bound, pick size of the probing interval
217 
218  // Override (for testing - will be chosen automatically in final
219  // implementation)
220  intervalSize = 0.1;
221 
222  if (intervalSize < COUENNE_AGGR_PROBING_MIN_INTERVAL){
223  intervalSize = COUENNE_AGGR_PROBING_MIN_INTERVAL;
224  }
225 
226  while ((fabs(currentBound) <= COUENNE_AGGR_PROBING_FINITE_BOUND) &&
227  ((CoinCpuTime() - startTime) < maxTime_) &&
228  (failedSteps < maxFailedSteps_) &&
229  (intervalSize >= COUENNE_AGGR_PROBING_MIN_INTERVAL) &&
230  iter < 100){
231 
233  if (probeLower){
234  tryBound = currentBound + intervalSize;
235  if (tryBound > initUpper){
236  // It does not make sense to use bounds larger than the initial
237  // ones
238  tryBound = initUpper;
239  }
240  if (lp->isInteger(index)){
241  tryBound = floor(tryBound);
242  }
243  // Relax bounds a little bit
244  lp->setColLower(index, currentBound - COUENNE_AGGR_PROBING_BND_RELAX);
245  problem->Lb()[index] = currentBound - COUENNE_AGGR_PROBING_BND_RELAX;
246  lp->setColUpper(index, tryBound + COUENNE_AGGR_PROBING_BND_RELAX);
247  problem->Ub()[index] = tryBound + COUENNE_AGGR_PROBING_BND_RELAX;
248  if (index < problem->nOrigVars()){
249  nlp->setColLower(index, currentBound - COUENNE_AGGR_PROBING_BND_RELAX);
250  nlp->setColUpper(index, tryBound + COUENNE_AGGR_PROBING_BND_RELAX);
251  }
252  }
253  else{
254  tryBound = currentBound - intervalSize;
255  if (tryBound < initLower){
256  // It does not make sense to use bounds larger than the initial
257  // ones
258  tryBound = initLower;
259  }
260  if (lp->isInteger(index)){
261  tryBound = ceil(tryBound);
262  }
263  // Relax bounds a little bit
264  lp->setColLower(index, tryBound - COUENNE_AGGR_PROBING_BND_RELAX);
265  problem->Lb()[index] = tryBound - COUENNE_AGGR_PROBING_BND_RELAX;
266  lp->setColUpper(index, currentBound + COUENNE_AGGR_PROBING_BND_RELAX);
267  problem->Ub()[index] = currentBound + COUENNE_AGGR_PROBING_BND_RELAX;
268  if (index < problem->nOrigVars()){
269  nlp->setColLower(index, tryBound - COUENNE_AGGR_PROBING_BND_RELAX);
270  nlp->setColUpper(index, currentBound + COUENNE_AGGR_PROBING_BND_RELAX);
271  }
272  }
273 
274  lp->resolve();
275  problem->domain()->push(numCols_, lp->getColSolution(),
276  lp->getColLower(), lp->getColUpper());
277 
280  CoinMin(maxTime_-(CoinCpuTime()-startTime),
281  maxTime_*0.5));
282 
283  if (restoreCutoff_){
284  if (indobj >= 0)
285  problem->Ub(indobj) = initCutoff;
286  problem->resetCutOff(initCutoff);
287  problem->installCutOff();
288  }
289 
290  std::cout << "Iteration " << iter << ", current bound " << currentBound
291  << ", try bound " << tryBound << std::endl;
292 
294  // Bonmin::Bab bb;
295  // bb.setUsingCouenne(true);
296 
297  CouenneBab bb;
298  bb(couenne_);
299 
300  problem->domain()->pop();
301 
302  double obj = 0.0;
304  bool intervalSearched = (bb.model().isProvenOptimal() ||
305  bb.model().isProvenInfeasible());
306 
307  if ((!intervalSearched) || // If the search is not complete
308  (restoreCutoff_ && // or we don't want to accept new solutions
309  problem->getCutOffSol() && // and we have a new feasible solution
310  problem->checkNLP(problem->getCutOffSol(), obj, true))){
312  if (lp->isInteger(index) && fabs(tryBound-currentBound) < 0.5){
314  failedSteps = maxFailedSteps_;
315  }
316  else{
317  intervalSize /= 2;
318  }
319  failedSteps++;
320  std::cout << "Probing failed; shrinking interval" << std::endl;
321  }
322  else{
327  if (lp->isInteger(index) && fabs(tryBound-currentBound) < 0.5){
330  intervalSize = 1.0;
331  }
332  else{
333  intervalSize *= 2;
334  }
335  currentBound = tryBound;
336  if (lp->isInteger(index)){
337  if (probeLower){
338  currentBound += 1.0;
339  }
340  else {
341  currentBound -= 1.0;
342  }
343  }
344  failedSteps = 0;
345  std::cout << "Probing succeeded; enlarging interval" << std::endl;
346  }
347 
348  // Check early termination condition: if we manage to fix the
349  // variable (unlikely), there is nothing more we can do
350  if ((probeLower && fabs(currentBound-initUpper) < COUENNE_EPS) ||
351  (!probeLower && fabs(currentBound-initLower) < COUENNE_EPS)){
352  failedSteps = maxFailedSteps_;
353  }
354 
355  // Reset cutoff
356  if (restoreCutoff_){
357  if (indobj >= 0)
358  problem->Ub(indobj) = initCutoff;
359  problem->resetCutOff(initCutoff);
360  problem->installCutOff();
361  }
362 
363  problem->domain()->pop();
364 
365  iter++;
366  }
367 
370  lp->setColLower(initLowerLp);
371  lp->setColUpper(initUpperLp);
372  nlp->setColLower(initLowerLp);
373  nlp->setColUpper(initUpperLp);
374  memcpy(problem->Lb(), initLowerLp, numCols_*sizeof(double));
375  memcpy(problem->Ub(), initUpperLp, numCols_*sizeof(double));
376 
378  problem->setCheckAuxBounds(false);
379  //couenne_->nodeComparisonMethod() = initNodeComparison;
380  couenne_->setNodeComparisonMethod(initNodeComparison);
384  couenne_->heuristics() = heuristics;
385 
387  if (restoreCutoff_){
388  problem->resetCutOff();
389  problem->setCutOff(initCutoff, initCutoffSol);
390  if (initCutoffSol){
391  delete[] initCutoffSol;
392  }
393  }
394 
395  delete[] initLowerLp;
396  delete[] initUpperLp;
397 
399  return currentBound;
400 
401 }
402 
403 double CouenneAggrProbing::probeVariable2(int index, bool probeLower){
404  // Does not work! It's impossible to get Maximization problems working...
405  // Adding extra variables doesn't seem to do the trick
406 
407  // Useful objects for easy access
408  OsiSolverInterface* lp = couenne_->continuousSolver();
409  CouenneProblem* problem = couenne_->couennePtr()->Problem();
410 
411  // Save initial bounds
412  double initUpper = lp->getColUpper()[index];
413  double initLower = lp->getColLower()[index];
414 
415  if (initUpper < initLower + COUENNE_EPS){
416  // Variable is fixed, so we can't probe anything
417  return ((probeLower) ? initLower : initUpper);
418  }
419 
423  Bonmin::BabSetupBase::NodeComparison initNodeComparison =
431 
436  couenne_->heuristics().clear();
437 
439  double* initLpObj = new double[numCols_];
440  memcpy(initLpObj, lp->getObjCoefficients(), numCols_*sizeof(double));
441  expression* initProbObj = problem->Obj(0)->Body();
442 
443  double* newLpObj = new double[numCols_];
444  memset(newLpObj, 0, numCols_*sizeof(double));
445 
446  //expression* exprObj = NULL; // PB: unused
447  expression* extraVar = NULL;
448 
449  lp->writeLp("before");
450 
451  if (probeLower){
452  std::cout << "Probing LOWER" << std::endl;
453  // Set LP objective function
454  newLpObj[index] = 1.0;
455  lp->setObjective(newLpObj);
456 
457  lp->writeLp("lower");
458 
459  // Set CouenneProblem objective
460  problem->Obj(0)->Body(problem->Variables()[index]);
461  // couenne_->setDoubleParameter(Bonmin::BabSetupBase::Cutoff, COIN_DBL_MAX);
462  // problem->setCutOff(initUpper, lp->getColUpper());
463  // problem->installCutOff();
464 
465  }
466  else{
467  // We cannot maximize an objective function in Couenne; so, we
468  // have to introduce an additional variable, equal to the opposite
469  // of the variable that we want to maximize, and minimize that.
470 
471  // Add one column and one row to the LP
472  int extraCol = numCols_;
473  lp->setObjective(newLpObj);
474  lp->addCol(0, NULL, NULL, -initUpper, -initLower, 1.0);
475 
476  // Create the row x_{extraCol} = -x_{index}
477  int rowIndices[2] = {index, extraCol};
478  double rowElements[2] = {1.0, 1.0};
479  lp->addRow(2, rowIndices, rowElements, 0.0, 0.0);
480  lp->resolve();
481 
482  // Create the expression x_{extraCol} = -x_{index} in CouenneProblem
483  extraVar = problem->addVariable(lp->isInteger(index), NULL);
484  // exprObj = new exprOpp(problem->Variables()[index]->clone());
485  // problem->addEQConstraint(extraVar, exprObj);
486  problem->Obj(0)->Body(extraVar);
487 
488  // couenne_->setDoubleParameter(Bonmin::BabSetupBase::Cutoff, COIN_DBL_MAX);
489  // problem->setCutOff(-initLower, lp->getColLower());
490  // problem->installCutOff();
491 
492  lp->writeLp("upper");
493  }
494 
498  maxTime_);
499 
500  // Bonmin::Bab bb;
501  // bb.setUsingCouenne(true);
502 
503  CouenneBab bb;
504  bb(couenne_);
505 
506  double bestBound = bb.model().getBestPossibleObjValue();
507 
508  std::cout << "Obtained bound: " << bb.model().getBestPossibleObjValue() << std::endl;
509 
510 
512  couenne_->setNodeComparisonMethod (initNodeComparison);
516  couenne_->heuristics() = heuristics;
517 
518  if (!probeLower){
519  int extra = lp->getNumCols()-1;
520  lp->deleteCols(1, &extra);
521  extra = lp->getNumRows()-1;
522  lp->deleteRows(1, &extra);
523  problem->Variables().pop_back();
524  delete extraVar;
525  // delete exprObj;
526  }
527 
528  lp->setObjective(initLpObj);
529  problem->Obj(0)->Body(initProbObj);
530 
531  delete[] initLpObj;
532  delete[] newLpObj;
533 
534  return ((probeLower) ? bestBound : -bestBound);
535 
536 }
int getIntParameter(const IntParameter &p) const
Return value of integer parameter.
CouNumber & Ub(int i) const
upper bound on
int nVars() const
Total number of variables.
CouenneCutGenerator * couennePtr() const
return pointer to cut generator (used to get pointer to problem)
void setCutOff(CouNumber cutoff, const CouNumber *sol=NULL) const
Set cutoff.
Definition: problem.cpp:326
expression * Body() const
get body
Cut Generator for aggressive BT; i.e., an aggressive probing.
CouNumber & Lb(int i) const
lower bound on
CouenneSetup * couenne_
Pointer to the CouenneProblem representation.
void setMaxNodes(int value)
Set/get maximum number of nodes to probe one variable.
CouenneProblem * Problem() const
return pointer to symbolic problem
Domain * domain() const
return current point &amp; bounds
CouNumber getCutOff() const
Get cutoff.
void setMaxTime(double value)
Set/get maximum time to probe one variable.
bool restoreCutoff_
Restore initial cutoff (value and solution)?
CouNumber * getCutOffSol() const
Get cutoff solution.
void setDoubleParameter(const DoubleParameter &p, const double val)
modify parameter (used for MaxTime)
limit on number of integer feasible solution.
virtual int Index() const
Return index of variable (only valid for exprVar and exprAux)
void fint fint fint real fint real real real real real real real real real fint real fint * lp
void pop()
restore previous point
Definition: domain.cpp:255
static void registerOptions(Ipopt::SmartPtr< Bonmin::RegisteredOptions > roptions)
Add list of options to be read from file.
expression * addVariable(bool isint=false, Domain *d=NULL)
Add original variable.
int numCols_
Number of columns (want to have this handy)
Class for MINLP problems with symbolic information.
double initCutoff_
Initial cutoff.
#define COUENNE_AGGR_PROBING_BND_RELAX
#define COUENNE_EPS
std::vector< exprVar * > & Variables()
Return vector of variables (symbolic representation)
double probeVariable(int index, bool probeLower)
Probe one variable (try to tigthen the lower or the upper bound, depending on the value of the second...
int maxFailedSteps_
Maximum number of failed iterations.
#define COUENNE_AGGR_PROBING_MIN_INTERVAL
double probeVariable2(int index, bool lower)
Alternative probing algorithm.
void installCutOff() const
Make cutoff known to the problem.
Definition: problem.cpp:363
NodeComparison & nodeComparisonMethod()
Method used to compare nodes.
void setNodeComparisonMethod(Bonmin::BabSetupBase::NodeComparison c)
NodeComparison
Strategies for comparing the nodes on the heap.
void resetCutOff(CouNumber value=COUENNE_INFINITY) const
Reset cutoff.
Definition: problem.cpp:348
#define COUENNE_INFINITY
void push(int dim, CouNumber *x, CouNumber *lb, CouNumber *ub, bool copy=true)
save current point and start using another
Definition: domain.cpp:166
OsiTMINLPInterface * nonlinearSolver()
Pointer to the non-linear solver used.
OsiSolverInterface * continuousSolver()
Pointer to the continuous solver to use for relaxations.
#define COUENNE_AGGR_PROBING_FINITE_BOUND
void setIntParameter(const IntParameter &p, const int v)
Return value of integer parameter.
bool checkNLP(const double *solution, double &obj, bool recompute=false) const
Check if solution is MINLP feasible.
Definition: checkNLP.cpp:23
HeuristicMethods & heuristics()
list of Heuristic methods to use.
Expression base class.
void setMaxFailedSteps(int value)
Set/get maximum number of failed steps.
CouenneObjective * Obj(int i) const
i-th objective
void setCheckAuxBounds(bool value)
set the value for checkAuxBounds.
std::list< HeuristicMethod > HeuristicMethods
int maxNodes_
Maximum number of nodes in probing.
const CbcModel & model() const
Get cbc model used to solve.
Definition: BonCbc.hpp:87
void setRestoreCutoff(bool value)
Set/get restoreCutoff parameter (should we restore the initial cutoff value after each probing run...
double maxTime_
Maximum time to probe one variable.
double getDoubleParameter(const DoubleParameter &p) const
modify parameter (used for MaxTime)
CouenneAggrProbing(CouenneSetup *couenne, const Ipopt::SmartPtr< Ipopt::OptionsList > options)
Constructor.