Re: Solver delivers wrong answer when 2 constraints are close

Looking through the code I can see that it seems to be a mismatch in 
checking for lower/upper bound in other places too (see bellow).

Also testing in AMPL it only returns 'x = 1' for 'min_bound = 1e-06', 
and it works for 'min_bound = 1e-05' and return 'x = 1.00001'.

===

int npp_implied_value(NPP *npp, NPPCOL *q, glp_double s)
{

...

       /* check current column lower bound */
       if (q->lb != -GLP_DBL_MAX)
       {  eps = (q->is_int ? GLP_NPP_EPS5 : GLP_NPP_EPS_5_8_fabs(q->lb));
...

       /* check current column upper bound */
       if (q->ub != +GLP_DBL_MAX)
       {  eps = (q->is_int ? GLP_NPP_EPS5 : GLP_NPP_EPS_5_8_fabs(q->ub));
...

}

...

int npp_implied_lower(NPP *npp, NPPCOL *q, glp_double l)
{     /* process implied column lower bound */
...

       /* check current column lower bound */
       if (q->lb != -GLP_DBL_MAX)
       {  eps = (q->is_int ? GLP_NPP_EPS : GLP_NPP_EPS_3_6_fabs(q->lb));
...

       /* check current column upper bound */
       if (q->ub != +GLP_DBL_MAX)
       {  eps = (q->is_int ? GLP_NPP_EPS5 : GLP_NPP_EPS_5_8_fabs(q->ub));
...

}

...

int npp_implied_upper(NPP *npp, NPPCOL *q, glp_double u)
{     int ret;
...

       /* check current column upper bound */
       if (q->ub != +GLP_DBL_MAX)
       {  eps = (q->is_int ? GLP_NPP_EPS : GLP_NPP_EPS_3_6_fabs(q->ub));
...

       /* check current column lower bound */
       if (q->lb != -GLP_DBL_MAX)
       {  eps = (q->is_int ? GLP_NPP_EPS5 : GLP_NPP_EPS_5_8_fabs(q->lb));
...

}

===

Cheers !

On 4/3/21 18:19, Thiago Neves wrote:
> Thank you so much for the explanation Antti.
> You are right, this seems more like a "design flaw" than a "bug". I'm 
> sorry for the misconception.
>
> Why does "eps" need to be so big (1e-3)? Is there a way to set a fixed 
> threshold in the command line (glpsol) or at least the 1e-3 part?
>
> In the case (1), wouldn't it be better to choose max(l[q], l'[q])?
>
> Att,
>
> *Thiago H. Neves*
> (31) 98608-0666
>
>
>
> ---------- Forwarded message ---------
> De: *Antti Lehtila* <[email protected] <mailto:[email protected]>>
> Date: qui., 4 de mar. de 2021 às 12:44
> Subject: Re: Solver delivers wrong answer when 2 constraints are close
> To: <[email protected] <mailto:[email protected]>>
>
>
> Hi,
>
> I think it works fully as documented, and so *per design*. Singleton
> rows are treated as column bounds by the preprocessor. See documentation
> for *npp_implied_lower*:
> *---------------
> *  Processing implied column lower bound l'[q] includes the following
> *  cases:
> *
> *  1) if l'[q] < l[q] + eps, implied lower bound is redundant;
> *
> *  2) if l[q] + eps <= l[q] <= u[q] + eps, current column lower bound
> *     l[q] can be strengthened by replacing it with l'[q]. If in this
> *     case new column lower bound becomes close to current column upper
> *     bound u[q], the column can be fixed on its upper bound;
> *
> *  3) if l'[q] > u[q] + eps, implied lower bound violates current
> *     column upper bound u[q], in which case the problem has no primal
> *     feasible solution. */
> *---------------
> The lower bound can have only a single value, but if you define multiple
> values for a column lower bound, they must of course be processed in
> some order.  In this case, the lower bound is first defined l(q)=0, then
> l'(q)=1, and finally l'(q)=1.0001. The third value for the bound is thus
> considered *redundant*, as per design, and so the second value l(q)=1
> remains in effect. This is because *eps* is defined as 1e-3 + 1e-6 *
> fabs(q->lb)).
>
> I guess it may be a considered a design flaw, but I think it should not
> be called a bug, as it is working as designed and documented. Besides,
> I think one should use the Bounds section for bounds, instead of using
> multiple constraints for defining a single lower bound.
>
> Best,
> Antti
> _______________________
> On 04.03.2021 11:38, Domingo Alvarez Duarte wrote:
> >
> > Testing this problem I discover that if we change the order of
> > constraint declarations it seems to give the expected answer as stated
> > by Thiago (what I think could be another bug).
> >
> > ====
> >
> > /param min_bound default 0;/
> > /var x >= 0;
> >
> > minimize y: x;/
> > *
> > */*s.t. PART_MIN_X: x >= 1 + min_bound;*
> > /
> > /*s.t. LIM_INF_X: x >= 1;
> > *
> > /
> > /solve;
> > display min_bound;
> > display x; # EXPECTED RESULT: X ==  1.0001
> >
> > data;
> > param min_bound := 1e-4;
> > end;/
> >
> > ====
> >
> > Output:
> >
> > ====
> >
> > x.val = 1.0001
> >
> > ====
> >
> > On 3/3/21 19:19, Thiago Neves wrote:
> >> Hi.
> >> I've found a strange behaviour in glpk which I don't know how to fix
> >> nor how to contour it. It seems like GLPK can't distinguish
> >> constraints that differs from about 1e-4.
> >>
> >> Follows simple examples that explain and reproduce the problem.**
> >> *
> >> *
> >> *The first model gives the desired answer (x = 1.0001):*
> >> /
> >> param min_bound default 0;/
> >> /var x >= 0;
> >>
> >> minimize y: x;/
> >> /*
> >> s.t. PART_MIN_X: x >= 1 + min_bound;*
> >>
> >> solve;
> >> display min_bound;
> >> display x; # EXPECTED RESULT: X ==  1.0001
> >>
> >> data;
> >> param min_bound := 1e-4;
> >> end;
> >> /
> >> /_____________________________________/
> >> /OUTPUT:/
> >> /x.val = 1.0001/
> >> /_____________________________________ /
> >>
> >> *Now, if I add a second constraint "close" to the first one, the
> >> solver will deliver an answer that is actually infeasible:*
> >>
> >> /param min_bound default 0;/
> >> /var x >= 0;
> >>
> >> minimize y: x;/
> >>
> >> *s.t. LIM_INF_X: x >= 1;
> >>
> >> */*s.t. PART_MIN_X: x >= 1 + min_bound;*
> >>
> >> solve;
> >> display min_bound;
> >> display x; # EXPECTED RESULT: X ==  1.0001
> >>
> >> data;
> >> param min_bound := 1e-4;
> >> end;/
> >> /_____________________________________/
> >> /OUTPUT:/
> >> x.val = 1
> >> /_____________________________________ /
> >>
> >> *If I change the "min_bound" parameter to 1e-2, the second model
> >> works as expected (x = 1.01):*
> >>
> >> /param min_bound default 0;/
> >> /
> >> /
> >> /var x >= 0;
> >>
> >> minimize y: x;/
> >> *
> >> *
> >> *s.t. LIM_INF_X: x >= 1;
> >>
> >> */*s.t. PART_MIN_X: x >= 1 + min_bound;*
> >>
> >> solve;/
> >> /
> >> display x; # EXPECTED RESULT: X ==  1.01
> >>
> >> data;
> >> param min_bound := 1e-2;
> >> end;/
> >> /_____________________________________/
> >> /OUTPUT:/
> >> x.val = 1.01
> >> /_____________________________________ /
> >> /
> >> /
> >>
> >> Att,
> >>
> >> *Thiago H. Neves*
> >> (31) 98608-0666
> >>
>
>