limp-cbc-0.3.2.0: cbits/coin/OsiChooseVariable.cpp
// Copyright (C) 2006, International Business Machines
// Corporation and others. All Rights Reserved.
// This code is licensed under the terms of the Eclipse Public License (EPL).
#include <string>
#include <cassert>
#include <cfloat>
#include <cmath>
#include "CoinPragma.hpp"
#include "OsiSolverInterface.hpp"
#include "OsiAuxInfo.hpp"
#include "OsiSolverBranch.hpp"
#include "CoinWarmStartBasis.hpp"
#include "CoinPackedMatrix.hpp"
#include "CoinTime.hpp"
#include "CoinSort.hpp"
#include "CoinFinite.hpp"
#include "OsiChooseVariable.hpp"
using namespace std;
OsiChooseVariable::OsiChooseVariable() :
goodObjectiveValue_(COIN_DBL_MAX),
upChange_(0.0),
downChange_(0.0),
goodSolution_(NULL),
list_(NULL),
useful_(NULL),
solver_(NULL),
status_(-1),
bestObjectIndex_(-1),
bestWhichWay_(-1),
firstForcedObjectIndex_(-1),
firstForcedWhichWay_(-1),
numberUnsatisfied_(0),
numberStrong_(0),
numberOnList_(0),
numberStrongDone_(0),
numberStrongIterations_(0),
numberStrongFixed_(0),
trustStrongForBound_(true),
trustStrongForSolution_(true)
{
}
OsiChooseVariable::OsiChooseVariable(const OsiSolverInterface * solver) :
goodObjectiveValue_(COIN_DBL_MAX),
upChange_(0.0),
downChange_(0.0),
goodSolution_(NULL),
solver_(solver),
status_(-1),
bestObjectIndex_(-1),
bestWhichWay_(-1),
firstForcedObjectIndex_(-1),
firstForcedWhichWay_(-1),
numberUnsatisfied_(0),
numberStrong_(0),
numberOnList_(0),
numberStrongDone_(0),
numberStrongIterations_(0),
numberStrongFixed_(0),
trustStrongForBound_(true),
trustStrongForSolution_(true)
{
// create useful arrays
int numberObjects = solver_->numberObjects();
list_ = new int [numberObjects];
useful_ = new double [numberObjects];
}
OsiChooseVariable::OsiChooseVariable(const OsiChooseVariable & rhs)
{
goodObjectiveValue_ = rhs.goodObjectiveValue_;
upChange_ = rhs.upChange_;
downChange_ = rhs.downChange_;
status_ = rhs.status_;
bestObjectIndex_ = rhs.bestObjectIndex_;
bestWhichWay_ = rhs.bestWhichWay_;
firstForcedObjectIndex_ = rhs.firstForcedObjectIndex_;
firstForcedWhichWay_ = rhs.firstForcedWhichWay_;
numberUnsatisfied_ = rhs.numberUnsatisfied_;
numberStrong_ = rhs.numberStrong_;
numberOnList_ = rhs.numberOnList_;
numberStrongDone_ = rhs.numberStrongDone_;
numberStrongIterations_ = rhs.numberStrongIterations_;
numberStrongFixed_ = rhs.numberStrongFixed_;
trustStrongForBound_ = rhs.trustStrongForBound_;
trustStrongForSolution_ = rhs.trustStrongForSolution_;
solver_ = rhs.solver_;
if (solver_) {
int numberObjects = solver_->numberObjects();
int numberColumns = solver_->getNumCols();
if (rhs.goodSolution_) {
goodSolution_ = CoinCopyOfArray(rhs.goodSolution_,numberColumns);
} else {
goodSolution_ = NULL;
}
list_ = CoinCopyOfArray(rhs.list_,numberObjects);
useful_ = CoinCopyOfArray(rhs.useful_,numberObjects);
} else {
goodSolution_ = NULL;
list_ = NULL;
useful_ = NULL;
}
}
OsiChooseVariable &
OsiChooseVariable::operator=(const OsiChooseVariable & rhs)
{
if (this != &rhs) {
delete [] goodSolution_;
delete [] list_;
delete [] useful_;
goodObjectiveValue_ = rhs.goodObjectiveValue_;
upChange_ = rhs.upChange_;
downChange_ = rhs.downChange_;
status_ = rhs.status_;
bestObjectIndex_ = rhs.bestObjectIndex_;
bestWhichWay_ = rhs.bestWhichWay_;
firstForcedObjectIndex_ = rhs.firstForcedObjectIndex_;
firstForcedWhichWay_ = rhs.firstForcedWhichWay_;
numberUnsatisfied_ = rhs.numberUnsatisfied_;
numberStrong_ = rhs.numberStrong_;
numberOnList_ = rhs.numberOnList_;
numberStrongDone_ = rhs.numberStrongDone_;
numberStrongIterations_ = rhs.numberStrongIterations_;
numberStrongFixed_ = rhs.numberStrongFixed_;
trustStrongForBound_ = rhs.trustStrongForBound_;
trustStrongForSolution_ = rhs.trustStrongForSolution_;
solver_ = rhs.solver_;
if (solver_) {
int numberObjects = solver_->numberObjects();
int numberColumns = solver_->getNumCols();
if (rhs.goodSolution_) {
goodSolution_ = CoinCopyOfArray(rhs.goodSolution_,numberColumns);
} else {
goodSolution_ = NULL;
}
list_ = CoinCopyOfArray(rhs.list_,numberObjects);
useful_ = CoinCopyOfArray(rhs.useful_,numberObjects);
} else {
goodSolution_ = NULL;
list_ = NULL;
useful_ = NULL;
}
}
return *this;
}
OsiChooseVariable::~OsiChooseVariable ()
{
delete [] goodSolution_;
delete [] list_;
delete [] useful_;
}
// Clone
OsiChooseVariable *
OsiChooseVariable::clone() const
{
return new OsiChooseVariable(*this);
}
// Set solver and redo arrays
void
OsiChooseVariable::setSolver (const OsiSolverInterface * solver)
{
solver_ = solver;
delete [] list_;
delete [] useful_;
// create useful arrays
int numberObjects = solver_->numberObjects();
list_ = new int [numberObjects];
useful_ = new double [numberObjects];
}
// Initialize
int
OsiChooseVariable::setupList ( OsiBranchingInformation *info, bool initialize)
{
if (initialize) {
status_=-2;
delete [] goodSolution_;
bestObjectIndex_=-1;
numberStrongDone_=0;
numberStrongIterations_ = 0;
numberStrongFixed_ = 0;
goodSolution_ = NULL;
goodObjectiveValue_ = COIN_DBL_MAX;
}
numberOnList_=0;
numberUnsatisfied_=0;
int numberObjects = solver_->numberObjects();
assert (numberObjects);
double check = 0.0;
int checkIndex=0;
int bestPriority=COIN_INT_MAX;
// pretend one strong even if none
int maximumStrong= numberStrong_ ? CoinMin(numberStrong_,numberObjects) : 1;
int putOther = numberObjects;
int i;
for (i=0;i<maximumStrong;i++) {
list_[i]=-1;
useful_[i]=0.0;
}
OsiObject ** object = info->solver_->objects();
// Say feasible
bool feasible = true;
for ( i=0;i<numberObjects;i++) {
int way;
double value = object[i]->infeasibility(info,way);
if (value>0.0) {
numberUnsatisfied_++;
if (value==COIN_DBL_MAX) {
// infeasible
feasible=false;
break;
}
int priorityLevel = object[i]->priority();
// Better priority? Flush choices.
if (priorityLevel<bestPriority) {
for (int j=0;j<maximumStrong;j++) {
if (list_[j]>=0) {
int iObject = list_[j];
list_[j]=-1;
useful_[j]=0.0;
list_[--putOther]=iObject;
}
}
bestPriority = priorityLevel;
check=0.0;
}
if (priorityLevel==bestPriority) {
if (value>check) {
//add to list
int iObject = list_[checkIndex];
if (iObject>=0)
list_[--putOther]=iObject; // to end
list_[checkIndex]=i;
useful_[checkIndex]=value;
// find worst
check=COIN_DBL_MAX;
for (int j=0;j<maximumStrong;j++) {
if (list_[j]>=0) {
if (useful_[j]<check) {
check=useful_[j];
checkIndex=j;
}
} else {
check=0.0;
checkIndex = j;
break;
}
}
} else {
// to end
list_[--putOther]=i;
}
} else {
// to end
list_[--putOther]=i;
}
}
}
// Get list
numberOnList_=0;
if (feasible) {
for (i=0;i<maximumStrong;i++) {
if (list_[i]>=0) {
list_[numberOnList_]=list_[i];
useful_[numberOnList_++]=-useful_[i];
}
}
if (numberOnList_) {
// Sort
CoinSort_2(useful_,useful_+numberOnList_,list_);
// move others
i = numberOnList_;
for (;putOther<numberObjects;putOther++)
list_[i++]=list_[putOther];
assert (i==numberUnsatisfied_);
if (!numberStrong_)
numberOnList_=0;
}
} else {
// not feasible
numberUnsatisfied_=-1;
}
return numberUnsatisfied_;
}
/* Choose a variable
Returns -
-1 Node is infeasible
0 Normal termination - we have a candidate
1 All looks satisfied - no candidate
2 We can change the bound on a variable - but we also have a strong branching candidate
3 We can change the bound on a variable - but we have a non-strong branching candidate
4 We can change the bound on a variable - no other candidates
We can pick up branch from whichObject() and whichWay()
We can pick up a forced branch (can change bound) from whichForcedObject() and whichForcedWay()
If we have a solution then we can pick up from goodObjectiveValue() and goodSolution()
*/
int
OsiChooseVariable::chooseVariable( OsiSolverInterface * solver, OsiBranchingInformation *, bool )
{
if (numberUnsatisfied_) {
bestObjectIndex_=list_[0];
bestWhichWay_ = solver->object(bestObjectIndex_)->whichWay();
firstForcedObjectIndex_ = -1;
firstForcedWhichWay_ =-1;
return 0;
} else {
return 1;
}
}
// Returns true if solution looks feasible against given objects
bool
OsiChooseVariable::feasibleSolution(const OsiBranchingInformation * info,
const double * solution,
int numberObjects,
const OsiObject ** objects)
{
bool satisfied=true;
const double * saveSolution = info->solution_;
info->solution_ = solution;
for (int i=0;i<numberObjects;i++) {
double value = objects[i]->checkInfeasibility(info);
if (value>0.0) {
satisfied=false;
break;
}
}
info->solution_ = saveSolution;
return satisfied;
}
// Saves a good solution
void
OsiChooseVariable::saveSolution(const OsiSolverInterface * solver)
{
delete [] goodSolution_;
int numberColumns = solver->getNumCols();
goodSolution_ = CoinCopyOfArray(solver->getColSolution(),numberColumns);
goodObjectiveValue_ = solver->getObjSense()*solver->getObjValue();
}
// Clears out good solution after use
void
OsiChooseVariable::clearGoodSolution()
{
delete [] goodSolution_;
goodSolution_ = NULL;
goodObjectiveValue_ = COIN_DBL_MAX;
}
/* This is a utility function which does strong branching on
a list of objects and stores the results in OsiHotInfo.objects.
On entry the object sequence is stored in the OsiHotInfo object
and maybe more.
It returns -
-1 - one branch was infeasible both ways
0 - all inspected - nothing can be fixed
1 - all inspected - some can be fixed (returnCriterion==0)
2 - may be returning early - one can be fixed (last one done) (returnCriterion==1)
3 - returning because max time
*/
int
OsiChooseStrong::doStrongBranching( OsiSolverInterface * solver,
OsiBranchingInformation *info,
int numberToDo, int returnCriterion)
{
// Might be faster to extend branch() to return bounds changed
double * saveLower = NULL;
double * saveUpper = NULL;
int numberColumns = solver->getNumCols();
solver->markHotStart();
const double * lower = info->lower_;
const double * upper = info->upper_;
saveLower = CoinCopyOfArray(info->lower_,numberColumns);
saveUpper = CoinCopyOfArray(info->upper_,numberColumns);
numResults_=0;
int returnCode=0;
double timeStart = CoinCpuTime();
for (int iDo=0;iDo<numberToDo;iDo++) {
OsiHotInfo * result = results_ + iDo;
// For now just 2 way
OsiBranchingObject * branch = result->branchingObject();
assert (branch->numberBranches()==2);
/*
Try the first direction. Each subsequent call to branch() performs the
specified branch and advances the branch object state to the next branch
alternative.)
*/
OsiSolverInterface * thisSolver = solver;
if (branch->boundBranch()) {
// ordinary
branch->branch(solver);
// maybe we should check bounds for stupidities here?
solver->solveFromHotStart() ;
} else {
// adding cuts or something
thisSolver = solver->clone();
branch->branch(thisSolver);
// set hot start iterations
int limit;
thisSolver->getIntParam(OsiMaxNumIterationHotStart,limit);
thisSolver->setIntParam(OsiMaxNumIteration,limit);
thisSolver->resolve();
}
// can check if we got solution
// status is 0 finished, 1 infeasible and 2 unfinished and 3 is solution
int status0 = result->updateInformation(thisSolver,info,this);
numberStrongIterations_ += thisSolver->getIterationCount();
if (status0==3) {
// new solution already saved
if (trustStrongForSolution_) {
info->cutoff_ = goodObjectiveValue_;
status0=0;
}
}
if (solver!=thisSolver)
delete thisSolver;
// Restore bounds
for (int j=0;j<numberColumns;j++) {
if (saveLower[j] != lower[j])
solver->setColLower(j,saveLower[j]);
if (saveUpper[j] != upper[j])
solver->setColUpper(j,saveUpper[j]);
}
/*
Try the next direction
*/
thisSolver = solver;
if (branch->boundBranch()) {
// ordinary
branch->branch(solver);
// maybe we should check bounds for stupidities here?
solver->solveFromHotStart() ;
} else {
// adding cuts or something
thisSolver = solver->clone();
branch->branch(thisSolver);
// set hot start iterations
int limit;
thisSolver->getIntParam(OsiMaxNumIterationHotStart,limit);
thisSolver->setIntParam(OsiMaxNumIteration,limit);
thisSolver->resolve();
}
// can check if we got solution
// status is 0 finished, 1 infeasible and 2 unfinished and 3 is solution
int status1 = result->updateInformation(thisSolver,info,this);
numberStrongDone_++;
numberStrongIterations_ += thisSolver->getIterationCount();
if (status1==3) {
// new solution already saved
if (trustStrongForSolution_) {
info->cutoff_ = goodObjectiveValue_;
status1=0;
}
}
if (solver!=thisSolver)
delete thisSolver;
// Restore bounds
for (int j=0;j<numberColumns;j++) {
if (saveLower[j] != lower[j])
solver->setColLower(j,saveLower[j]);
if (saveUpper[j] != upper[j])
solver->setColUpper(j,saveUpper[j]);
}
/*
End of evaluation for this candidate variable. Possibilities are:
* Both sides below cutoff; this variable is a candidate for branching.
* Both sides infeasible or above the objective cutoff: no further action
here. Break from the evaluation loop and assume the node will be purged
by the caller.
* One side below cutoff: Install the branch (i.e., fix the variable). Possibly break
from the evaluation loop and assume the node will be reoptimised by the
caller.
*/
numResults_++;
if (status0==1&&status1==1) {
// infeasible
returnCode=-1;
break; // exit loop
} else if (status0==1||status1==1) {
numberStrongFixed_++;
if (!returnCriterion) {
returnCode=1;
} else {
returnCode=2;
break;
}
}
bool hitMaxTime = ( CoinCpuTime()-timeStart > info->timeRemaining_);
if (hitMaxTime) {
returnCode=3;
break;
}
}
delete [] saveLower;
delete [] saveUpper;
// Delete the snapshot
solver->unmarkHotStart();
return returnCode;
}
// Given a candidate fill in useful information e.g. estimates
void
OsiChooseVariable::updateInformation(const OsiBranchingInformation *info,
int , OsiHotInfo * hotInfo)
{
int index = hotInfo->whichObject();
assert (index<solver_->numberObjects());
//assert (branch<2);
OsiObject ** object = info->solver_->objects();
upChange_ = object[index]->upEstimate();
downChange_ = object[index]->downEstimate();
}
#if 1
// Given a branch fill in useful information e.g. estimates
void
OsiChooseVariable::updateInformation( int index, int branch,
double , double ,
int )
{
assert (index<solver_->numberObjects());
assert (branch<2);
OsiObject ** object = solver_->objects();
if (branch)
upChange_ = object[index]->upEstimate();
else
downChange_ = object[index]->downEstimate();
}
#endif
//##############################################################################
void
OsiPseudoCosts::gutsOfDelete()
{
if (numberObjects_ > 0) {
numberObjects_ = 0;
numberBeforeTrusted_ = 0;
delete[] upTotalChange_; upTotalChange_ = NULL;
delete[] downTotalChange_; downTotalChange_ = NULL;
delete[] upNumber_; upNumber_ = NULL;
delete[] downNumber_; downNumber_ = NULL;
}
}
void
OsiPseudoCosts::gutsOfCopy(const OsiPseudoCosts& rhs)
{
numberObjects_ = rhs.numberObjects_;
numberBeforeTrusted_ = rhs.numberBeforeTrusted_;
if (numberObjects_ > 0) {
upTotalChange_ = CoinCopyOfArray(rhs.upTotalChange_,numberObjects_);
downTotalChange_ = CoinCopyOfArray(rhs.downTotalChange_,numberObjects_);
upNumber_ = CoinCopyOfArray(rhs.upNumber_,numberObjects_);
downNumber_ = CoinCopyOfArray(rhs.downNumber_,numberObjects_);
}
}
OsiPseudoCosts::OsiPseudoCosts() :
upTotalChange_(NULL),
downTotalChange_(NULL),
upNumber_(NULL),
downNumber_(NULL),
numberObjects_(0),
numberBeforeTrusted_(0)
{
}
OsiPseudoCosts::~OsiPseudoCosts()
{
gutsOfDelete();
}
OsiPseudoCosts::OsiPseudoCosts(const OsiPseudoCosts& rhs) :
upTotalChange_(NULL),
downTotalChange_(NULL),
upNumber_(NULL),
downNumber_(NULL),
numberObjects_(0),
numberBeforeTrusted_(0)
{
gutsOfCopy(rhs);
}
OsiPseudoCosts&
OsiPseudoCosts::operator=(const OsiPseudoCosts& rhs)
{
if (this != &rhs) {
gutsOfDelete();
gutsOfCopy(rhs);
}
return *this;
}
void
OsiPseudoCosts::initialize(int n)
{
gutsOfDelete();
numberObjects_ = n;
if (numberObjects_ > 0) {
upTotalChange_ = new double [numberObjects_];
downTotalChange_ = new double [numberObjects_];
upNumber_ = new int [numberObjects_];
downNumber_ = new int [numberObjects_];
CoinZeroN(upTotalChange_,numberObjects_);
CoinZeroN(downTotalChange_,numberObjects_);
CoinZeroN(upNumber_,numberObjects_);
CoinZeroN(downNumber_,numberObjects_);
}
}
//##############################################################################
OsiChooseStrong::OsiChooseStrong() :
OsiChooseVariable(),
shadowPriceMode_(0),
pseudoCosts_(),
results_(NULL),
numResults_(0)
{
}
OsiChooseStrong::OsiChooseStrong(const OsiSolverInterface * solver) :
OsiChooseVariable(solver),
shadowPriceMode_(0),
pseudoCosts_(),
results_(NULL),
numResults_(0)
{
// create useful arrays
pseudoCosts_.initialize(solver_->numberObjects());
}
OsiChooseStrong::OsiChooseStrong(const OsiChooseStrong & rhs) :
OsiChooseVariable(rhs),
shadowPriceMode_(rhs.shadowPriceMode_),
pseudoCosts_(rhs.pseudoCosts_),
results_(NULL),
numResults_(0)
{
}
OsiChooseStrong &
OsiChooseStrong::operator=(const OsiChooseStrong & rhs)
{
if (this != &rhs) {
OsiChooseVariable::operator=(rhs);
shadowPriceMode_ = rhs.shadowPriceMode_;
pseudoCosts_ = rhs.pseudoCosts_;
delete[] results_;
results_ = NULL;
numResults_ = 0;
}
return *this;
}
OsiChooseStrong::~OsiChooseStrong ()
{
delete[] results_;
}
// Clone
OsiChooseVariable *
OsiChooseStrong::clone() const
{
return new OsiChooseStrong(*this);
}
#define MAXMIN_CRITERION 0.85
// Initialize
int
OsiChooseStrong::setupList ( OsiBranchingInformation *info, bool initialize)
{
if (initialize) {
status_=-2;
delete [] goodSolution_;
bestObjectIndex_=-1;
numberStrongDone_=0;
numberStrongIterations_ = 0;
numberStrongFixed_ = 0;
goodSolution_ = NULL;
goodObjectiveValue_ = COIN_DBL_MAX;
}
numberOnList_=0;
numberUnsatisfied_=0;
int numberObjects = solver_->numberObjects();
if (numberObjects>pseudoCosts_.numberObjects()) {
// redo useful arrays
pseudoCosts_.initialize(numberObjects);
}
double check = -COIN_DBL_MAX;
int checkIndex=0;
int bestPriority=COIN_INT_MAX;
int maximumStrong= CoinMin(numberStrong_,numberObjects) ;
int putOther = numberObjects;
int i;
for (i=0;i<numberObjects;i++) {
list_[i]=-1;
useful_[i]=0.0;
}
OsiObject ** object = info->solver_->objects();
// Get average pseudo costs and see if pseudo shadow prices possible
int shadowPossible=shadowPriceMode_;
if (shadowPossible) {
for ( i=0;i<numberObjects;i++) {
if ( !object[i]->canHandleShadowPrices()) {
shadowPossible=0;
break;
}
}
if (shadowPossible) {
int numberRows = solver_->getNumRows();
const double * pi = info->pi_;
double sumPi=0.0;
for (i=0;i<numberRows;i++)
sumPi += fabs(pi[i]);
sumPi /= static_cast<double> (numberRows);
// and scale back
sumPi *= 0.01;
info->defaultDual_ = sumPi; // switch on
int numberColumns = solver_->getNumCols();
int size = CoinMax(numberColumns,2*numberRows);
info->usefulRegion_ = new double [size];
CoinZeroN(info->usefulRegion_,size);
info->indexRegion_ = new int [size];
}
}
double sumUp=0.0;
double numberUp=0.0;
double sumDown=0.0;
double numberDown=0.0;
const double* upTotalChange = pseudoCosts_.upTotalChange();
const double* downTotalChange = pseudoCosts_.downTotalChange();
const int* upNumber = pseudoCosts_.upNumber();
const int* downNumber = pseudoCosts_.downNumber();
const int numberBeforeTrusted = pseudoCosts_.numberBeforeTrusted();
for ( i=0;i<numberObjects;i++) {
sumUp += upTotalChange[i];
numberUp += upNumber[i];
sumDown += downTotalChange[i];
numberDown += downNumber[i];
}
double upMultiplier=(1.0+sumUp)/(1.0+numberUp);
double downMultiplier=(1.0+sumDown)/(1.0+numberDown);
// Say feasible
bool feasible = true;
#if 0
int pri[]={10,1000,10000};
int priCount[]={0,0,0};
#endif
for ( i=0;i<numberObjects;i++) {
int way;
double value = object[i]->infeasibility(info,way);
if (value>0.0) {
numberUnsatisfied_++;
if (value==COIN_DBL_MAX) {
// infeasible
feasible=false;
break;
}
int priorityLevel = object[i]->priority();
#if 0
for (int k=0;k<3;k++) {
if (priorityLevel==pri[k])
priCount[k]++;
}
#endif
// Better priority? Flush choices.
if (priorityLevel<bestPriority) {
for (int j=maximumStrong-1;j>=0;j--) {
if (list_[j]>=0) {
int iObject = list_[j];
list_[j]=-1;
useful_[j]=0.0;
list_[--putOther]=iObject;
}
}
maximumStrong = CoinMin(maximumStrong,putOther);
bestPriority = priorityLevel;
check=-COIN_DBL_MAX;
checkIndex=0;
}
if (priorityLevel==bestPriority) {
// Modify value
sumUp = upTotalChange[i]+1.0e-30;
numberUp = upNumber[i];
sumDown = downTotalChange[i]+1.0e-30;
numberDown = downNumber[i];
double upEstimate = object[i]->upEstimate();
double downEstimate = object[i]->downEstimate();
if (shadowPossible<2) {
upEstimate = numberUp ? ((upEstimate*sumUp)/numberUp) : (upEstimate*upMultiplier);
if (numberUp<numberBeforeTrusted)
upEstimate *= (numberBeforeTrusted+1.0)/(numberUp+1.0);
downEstimate = numberDown ? ((downEstimate*sumDown)/numberDown) : (downEstimate*downMultiplier);
if (numberDown<numberBeforeTrusted)
downEstimate *= (numberBeforeTrusted+1.0)/(numberDown+1.0);
} else {
// use shadow prices always
}
value = MAXMIN_CRITERION*CoinMin(upEstimate,downEstimate) + (1.0-MAXMIN_CRITERION)*CoinMax(upEstimate,downEstimate);
if (value>check) {
//add to list
int iObject = list_[checkIndex];
if (iObject>=0) {
assert (list_[putOther-1]<0);
list_[--putOther]=iObject; // to end
}
list_[checkIndex]=i;
assert (checkIndex<putOther);
useful_[checkIndex]=value;
// find worst
check=COIN_DBL_MAX;
maximumStrong = CoinMin(maximumStrong,putOther);
for (int j=0;j<maximumStrong;j++) {
if (list_[j]>=0) {
if (useful_[j]<check) {
check=useful_[j];
checkIndex=j;
}
} else {
check=0.0;
checkIndex = j;
break;
}
}
} else {
// to end
assert (list_[putOther-1]<0);
list_[--putOther]=i;
maximumStrong = CoinMin(maximumStrong,putOther);
}
} else {
// worse priority
// to end
assert (list_[putOther-1]<0);
list_[--putOther]=i;
maximumStrong = CoinMin(maximumStrong,putOther);
}
}
}
#if 0
printf("%d at %d, %d at %d and %d at %d\n",priCount[0],pri[0],
priCount[1],pri[1],priCount[2],pri[2]);
#endif
// Get list
numberOnList_=0;
if (feasible) {
for (i=0;i<CoinMin(maximumStrong,putOther);i++) {
if (list_[i]>=0) {
list_[numberOnList_]=list_[i];
useful_[numberOnList_++]=-useful_[i];
}
}
if (numberOnList_) {
// Sort
CoinSort_2(useful_,useful_+numberOnList_,list_);
// move others
i = numberOnList_;
for (;putOther<numberObjects;putOther++)
list_[i++]=list_[putOther];
assert (i==numberUnsatisfied_);
if (!numberStrong_)
numberOnList_=0;
}
} else {
// not feasible
numberUnsatisfied_=-1;
}
// Get rid of any shadow prices info
info->defaultDual_ = -1.0; // switch off
delete [] info->usefulRegion_;
delete [] info->indexRegion_;
return numberUnsatisfied_;
}
void
OsiChooseStrong::resetResults(int num)
{
delete[] results_;
numResults_ = 0;
results_ = new OsiHotInfo[num];
}
/* Choose a variable
Returns -
-1 Node is infeasible
0 Normal termination - we have a candidate
1 All looks satisfied - no candidate
2 We can change the bound on a variable - but we also have a strong branching candidate
3 We can change the bound on a variable - but we have a non-strong branching candidate
4 We can change the bound on a variable - no other candidates
We can pick up branch from whichObject() and whichWay()
We can pick up a forced branch (can change bound) from whichForcedObject() and whichForcedWay()
If we have a solution then we can pick up from goodObjectiveValue() and goodSolution()
*/
int
OsiChooseStrong::chooseVariable( OsiSolverInterface * solver, OsiBranchingInformation *info, bool fixVariables)
{
if (numberUnsatisfied_) {
const double* upTotalChange = pseudoCosts_.upTotalChange();
const double* downTotalChange = pseudoCosts_.downTotalChange();
const int* upNumber = pseudoCosts_.upNumber();
const int* downNumber = pseudoCosts_.downNumber();
int numberBeforeTrusted = pseudoCosts_.numberBeforeTrusted();
// Somehow we can get here with it 0 !
if (!numberBeforeTrusted) {
numberBeforeTrusted=5;
pseudoCosts_.setNumberBeforeTrusted(numberBeforeTrusted);
}
int numberLeft = CoinMin(numberStrong_-numberStrongDone_,numberUnsatisfied_);
int numberToDo=0;
resetResults(numberLeft);
int returnCode=0;
bestObjectIndex_ = -1;
bestWhichWay_ = -1;
firstForcedObjectIndex_ = -1;
firstForcedWhichWay_ =-1;
double bestTrusted=-COIN_DBL_MAX;
for (int i=0;i<numberLeft;i++) {
int iObject = list_[i];
if (upNumber[iObject]<numberBeforeTrusted||downNumber[iObject]<numberBeforeTrusted) {
results_[numberToDo++] = OsiHotInfo(solver, info,
solver->objects(), iObject);
} else {
const OsiObject * obj = solver->object(iObject);
double upEstimate = (upTotalChange[iObject]*obj->upEstimate())/upNumber[iObject];
double downEstimate = (downTotalChange[iObject]*obj->downEstimate())/downNumber[iObject];
double value = MAXMIN_CRITERION*CoinMin(upEstimate,downEstimate) + (1.0-MAXMIN_CRITERION)*CoinMax(upEstimate,downEstimate);
if (value > bestTrusted) {
bestObjectIndex_=iObject;
bestWhichWay_ = upEstimate>downEstimate ? 0 : 1;
bestTrusted = value;
}
}
}
int numberFixed=0;
if (numberToDo) {
returnCode = doStrongBranching(solver, info, numberToDo, 1);
if (returnCode>=0&&returnCode<=2) {
if (returnCode) {
returnCode=4;
if (bestObjectIndex_>=0)
returnCode=3;
}
for (int i=0;i<numResults_;i++) {
int iObject = results_[i].whichObject();
double upEstimate;
if (results_[i].upStatus()!=1) {
assert (results_[i].upStatus()>=0);
upEstimate = results_[i].upChange();
} else {
// infeasible - just say expensive
if (info->cutoff_<1.0e50)
upEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
else
upEstimate = 2.0*fabs(info->objectiveValue_);
if (firstForcedObjectIndex_ <0) {
firstForcedObjectIndex_ = iObject;
firstForcedWhichWay_ =0;
}
numberFixed++;
if (fixVariables) {
const OsiObject * obj = solver->object(iObject);
OsiBranchingObject * branch = obj->createBranch(solver,info,0);
branch->branch(solver);
delete branch;
}
}
double downEstimate;
if (results_[i].downStatus()!=1) {
assert (results_[i].downStatus()>=0);
downEstimate = results_[i].downChange();
} else {
// infeasible - just say expensive
if (info->cutoff_<1.0e50)
downEstimate = 2.0*(info->cutoff_-info->objectiveValue_);
else
downEstimate = 2.0*fabs(info->objectiveValue_);
if (firstForcedObjectIndex_ <0) {
firstForcedObjectIndex_ = iObject;
firstForcedWhichWay_ =1;
}
numberFixed++;
if (fixVariables) {
const OsiObject * obj = solver->object(iObject);
OsiBranchingObject * branch = obj->createBranch(solver,info,1);
branch->branch(solver);
delete branch;
}
}
double value = MAXMIN_CRITERION*CoinMin(upEstimate,downEstimate) + (1.0-MAXMIN_CRITERION)*CoinMax(upEstimate,downEstimate);
if (value>bestTrusted) {
bestTrusted = value;
bestObjectIndex_ = iObject;
bestWhichWay_ = upEstimate>downEstimate ? 0 : 1;
// but override if there is a preferred way
const OsiObject * obj = solver->object(iObject);
if (obj->preferredWay()>=0&&obj->infeasibility())
bestWhichWay_ = obj->preferredWay();
if (returnCode)
returnCode=2;
}
}
} else if (returnCode==3) {
// max time - just choose one
bestObjectIndex_ = list_[0];
bestWhichWay_ = 0;
returnCode=0;
}
} else {
bestObjectIndex_=list_[0];
}
if ( bestObjectIndex_ >=0 ) {
OsiObject * obj = solver->objects()[bestObjectIndex_];
obj->setWhichWay( bestWhichWay_);
}
if (numberFixed==numberUnsatisfied_&&numberFixed)
returnCode=4;
return returnCode;
} else {
return 1;
}
}
// Given a candidate fill in useful information e.g. estimates
void
OsiPseudoCosts::updateInformation(const OsiBranchingInformation *info,
int branch, OsiHotInfo * hotInfo)
{
int index = hotInfo->whichObject();
assert (index<info->solver_->numberObjects());
const OsiObject * object = info->solver_->object(index);
assert (object->upEstimate()>0.0&&object->downEstimate()>0.0);
assert (branch<2);
if (branch) {
if (hotInfo->upStatus()!=1) {
assert (hotInfo->upStatus()>=0);
upTotalChange_[index] += hotInfo->upChange()/object->upEstimate();
upNumber_[index]++;
} else {
#if 0
// infeasible - just say expensive
if (info->cutoff_<1.0e50)
upTotalChange_[index] += 2.0*(info->cutoff_-info->objectiveValue_)/object->upEstimate();
else
upTotalChange_[index] += 2.0*fabs(info->objectiveValue_)/object->upEstimate();
#endif
}
} else {
if (hotInfo->downStatus()!=1) {
assert (hotInfo->downStatus()>=0);
downTotalChange_[index] += hotInfo->downChange()/object->downEstimate();
downNumber_[index]++;
} else {
#if 0
// infeasible - just say expensive
if (info->cutoff_<1.0e50)
downTotalChange_[index] += 2.0*(info->cutoff_-info->objectiveValue_)/object->downEstimate();
else
downTotalChange_[index] += 2.0*fabs(info->objectiveValue_)/object->downEstimate();
#endif
}
}
}
#if 1
// Given a branch fill in useful information e.g. estimates
void
OsiPseudoCosts::updateInformation(int index, int branch,
double changeInObjective,
double changeInValue,
int status)
{
//assert (index<solver_->numberObjects());
assert (branch<2);
assert (changeInValue>0.0);
assert (branch<2);
if (branch) {
if (status!=1) {
assert (status>=0);
upTotalChange_[index] += changeInObjective/changeInValue;
upNumber_[index]++;
}
} else {
if (status!=1) {
assert (status>=0);
downTotalChange_[index] += changeInObjective/changeInValue;
downNumber_[index]++;
}
}
}
#endif
OsiHotInfo::OsiHotInfo() :
originalObjectiveValue_(COIN_DBL_MAX),
changes_(NULL),
iterationCounts_(NULL),
statuses_(NULL),
branchingObject_(NULL),
whichObject_(-1)
{
}
OsiHotInfo::OsiHotInfo(OsiSolverInterface * solver,
const OsiBranchingInformation * info,
const OsiObject * const * objects,
int whichObject) :
originalObjectiveValue_(COIN_DBL_MAX),
whichObject_(whichObject)
{
originalObjectiveValue_ = info->objectiveValue_;
const OsiObject * object = objects[whichObject_];
// create object - "down" first
branchingObject_ = object->createBranch(solver,info,0);
// create arrays
int numberBranches = branchingObject_->numberBranches();
changes_ = new double [numberBranches];
iterationCounts_ = new int [numberBranches];
statuses_ = new int [numberBranches];
CoinZeroN(changes_,numberBranches);
CoinZeroN(iterationCounts_,numberBranches);
CoinFillN(statuses_,numberBranches,-1);
}
OsiHotInfo::OsiHotInfo(const OsiHotInfo & rhs)
{
originalObjectiveValue_ = rhs.originalObjectiveValue_;
whichObject_ = rhs.whichObject_;
if (rhs.branchingObject_) {
branchingObject_ = rhs.branchingObject_->clone();
int numberBranches = branchingObject_->numberBranches();
changes_ = CoinCopyOfArray(rhs.changes_,numberBranches);
iterationCounts_ = CoinCopyOfArray(rhs.iterationCounts_,numberBranches);
statuses_ = CoinCopyOfArray(rhs.statuses_,numberBranches);
} else {
branchingObject_ = NULL;
changes_ = NULL;
iterationCounts_ = NULL;
statuses_ = NULL;
}
}
OsiHotInfo &
OsiHotInfo::operator=(const OsiHotInfo & rhs)
{
if (this != &rhs) {
delete branchingObject_;
delete [] changes_;
delete [] iterationCounts_;
delete [] statuses_;
originalObjectiveValue_ = rhs.originalObjectiveValue_;
whichObject_ = rhs.whichObject_;
if (rhs.branchingObject_) {
branchingObject_ = rhs.branchingObject_->clone();
int numberBranches = branchingObject_->numberBranches();
changes_ = CoinCopyOfArray(rhs.changes_,numberBranches);
iterationCounts_ = CoinCopyOfArray(rhs.iterationCounts_,numberBranches);
statuses_ = CoinCopyOfArray(rhs.statuses_,numberBranches);
} else {
branchingObject_ = NULL;
changes_ = NULL;
iterationCounts_ = NULL;
statuses_ = NULL;
}
}
return *this;
}
OsiHotInfo::~OsiHotInfo ()
{
delete branchingObject_;
delete [] changes_;
delete [] iterationCounts_;
delete [] statuses_;
}
// Clone
OsiHotInfo *
OsiHotInfo::clone() const
{
return new OsiHotInfo(*this);
}
/* Fill in useful information after strong branch
*/
int OsiHotInfo::updateInformation( const OsiSolverInterface * solver, const OsiBranchingInformation * info,
OsiChooseVariable * choose)
{
int iBranch = branchingObject_->branchIndex()-1;
assert (iBranch>=0&&iBranch<branchingObject_->numberBranches());
iterationCounts_[iBranch] += solver->getIterationCount();
int status;
if (solver->isProvenOptimal())
status=0; // optimal
else if (solver->isIterationLimitReached()
&&!solver->isDualObjectiveLimitReached())
status=2; // unknown
else
status=1; // infeasible
// Could do something different if we can't trust
double newObjectiveValue = solver->getObjSense()*solver->getObjValue();
changes_[iBranch] =CoinMax(0.0,newObjectiveValue-originalObjectiveValue_);
// we might have got here by primal
if (choose->trustStrongForBound()) {
if (!status&&newObjectiveValue>=info->cutoff_) {
status=1; // infeasible
changes_[iBranch] = 1.0e100;
}
}
statuses_[iBranch] = status;
if (!status&&choose->trustStrongForSolution()&&newObjectiveValue<choose->goodObjectiveValue()) {
// check if solution
const OsiSolverInterface * saveSolver = info->solver_;
info->solver_=solver;
const double * saveLower = info->lower_;
info->lower_ = solver->getColLower();
const double * saveUpper = info->upper_;
info->upper_ = solver->getColUpper();
// also need to make sure bounds OK as may not be info solver
#if 0
if (saveSolver->getMatrixByCol()) {
const CoinBigIndex * columnStart = info->columnStart_;
assert (saveSolver->getMatrixByCol()->getVectorStarts()==columnStart);
}
#endif
if (choose->feasibleSolution(info,solver->getColSolution(),solver->numberObjects(),
const_cast<const OsiObject **> (solver->objects()))) {
// put solution somewhere
choose->saveSolution(solver);
status=3;
}
info->solver_=saveSolver;
info->lower_ = saveLower;
info->upper_ = saveUpper;
}
// Now update - possible strong branching info
choose->updateInformation( info,iBranch,this);
return status;
}