Re: [Fwd: Bug in glpk integer programmer]

[Heinrich Schuchardt]

On 15.10.20 16:03, Andrew Makhorin wrote:
> -------- Forwarded Message --------
> From: s meagher <sjmeagher@gmail.com>
> To: bug-glpk@gnu.org
> Subject: Bug in glpk integer programmer
> Date: Fri, 16 Oct 2020 00:19:14 +1100
>
>> Hi
>>
>> The following code demonstrates a possible bug in glpk. I fully
>> appreciate that in fact the bug may be that I do understand glpk
>> properly and have not configured it correctly. If this is the case,
>> please accept my apologies for taking up your time.
>>
>> Kind regards
>>
>> Steve Meagher
>>
>> /*   Possible integer programming bug in glpk v4.65
>>
>>      (c) Stephen Meagher 2020
>>
>>      The following code was written to demonstrate a bug I believed I
>> found for integer programming in glpk. Instead, I seem to have

You did not define any variable as integer.

>> discovered a different issue - which may be just because I'm using
>> glpk incorrectly.
>>
>>      The following integer programme
>>
>>              objective function: f(x_1,x_2) = 0
>>              row 1                           x_1 + x_2 = 3
>>              row 2                   0<=     x_1 - x_2 <= 2
>>
>> clearly has feasible solution x1 = 2, x2 = 1, but if I set it up as
>>
>>              objective function: f(x_1,x_2) = 0
>>              row 1                           x_1 + x_2 = 3
>>              row 2                   0<=     x_1 - x_2 <=
>> 2*(1.05)
>>
>> then glpk tells me the solution is x1 = 2.55 and x2 = 0.45, which is
>> wrong for two reasons. (Note, even if I set it up with 2 instead of
>> 2*(1.05), glpk gives the wrong answer when configured as follows).
>>
>>      I have compiled this programme on a MacBook Pro running macOS
>> Catalina 10.15.6 using the command: 
>>
>>      gcc -std=c18 glpkipbug.c -lglpk
>>
>>      I have the gnu big num library installed and use it to configure
>> glpk when I installed it.
>>
>>      Here is the fulloutput:
>>
>>      GLPK Simplex Optimizer, v4.65
>> 2 rows, 2 columns, 4 non-zeros
>>       0: obj =   0.000000000e+00 inf =   3.000e+00 (1)
>>       2: obj =   0.000000000e+00 inf =   0.000e+00 (0)
>> OPTIMAL LP SOLUTION FOUND
>> GLPK Integer Optimizer, v4.65
>> 2 rows, 2 columns, 4 non-zeros
>> 2 integer variables, none of which are binary
>> Preprocessing...
>> 2 rows, 2 columns, 4 non-zeros
>> 2 integer variables, none of which are binary
>> Scaling...
>>  A: min|aij| =  1.000e+00  max|aij| =  1.000e+00  ratio =  1.000e+00
>> Problem data seem to be well scaled
>> Constructing initial basis...
>> Size of triangular part is 2
>> Solving LP relaxation...
>> GLPK Simplex Optimizer, v4.65
>> 2 rows, 2 columns, 4 non-zeros
>>       2: obj =   0.000000000e+00 inf =   3.000e+00 (1)
>>       3: obj =   0.000000000e+00 inf =   0.000e+00 (0)
>> OPTIMAL LP SOLUTION FOUND
>> Integer optimization begins...
>> Long-step dual simplex will be used
>> +     3: mip =     not found yet >=              -inf        (1; 0)
>> +     4: >>>>>   0.000000000e+00 >=   0.000000000e+00   0.0% (1; 0)
>> +     4: mip =   0.000000000e+00 >=     tree is empty   0.0% (0; 1)
>> INTEGER OPTIMAL SOLUTION FOUND
>> Intopt success
>> Solution as x1 = 2.550000, x2 = 0.450000
>>
>>
>> (If you're still reading: The original bug occurred when I set certain
>> row bounds to be inequalities. Then the rows which were equalities
>> were no longer satisfied. However, this bug disappeared if I set the
>> bounds to be integers instead of real numbers. Furthermore, I wrote
>> that programme in C++ and not C.)
>>
>> */
>>
>>
>> #include <glpk.h>
>> #include <stdio.h>
>>
>> /* Defines and demontrates the bug */
>> void bug() {
>>      glp_prob *lp;
>>      glp_iocp parm;
>>      int ia[5];
>>      int ja[5];
>>      double arr[5];
>>      lp = glp_create_prob();
>>      glp_set_obj_dir(lp,GLP_MIN);
>>      glp_init_iocp(&parm);
>>      parm.presolve = GLP_ON;
>>      glp_add_rows(lp,2);
>>      glp_set_row_name(lp,1,"r1");
>>      glp_set_row_bnds(lp,1,GLP_FX,3.0,3.0);
>>      glp_set_row_name(lp,2,"r2");
>>
>>      glp_set_row_bnds(lp,2,GLP_DB, 0,2*(1.05));      // non-integer
>> bounds
>>
>>      glp_add_cols(lp,2);
>>      glp_set_col_name(lp,1,"x1");
>>      glp_set_obj_coef(lp,1,0.0);
>>      glp_set_col_bnds(lp,1,GLP_LO,0.0,0.0);
>>      glp_set_col_kind(lp,1,GLP_IV);
>>      glp_set_col_name(lp,2,"x2");
>>      glp_set_obj_coef(lp,2,0.0);
>>      glp_set_col_bnds(lp,2,GLP_LO,0.0,0.0);
>>      glp_set_col_kind(lp,2,GLP_IV);
>>      ia[1] = 1; ja[1] = 1; arr[1] = 1; ia[2] = 1; ja[2] = 2; arr[2] =
>> 1;
>>      ia[3] = 2; ja[3] = 1; arr[3] = 1; ia[4] = 2; ja[4] = 2; arr[4] =
>> -1;
>>      glp_load_matrix(lp,4,ia,ja,arr);
>>      glp_simplex(lp,NULL);
>>      if (glp_intopt(lp, &parm) != 0) {
>>              printf("Intopt error\n");
>>      } else {
>>              printf("Intopt success\n");
>>      }
>>      double x1 = glp_get_col_prim(lp,1);
>>      double x2 = glp_get_col_prim(lp,2);

You should use glp_mip_col_val() here to get

x1 = 2.000000, x2 = 1.000000

Best regards

Heinrich

>>      printf("Solution as x1 = %lf, x2 = %lf\n",x1,x2);
>>      glp_delete_prob(lp);
>> }
>>
>> int main() {
>>      bug();
>>      return 0;
>> }
>>
>>
>>
>