[Dip] Dip/py code issues, suggestions for additional functionality and question about branching decisions
Marcus Kaiser
marcus.kaiser at mytum.de
Tue Mar 29 09:18:15 EDT 2016
Hello Ted,
this is a revision of my previous mail with added code snippets. They
are excerpts of a copy of DIP source code version 0.92.2 to which I
added some changes, and they show the fixes I made to address the
mentioned issues. It goes without saying that I do not claim the code to
be flawless. Nevertheless, I hope it is helpful to you.
Best regards
Marcus
------------------------------------------------------------------------
/
[///Implementation/]/ I found Dippy to suffer from *memory leaking*.
Since it interfaces with Python, it is responsible for maintaning the
reference counters of the Python objects it deals with. This is done in
some parts of the code, but by far not all. A crucial situation which
came to my attention is the method DippyDecompApp::solveRelaxed. The
retrieved columns are converted to objects of the class DecompVar, yet
the reference counter on the original Python object is not decreased,
which prevents them from being deleted. Hence, for large problem
instances the memory floods.
_/DippyDecompApp.cpp:317
/__/
/_DecompSolverStatus DippyDecompApp::solveRelaxed(const int whichBlock,
const double* redCostX,
const double convexDual,
DecompVarList& varList)
{
if (!m_pySolveRelaxed) {
return DecompSolStatNoSolution;
}
PyObject* pRelaxKey = PyList_GetItem(m_relaxedKeys, whichBlock);
PyObject* pRedCostList = pyTupleList_FromDoubleArray(redCostX,
m_colList);
PyObject* pConvexDual = PyFloat_FromDouble(convexDual);
// call solveRelaxed on DipProblem
PyObject* pStatandVarList = PyObject_CallMethod(m_pProb,
"solveRelaxed", "OOd",
pRelaxKey,
pRedCostList,
pConvexDual);
*
** Py_DECREF(pRedCostList);**
** Py_DECREF(pConvexDual);*
if ( (pStatandVarList == NULL) || (pStatandVarList == Py_None) ){
throw UtilException("Error calling method prob.solveRelaxed()",
"solveRelaxed", "DippyDecompApp");
}
// [status, varList] = relaxed_solver(...)
PyObject * pStatus = PyTuple_GetItem(pStatandVarList, 0);
int cStatus = PyInt_AsLong(pStatus);
DecompSolverStatus status = (DecompSolverStatus)cStatus;
PyObject * pVarList = PyTuple_GetItem(pStatandVarList, 1);
int nVars = PyObject_Length(pVarList);
// In the new design, we need to allow the possibility that the user
will solve
// the problem exactly, but not find any solutions with reduced
costs zero
// The below is is commented out and left in the source for posterity
// tkr 11/11/15
//if (nVars == 0) {
// throw UtilException("Empty variable list", "solveRelaxed",
"DippyDecompApp");
//}
// solveRelaxed returns 3-tuples (cost, reduced cost, dictionary of
(variable, value) pairs)
// We can use these to construct a C++ DecompVar objects
double cost, rc;
PyObject* pTuple, *pDict, *pKeys, *pCol;
string name;
double value;
for (int j = 0; j < nVars; j++) {
pTuple = PySequence_GetItem(pVarList, j);
cost = PyFloat_AsDouble(PyTuple_GetItem(pTuple, 0));
rc = PyFloat_AsDouble(PyTuple_GetItem(pTuple, 1));
pDict = PyTuple_GetItem(pTuple, 2);
pKeys = PyDict_Keys(pDict);
vector<int> varInds;
vector<double> varVals;
for (int n = 0; n < PyObject_Length(pDict); n++) {
pCol = PyList_GetItem(pKeys, n);
value = PyFloat_AsDouble(PyDict_GetItem(pDict, pCol));
varInds.push_back(m_colIndices[pCol]);
varVals.push_back(value);
}
**
** Py_DECREF(****pKeys);*****
** Py_DECREF(pTuple);*****
DecompVar* var = new DecompVar(varInds, varVals, rc, cost);
var->setBlockId(whichBlock);
varList.push_back(var);
}
* Py_DECREF(pStatandVarList);*
return status;
}_/
/_
------------------------------------------------------------------------
/[Implementation]/ The method AlpsDecompTreeNode::getBranchedVar seems
to be meant to return the *variable which is branched on*. However, only
one of the four vectors potentially holding branching decisions is
checked. I think both, lower bounds and upper bounds, need to be checked
for at least one branch, e.g. downBranchLB_ and downBranchUB_.
/_AlpsDecompTreeNode.h:76_
/ int getBranchedVar() const {
if (!downBranchLB_.empty()) {
return downBranchLB_[0].first;*
** } else if (!downBranchUB_.empty()) {**
** return downBranchUB_[0].first;*
*// Not sure whether the following is necessary. Seems to be hard to
decide as the getBranchedVar is not "natural", i.e. it does only
consi**der one variable when there could be more.**
**//** } else if (!upBranchLB_.empty()) {**
**// return upBranchLB_[0].first;**
**// } else if (!upBranchUB_.empty()) {**
**// return upBranchUB_[0].first;*
} else {
return -1;
}
}
------------------------------------------------------------------------
/[//Implementation/] In DecompAlgo::processNode, line 1855 the
processing of a node is terminated if the *lower bound meets the global
upper bound*. It takes into account numeric inexactness. The calling
method AlpsDecompTreeNode::process repeats this check in line 309 but in
an exact fashion why the node might not be fathomed when it should be.
In my oppinion this repeated condition should match the original one. Or
maybe even better, DecompAlgo::getStopCriteria could be used?
/_AlpsDecompTreeNode.cpp:309_
/
if (quality_ >= currentUB*- **DecompEpsilon*) {
doFathom = true;
UTIL_DEBUG(param.msgLevel, 3,
cout << "Fathom since thisQuality= "
<< setw(10) << UtilDblToStr(thisQuality)
<< " quality_= " << setw(10) << UtilDblToStr(quality_)
<< " currentUB = " << setw(10) <<
UtilDblToStr(currentUB)
<< " gap = " << setw(10) << UtilDblToStr(gap)
<< endl;
);
}
------------------------------------------------------------------------
/[Implementation]/ The member DecompAlgo::m_relGap is set in
DecompAlgo::updateObjBound and used in DecompAlgo::isGapTight.
DecompAlgo::m_relGap, however, is not reset when entering
DecompAlgo::processNode. Therefore, it has an *invalid value* based on a
node processed before - probably representing a tight *gap*. I think,
this might lead to stopping the processing of a node immediately as it
is mistakenly believed to have a tight gap, cf. DecompAlgo::phaseUpdate,
line 4274. I suggest to reset DecompAlgo::m_relGap appropriately or
replace it completely by calls to DecompAlgo::getNodeLPGap.
/_DecompAlgo.h:387_
/
bool isGapTight() {
//TODO: make param
double tightGap = m_param.MasterGapLimit;
//printf("isGapTight m_relGap = %g\n", m_relGap);
if (m_param.LogDebugLevel >= 2) {
(*m_osLog) << "DW GAP = " << UtilDblToStr(m_relGap)
<< " isTight = " << (m_relGap <= tightGap)
<< "\n";
}
if (*getNodeLPGap()* <= tightGap) {
return true;
} else {
return false;
}
}
------------------------------------------------------------------------
/[Functionality/Performance]/ For the problem I consider, the MILP
formulation of the subproblems is quite huge. Therefore, I use a dummy
formulation and prevent DIP from looking at it by always providing
columns with negative reduced cost via my implementation of
DecompApp::solveRelaxed if they exist. I do so in a two-step approach.
In a first step the subproblem is solved heuristically. If this does not
lead to new columns an exact method is applied. Since the exact method
is costly, I would like to stick to the heuristic as long as there are
new columns for /some/ subproblem. This stands in contrast to the
current interface of DIP as it tries to solve subproblems without new
columns using the MILP formulation.
Furthermore, the subproblems resemble each other for my problem. This
would allow to solve them in an accumulated fashion as far as the
reduced cost and branching decisions allow so. The latter is certainly
true for the root node of the branch-and-bound tree. Hence, I would
suggest to redesign the interface of DIP to enable a solution process
for *all the subproblems at once* if the user provides it. Maybe the
current treatment could act as a fallback.
My approach here was to introduce the method solveRelaxedAll for
DecompApp and DecompAlgo. DecompAlgo::solveRelaxedAll calls
DecompApp::solveRelaxedAll and decides on the return value whether it is
actually implemented by the user (true) or not (false). This result is
propagated back to DecompAlgo::generateVars which falls back on the
original method if there is no implementation. This design is not
perfect at all since there are large parts of code in
DecompAlgo::solveRelaxedAll copied from DecompAlgo::solveRelaxed.
_/DecompApp.h/_
virtual bool solveRelaxedAll(const double* redCost, const double*
convexDuals, DecompVarList& varList, std::vector<DecompSolverStatus>&
states) {
return false;
}
_/DecompAlgo.cpp/_
bool DecompAlgo::solveRelaxedAll(const double* redCost,
const double* origCost,
const double* convexDuals,
const int n_origCols,
DecompSolverResult* solveResult,
DecompVarList& vars,
double timeLimit
)
{
UtilPrintFuncBegin(m_osLog, m_classTag, "solveRelaxedAll()",
m_param.LogDebugLevel, 2);
if (m_param.SubProbParallel) {
m_stats.timerOther1.reset();
} else {
m_stats.timerOther2.reset();
}
* DecompVarList userVars;**
** vector<DecompSolverStatus> solverStatus(m_numConvexCon,
DecompSolStatNoSolution);**
**
** if (m_param.SolveRelaxAsIp != 1) {**
**
** if (!m_app->solveRelaxedAll(redCost, convexDuals, userVars,
solverStatus))**
** return false;**
**
** for (auto var : userVars)**
** if (var->getVarType() == DecompVar_Point)**
**var->setReducedCost(var->getReducedCost() -
convexDuals[var->getBlockId()]);**
** }**
**
** m_isColGenExact = true;**
** for (int subprobIndex = 0; subprobIndex < m_numConvexCon;
subprobIndex++)**
** m_isColGenExact &= solverStatus[subprobIndex] ==
DecompSolStatOptimal;**
**
** if ((m_isColGenExact || userVars.size() > 0) &&
(m_param.SolveRelaxAsIp != 2)) {**
** vars.splice(vars.end(), userVars);**
** } else {*
for (int subprobIndex = 0; subprobIndex < m_numConvexCon;
subprobIndex++) {
DecompSubModel& subModel = getModelRelax(subprobIndex);
OsiSolverInterface* subprobSI = subModel.getOsi();
int whichBlock = subModel.getBlockId();
bool isRoot = getNodeIndex() ? false : true;
DecompConstraintSet* model = subModel.getModel();
bool doCutoff = m_param.SubProbUseCutoff;
bool doExact = m_function == DecompFuncGenerateInitVars ?
false : true;
assert(subprobSI);
subModel.setOsiObjCoeff(redCost);
if (m_param.BranchEnforceInSubProb) {
subModel.setActiveColBounds(m_colLBNode, m_colUBNode);
}
if (m_param.LogDumpModel > 1) {
string baseName = "subProb_" + subModel.getModelName();
if (m_isStrongBranch)
baseName += "_SB";
std::cout << "problem name is " << baseName <<
m_nodeStats.nodeIndex << m_nodeStats.cutCallsTotal <<
m_nodeStats.priceCallsTotal << whichBlock << std::endl;
printCurrentProblem(subprobSI, baseName,
m_nodeStats.nodeIndex, m_nodeStats.cutCallsTotal,
m_nodeStats.priceCallsTotal, whichBlock);
}
subModel.solveAsMIP(solveResult, m_param, doExact,
doCutoff, isRoot, convexDuals[subprobIndex] - DecompEpsilon, timeLimit);
m_isColGenExact = solveResult->m_isOptimal;
if (solveResult->m_nSolutions) {
int k;
int nSol = std::min<int>(solveResult->m_nSolutions,
m_param.SubProbNumSolLimit);
for (k = 0; k < nSol; k++) {
const double* milpSolution =
solveResult->getSolution(k);
vector<int> ind;
vector<double> els;
int i, c;
double varRedCost = 0.0;
double varOrigCost = 0.0;
DecompVarType varType = !solveResult->m_isUnbounded
? DecompVar_Point : DecompVar_Ray;
if (model->isSparse()) {
map<int, int>::const_iterator mcit;
const map<int, int>& sparseToOrig =
model->getMapSparseToOrig();
for (mcit = sparseToOrig.begin();
mcit != sparseToOrig.end(); mcit++) {
i = mcit->first; //sparse-index
c = mcit->second; //original-index
if (!UtilIsZero(milpSolution[i],
m_app->m_param.TolZero)) {
ind.push_back(c);
els.push_back(milpSolution[i]);
varRedCost += redCost[c] * milpSolution[i];
varOrigCost += origCost[c] *
milpSolution[i];
}
}
}
else {
for (c = 0; c < n_origCols; c++) {
if (!UtilIsZero(milpSolution[c],
m_app->m_param.TolZero)) {
ind.push_back(c);
els.push_back(milpSolution[c]);
varRedCost += redCost[c] * milpSolution[c];
varOrigCost += origCost[c] *
milpSolution[c];
}
}
}
if (varType == DecompVar_Point) {
varRedCost -= convexDuals[subprobIndex];
}
DecompVar* var = new DecompVar(ind, els,
varRedCost, varOrigCost, varType);
var->setBlockId(whichBlock);
vars.push_back(var);
}
}
}
* }*
if (!m_param.SubProbParallel) {
m_stats.thisSolveRelax.push_back(m_stats.timerOther1.getRealTime());
}
UtilPrintFuncEnd(m_osLog, m_classTag, "solveRelaxedAll()",
m_param.LogDebugLevel, 2);
return true;
}
_/
Decom/__/pAlgo.cpp:4763/_
* timeLimit = max(m_param.SubProbTimeLimitExact -
m_stats.timerOverall.getRealTime(), 0.0);**
** if (!solveRelaxedAll(redCostX, origObjective, convexDuals,
nCoreCols, &solveResult, potentialVars, timeLimit)) {*
#ifdef _OPENMP
UTIL_DEBUG(m_app->m_param.LogDebugLevel, 3,
(*m_osLog) << "===== START Threaded solve of
subproblems. =====\n";
);
if (m_param.SubProbParallel) {
omp_set_num_threads(min(m_param.NumConcurrentThreadsSubProb,
m_numConvexCon));
}
else {
omp_set_num_threads(1);
}
#endif
DecompVarList* potentialVarsT = new
DecompVarList[m_numConvexCon];
CoinAssertHint(potentialVarsT, "Error: Out of Memory");
#pragma omp parallel for schedule(dynamic,
m_param.SubProbParallelChunksize)
for (int subprobIndex = 0; subprobIndex < m_numConvexCon;
subprobIndex++) {
DecompSubModel& subModel = getModelRelax(subprobIndex);
double alpha = u[nBaseCoreRows + subprobIndex];
DecompSolverResult solveResult(m_infinity);
#ifdef _OPENMP
UTIL_DEBUG(m_app->m_param.LogDebugLevel, 4,
(*m_osLog) << "THREAD " << omp_get_thread_num() <<
" solving subproblem " << subprobIndex << "\n";);
#else
UTIL_DEBUG(m_app->m_param.LogDebugLevel, 4,
(*m_osLog) << "solve relaxed model = " <<
subModel.getModelName() << endl;
);
#endif
timeLimit = max(m_param.SubProbTimeLimitExact -
m_stats.timerOverall.getRealTime(), 0.0);
solveRelaxed(redCostX, origObjective, alpha, nCoreCols,
false, subModel, &solveResult, potentialVarsT[subprobIndex], timeLimit);
if (solveResult.m_isCutoff) {
mostNegRCvec[subprobIndex] =
min(mostNegRCvec[subprobIndex], 0.0);
}
}
for (int subprobIndex = 0; subprobIndex < m_numConvexCon;
subprobIndex++) {
/* printf("arg[%d].vars size=%d\n",
t, static_cast<int>(arg[t].vars->size()));
*/
for (it = potentialVarsT[subprobIndex].begin();
it != potentialVarsT[subprobIndex].end(); it++) {
varRedCost = (*it)->getReducedCost();
whichBlock = (*it)->getBlockId();
if ((*it)->getVarType() == DecompVar_Point) {
alpha = u[nBaseCoreRows + whichBlock];
} else if ( (*it)->getVarType() == DecompVar_Ray) {
alpha = 0;
}
UTIL_DEBUG(m_app->m_param.LogDebugLevel, 3,
(*m_osLog)
<< "alpha[block=" << whichBlock << "]:" << alpha
<< " varRedCost: " << varRedCost << "\n";
);
}
}
#ifdef _OPENMP
UTIL_DEBUG(m_app->m_param.LogDebugLevel, 3,
(*m_osLog)
<< "===== END Threaded solve of subproblems. =====\n";);
#endif
for (int subprobIndex = 0; subprobIndex < m_numConvexCon;
subprobIndex++) {
for (it = potentialVarsT[subprobIndex].begin(); it !=
potentialVarsT[subprobIndex].end(); it++) {
potentialVars.push_back(*it);
}
}
UTIL_DELARR(potentialVarsT);
* }**
**
** UTIL_DELARR(convexDuals);*
------------------------------------------------------------------------
[/Functionality/] Based on the fact that my algorithm for the subproblem
is partly heuristic (see previous paragraphs), it is not until the last
iterations for solving a single node that it provides not only feasible
(DecompSolStatFeasible) but optimal solutions (DecompSolStatOptimal).
Thus, a lower bound for the node can only be computed at the end of the
solution process of each node. This prevents me from using the
*tailing-off mechanism* provided by DIP (DecompAlgo::isTailoffLB) as it
is based on the lower bound. For that reason, I think it would be useful
to introduce an alternative tailing-off control based on the progression
of the upper bound for the relaxed problem. Would this be a reasonable
approach?
------------------------------------------------------------------------
/[Implementation/Functionality]/ If a node is solved to optimality in
the pricing phase (PHASE_PRICE2) no more columns are generated and the
algorithm switches to PHASE_CUT (cf. DecompAlgo::phaseUpdate, line
4215). It prevents stopping on a tight gap as checked in
DecompAlgo::phaseUpdate, line 4274. Thus, the switch to the cutting
phase it carried out. However, there is a parameter called PCStrategy.
Setting it tofavorPrice will make DecompAlgo::phaseUpdate to immediately
switch the phase from PHASE_CUT to PHASE_PRICE2 back again in line 4153
as long as the limit on the number of price calls is not reached. This
can result in alternation of the two phases and hence an*infinite loop*.
I found the remedy of setting the RoundCutItersLimit to 0, which
probably suits my intention better. Yet, I wonder what the actual use of
the PCStrategy parameter is. At the moment it seems to be redundant as
the RoundPriceItersLimit and RoundCutItersLimit are the controlling
paramerters.
------------------------------------------------------------------------
/[Performance/Implementation]/ I found that the *checks for duplicate
and parallel columns* in DecompApp::addVarsToPool, lines 5565 sqq. are
quite expensive. First of all, I believe that the check for parallel
columns in line 5609 is redundant if the parameter ParallelColsLimit is
1.0 (as pointed out in the comment preceeding that line of code). Since
this is the default value of the parameter, I recommend checking for
parallel columns only if the parameter is smaller than 1.0. Apart from
that in my understanding of column generation, columns with negative
reduced cost cannot have been included before. This would render it
unnecessary to check for duplicates in the existing columns. Maybe this
is not true for some configurations like DualStab? As mentioned in the
comments in the code, the hashing of the columns is not efficient at
all. A comparison based on this hashes is even more expensive than a
direct comparison of the sparse representation. A better hashing and
conditional exact comparison would result in a noticable speed-up, I
suppose.
_/DecompAlgo.cpp:5609
/_ if (foundGoodCol *&& m_param.ParallelColsLimit < 1.0*&&
m_varpool.isParallel(m_vars, waitingCol,
m_param.ParallelColsLimit)) {
UTIL_DEBUG(m_app->m_param.LogDebugLevel, 3,
(*m_osLog) << "Parallel variable, already in vars.\n";
);
waitingCol.deleteVar();
waitingCol.deleteCol();
if (m_algo != RELAX_AND_CUT) { //??
m_nodeStats.varsThisCall--;
m_nodeStats.varsThisRound--;
}
continue;
}
_/DecompVarPool.cp/__/p:/__/150/_
bool DecompVarPool::isDuplicate(const DecompVarList& vars,
const DecompWaitingCol& wcol)
{
const DecompVar* var = wcol.getVarPtr();
const int block = var->getBlockId();
const int len = var->m_s.getNumElements();
const int* indices = var->m_s.getIndices();
const double* values = var->m_s.getElements();
// const string hash = var->getStrHash();
for (DecompVarList::const_iterator vi = vars.begin(); vi !=
vars.end(); vi++) {
if ((*vi)->getBlockId() != block)
continue;
// Could be reasonalbe if hashing is lightweight and yields an adequate
number of "buckets"
// if ((*vi)->getStrHash() != hash)
// continue;
if ((*vi)->m_s.getNumElements() != len)
continue;
const int* other_indices = (*vi)->m_s.getIndices();
const double* other_values = (*vi)->m_s.getElements();
for (int i = 0; i < len; i++)
if (other_indices[i] != indices[i] || other_values[i] !=
values[i])
goto next;
return true;
next: {}
}
return false;
}
bool DecompVarPool::isDuplicate(const DecompWaitingCol& wcol)
{
// Repetition of code from above. Refactor?
const DecompVar* var = wcol.getVarPtr();
const int block = var->getBlockId();
const int len = var->m_s.getNumElements();
const int* indices = var->m_s.getIndices();
const double* values = var->m_s.getElements();
// const string hash = var->getStrHash();
for (auto vi = begin(); vi != end(); vi++) {
const DecompVar* other = (*vi).getVarPtr();
if (other->getBlockId() != block)
continue;
// Could be reasonalbe if hashing is lightweight and yields an adequate
number of "buckets"
// if (other->getStrHash() != hash)
// continue;
if (other->m_s.getNumElements() != len)
continue;
const int* other_indices = other->m_s.getIndices();
const double* other_values = other->m_s.getElements();
for (int i = 0; i < len; i++)
if (other_indices[i] != indices[i] || other_values[i] !=
values[i])
goto next;
return true;
next: {}
}
return false;
}
------------------------------------------------------------------------
/[Question]/ Finally, I did not understand the usage of the parameters
BranchEnforceInMaster and BranchEnforceInSubProb and would be grateful
if you could explain the basic meaning of them to me. They seem to
control the treatment of the branching decisions. BranchEnforceInMaster
suggests to include the decisions as new rows in the master problem and
enforce them via reduced costs, when BranchEnforceInSubProb suggests to
make it the subproblems' task to deal with them. In the latter case,
what happens to conditions which include master-only variables and
therefore cannot be treated in the subproblems? What is the best way for
getting the branching decisions for the current node when using
BranchEnforceInSubProb?
------------------------------------------------------------------------
Am 27.03.2016 um 00:53 schrieb Ted Ralphs:
> Hi Marcus,
>
> Thanks very much for the detailed feedback! Sorry for the delay in
> responding. I had a quick look over your comments and they all seem
> reasonable. The implementation issues look like they can be fixed
> fairly quickly and easily. For the others, we need to have a more
> detailed look to see what can be done in the short and long run. In
> the next couple of days, I'll try to respond with detailed comments on
> each of your issues. Of course, if you want to provide patches for any
> of these, that would be most welcome :).
>
> Cheers,
>
> Ted
>
>
>
> On Mon, Mar 21, 2016 at 9:59 AM, Marcus Kaiser <marcus.kaiser at mytum.de
> <mailto:marcus.kaiser at mytum.de>> wrote:
>
> Hello Dip/py community,
> in my master's thesis im dealing with a Dantzig-Wolfe
> decomposition and the resulting column generation. For the
> implementation I use the DIP framework, version 0.92.2 on Windows
> 7 with MSVC 2013. I appreciate that you provide such an extensive
> framework. While working with it and debbuging my code the
> following things came to my mind. Some of them are performance and
> implementation related others address the functionality of the
> algorithms in DIP. Additionally, I have a question about the
> treatment of branching decisions.
>
> /[///Implementation/]/ I found Dippy to suffer from *memory
> leaking*. Since it interfaces with Python, it is responsible for
> maintaning the reference counters of the Python objects it deals
> with. This is done in some parts of the code, but by far not all.
> A crucial situation which came to my attention is the method
> DippyDecompApp::solveRelaxed. The retrieved columns are converted
> to objects of the class DecompVar, yet the reference counter on
> the original Python object is not decreased, which prevents them
> from being deleted. Hence, for large problem instances the memory
> floods.
>
> /[Implementation]/ The method AlpsDecompTreeNode::getBranchedVar
> seems to be meant to return the *variable which is branched on*.
> However, only one of the four vectors potentially holding
> branching decisions is checked. I think both, lower bounds and
> upper bounds, need to be checked for at least one branch, e.g.
> downBranchLB_ and downBranchUB_.
>
> /[//Implementation/] In DecompAlgo::processNode, line 1855 the
> processing of a node is terminated if the *lower bound meets the
> global upper bound*. It takes into account numeric inexactness.
> The calling method AlpsDecompTreeNode::process repeats this check
> in line 309 but in an exact fashion why the node might not be
> fathomed when it should be. In my oppinion this repeated condition
> should match the original one. Or maybe even better,
> DecompAlgo::getStopCriteria could be used?
>
> /[Implementation]/ The member DecompAlgo::m_relGap is set in
> DecompAlgo::updateObjBound and used in DecompAlgo::isGapTight.
> DecompAlgo::m_relGap, however, is not reset when entering
> DecompAlgo::processNode. Therefore, it has an *invalid value*
> based on a node processed before - probably representing a tight
> *gap*. I think, this might lead to stopping the processing of a
> node immediately as it is mistakenly believed to have a tight gap,
> cf. DecompAlgo::phaseUpdate, line 4274. I suggest to reset
> DecompAlgo::m_relGap appropriately or replace it completely by
> calls to DecompAlgo::getNodeLPGap.
>
> /[Functionality/Performance]/ For the problem I consider, the MILP
> formulation of the subproblems is quite huge. Therefore, I use a
> dummy formulation and prevent DIP from looking at it by always
> providing columns with negative reduced cost via my implementation
> of DecompApp::solveRelaxed if they exist. I do so in a two-step
> approach. In a first step the subproblem is solved heuristically.
> If this does not lead to new columns an exact method is applied.
> Since the exact method is costly, I would like to stick to the
> heuristic as long as there are new columns for /some/ subproblem.
> This stands in contrast to the current interface of DIP as it
> tries to solve subproblems without new columns using the MILP
> formulation.
> Furthermore, the subproblems resemble each other for my problem.
> This would allow to solve them in an accumulated fashion as far as
> the reduced cost and branching decisions allow so. The latter is
> certainly true for the root node of the branch-and-bound tree.
> Hence, I would suggest to redesign the interface of DIP to enable
> a solution process for *all the subproblems at once* if the user
> provides it. Maybe the current treatment could act as a fallback.
>
> [/Functionality/] Based on the fact that my algorithm for the
> subproblem is partly heuristic (see previous paragraphs), it is
> not until the last iterations for solving a single node that it
> provides not only feasible (DecompSolStatFeasible) but optimal
> solutions (DecompSolStatOptimal). Thus, a lower bound for the node
> can only be computed at the end of the solution process of each
> node. This prevents me from using the *tailing-off mechanism*
> provided by DIP (DecompAlgo::isTailoffLB) as it is based on the
> lower bound. For that reason, I think it would be useful to
> introduce an alternative tailing-off control based on the
> progression of the upper bound for the relaxed problem. Would this
> be a reasonable approach?
>
> /[Implementation/Functionality]/ If a node is solved to optimality
> in the pricing phase (PHASE_PRICE2) no more columns are generated
> and the algorithm switches to PHASE_CUT (cf.
> DecompAlgo::phaseUpdate, line 4215). It prevents stopping on a
> tight gap as checked in DecompAlgo::phaseUpdate, line 4274. Thus,
> the switch to the cutting phase it carried out. However, there is
> a parameter called PCStrategy. Setting it tofavorPrice will make
> DecompAlgo::phaseUpdate to immediately switch the phase from
> PHASE_CUT to PHASE_PRICE2 back again in line 4153 as long as the
> limit on the number of price calls is not reached. This can result
> in alternation of the two phases and hence an*infinite loop*. I
> found the remedy of setting the RoundCutItersLimit to 0, which
> probably suits my intention better. Yet, I wonder what the actual
> use of the PCStrategy parameter is. At the moment it seems to be
> redundant as the RoundPriceItersLimit and RoundCutItersLimit are
> the controlling paramerters.
>
> /[Performance/Implementation]/ I found that the *checks for
> duplicate and parallel columns* in DecompApp::addVarsToPool, lines
> 5565 sqq. are quite expensive. First of all, I believe that the
> check for parallel columns in line 5609 is redundant if the
> parameter ParallelColsLimit is 1.0 (as pointed out in the comment
> preceeding that line of code). Since this is the default value of
> the parameter, I recommend checking for parallel columns only if
> the parameter is smaller than 1.0. Apart from that in my
> understanding of column generation, columns with negative reduced
> cost cannot have been included before. This would render it
> unnecessary to check for duplicates in the existing columns. Maybe
> this is not true for some configurations like DualStab? As
> mentioned in the comments in the code, the hashing of the columns
> is not efficient at all. A comparison based on this hashes is even
> more expensive than a direct comparison of the sparse
> representation. A better hashing and conditional exact comparison
> would result in a noticable speed-up, I suppose.
>
> /[Question]/ Finally, I did not understand the usage of the
> parameters BranchEnforceInMaster and BranchEnforceInSubProb and
> would be grateful if you could explain the basic meaning of them
> to me. They seem to control the treatment of the branching
> decisions. BranchEnforceInMaster suggests to include the decisions
> as new rows in the master problem and enforce them via reduced
> costs, when BranchEnforceInSubProb suggests to make it the
> subproblems' task to deal with them. In the latter case, what
> happens to conditions which include master-only variables and
> therefore cannot be treated in the subproblems? What is the best
> way for getting the branching decisions for the current node when
> using BranchEnforceInSubProb?
>
> Thank you in advance,
> Marcus
>
> _______________________________________________
> Dip mailing list
> Dip at list.coin-or.org <mailto:Dip at list.coin-or.org>
> http://list.coin-or.org/mailman/listinfo/dip
>
>
>
>
> --
> Dr. Ted Ralphs
> Professor, Lehigh University
> (610) 628-1280
> ted 'at' lehigh 'dot' edu
> coral.ie.lehigh.edu/~ted <http://coral.ie.lehigh.edu/%7Eted>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://list.coin-or.org/pipermail/dip/attachments/20160329/31e7f2a6/attachment.html>
More information about the Dip
mailing list