bug-glpk
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[Bug-glpk] bug detected in the cut generator module (glpios03)


From: Andrew Makhorin
Subject: [Bug-glpk] bug detected in the cut generator module (glpios03)
Date: Tue, 5 Feb 2008 12:30:53 +0300

A bug has been detected in the cut generator module (glpios03).
The bug appears if the current subproblem to be solved has no integer
feasible solution, and the cut generator generates a cut, which has
zero coefficients and non-zero rhs. This leads to zero divide error
on normalizing the coefficient vector, since its euclidean norm is zero
in such cases.

The bug will be fixed in a next release of the package.

To fix the bug right now please replace internal routines efficacy and
parallel (file glpios03.c, lines 1558-1594) by the following ones:

static double efficacy(glp_tree *tree, IOSCUT *cut)
{     glp_prob *mip = tree->mip;
      IOSAIJ *aij;
      double s = 0.0, t = 0.0, temp;
      for (aij = cut->ptr; aij != NULL; aij = aij->next)
      {  xassert(1 <= aij->j && aij->j <= mip->n);
         s += aij->val * mip->col[aij->j]->prim;
         t += aij->val * aij->val;
      }
      temp = sqrt(t);
      if (temp < DBL_EPSILON) temp = DBL_EPSILON;
      switch (cut->type)
      {  case GLP_LO:
            temp = (s >= cut->rhs ? 0.0 : (cut->rhs - s) / temp);
            break;
         case GLP_UP:
            temp = (s <= cut->rhs ? 0.0 : (s - cut->rhs) / temp);
            break;
         default:
            xassert(cut != cut);
      }
      return temp;
}

static double parallel(IOSCUT *a, IOSCUT *b, double work[])
{     IOSAIJ *aij;
      double s = 0.0, sa = 0.0, sb = 0.0, temp;
      for (aij = a->ptr; aij != NULL; aij = aij->next)
      {  work[aij->j] = aij->val;
         sa += aij->val * aij->val;
      }
      for (aij = b->ptr; aij != NULL; aij = aij->next)
      {  s += work[aij->j] * aij->val;
         sb += aij->val * aij->val;
      }
      for (aij = a->ptr; aij != NULL; aij = aij->next)
         work[aij->j] = 0.0;
      temp = sqrt(sa * sb);
      if (temp < DBL_EPSILON) temp = DBL_EPSILON;
      return s / temp;
}






reply via email to

[Prev in Thread] Current Thread [Next in Thread]