Re: [Help-glpk] non-official updated version of glpk (4.63 pre-release)

[Andrew Makhorin]

Hi Chris,

> 
>         A "smart" LP perturbation was also implemented in the dual
>         simplex
>         solver. This feature is similar to the one implemented in the
>         primal
>         simplex solver (see below).
> 
> 
> I did some testing of this and it seems like the dual simplex is less
> stable than the previous version.
> 
> An example to see this is problem sp98ir from miplib 2010. Trying to
> solve this with --pcost I see lots of warnings like:
> Warning: numerical instability (dual simplex, phase II)
> or:
> Warning: numerical instability (dual simplex, phase I)
> and finally:
> Warning: dual simplex failed due to excessive numerical instability
> 
> 
> and later on it seems to get stuck in a sequence like:
> *279500: obj =   2.196544071e+08 inf =   1.239e-13 (2)
> *280000: obj =   2.196544071e+08 inf =   1.239e-13 (2)
> *280500: obj =   2.196544071e+08 inf =   1.239e-13 (2)
> ...
> *9366000: obj =   2.196544071e+08 inf =   1.239e-13 (2)
> *9366500: obj =   2.196544071e+08 inf =   1.239e-13 (2)
> *9367000: obj =   2.196544071e+08 inf =   1.239e-13 (2)
> and so on.
> 
> 
> 
> With version 4.62 this problem is solved rather quickly with only a
> few occasional:
> Warning: numerical instability (dual simplex, phase II)
> 
> messages.
> 

Thank you for evaluation.

I think you could solve sp98ir with 4.62 due to a happy chance. The only
difference between 4.62 and the most recent 4.63 is that in the latter
version the LP perturbation is activated later, i.e. not at the very
beginning.

In case of sp98ir (as well as in cases of similar instances) numerical
instability in the dual simplex solver happens because that instance has
quite large objective coefficients (about 1e8). For the same reason the
primal simplex solver (called if the dual simplex solver fails) gets
into an infinite loop; namely, at the optimum the reduced costs remain
quite large, so the termination test doesn't pass. Obviously,
multiplying the objective by a constant leads to proportional changing
of the reduced costs, so when I reduced the original objective
coefficients by multiplying them by 1e-6, the numerical instability
disappeared.

This kind of numerical difficulties can be avoid by scaling the
objective (internally, within the solver). However, I didn't implement
this feature yet, because simple scaling may not help. Ideally, I'd like
to replace the original objective by an equivalent one, which is
something like a normalized projection of the original objective onto
the affine subspace generated by the equality constraints, but I still
have no idea how to efficiently compute it.

Best regards,

Andrew Makhorin