[Help-glpk] Re: Numerical Instability (Dual and Primal ..)

[Andrew Makhorin]

> I have a "numerical instability" problem in my formulation. I do
> not know
> what is the reason.

A similar problem was already reported earlier; please see:

http://mail.gnu.org/archive/html/help-glpk/2003-10/msg00004.html

However, in your case numerical difficulties are caused due to "big M"
formulation.

I believe that the geometric mean scaling (used by default) is not
appropriate for badly scaled mip models, and the equilibration scaling
would be much more appropriate. Unfortunately, you cannot change the
scaling method from the command line, so please replace line 46 in the
file source/glplpx3.c, routine lpx_reset_parms (glpk 4.1):

   lp->scale = 3;

by the following line:

   lp->scale = 1;

to use the equilibration scaling by default. (I will made this change
in the next version of the package.)

I tried to solve your example with lp->scale = 1:

lpt_read_prob: reading LP data from `CM162.0.cplex'...
lpt_read_prob: 372 variables, 676 constraints
lpt_read_prob: 955 lines were read
lpx_simplex: original LP has 676 rows, 372 columns, 1856 non-zeros
lpx_simplex: presolved LP has 647 rows, 254 columns, 1678 non-zeros
lpx_adv_basis: size of triangular part = 647
      0:   objval =   0.000000000e+00   infeas =   1.000000000e+00 (0)
     63:   objval =   0.000000000e+00   infeas =   0.000000000e+00 (0)
*    63:   objval =   0.000000000e+00   infeas =   0.000000000e+00 (0)
*    90:   objval =   5.000000000e+00   infeas =   0.000000000e+00 (0)
OPTIMAL SOLUTION FOUND
Integer optimization begins...
Objective function is integral
+    90: mip =     not found yet; lp =   5.000000000e+00 (1, 0)
lpx_simplex: numerical instability (dual simplex)
lpx_simplex: numerical instability (dual simplex)
lpx_simplex: numerical instability (dual simplex)
lpx_simplex: numerical instability (dual simplex)
+   428: mip =   4.999999000e+00; lp =   5.000000000e+00 (24, 35)
+   451: mip =   5.000000000e+00; lp =   5.000000000e+00 (30, 43)
+   451: mip =   5.000000000e+00; lp =     tree is empty (0, 73)
INTEGER OPTIMAL SOLUTION FOUND
Time used:   0.4 secs
Memory used: 0.9M (902432 bytes)

Note that the message about numerical instability is a warning. If the
solver is able to reach the end, it reports correct results.

Another, more radical way to avoid numerical difficulties might be
reformulating your problem.

Andrew Makhorin