x-gurobi: count constraint violated in optimal solution

[4er4er4er]

I modified the book example multmip3.mod (see attached files) to a version called x-multmip3.mod that uses various logical operators and functions instead of binary variables (see Total_Cost, Ship_Range, and Max_Serve):

set ORIG;   # origins
set DEST;   # destinations
set PROD;   # products

param supply {ORIG,PROD} >= 0;  # amounts available at origins
param demand {DEST,PROD} >= 0;  # amounts required at destinations

   check {p in PROD}:
      sum {i in ORIG} supply[i,p] = sum {j in DEST} demand[j,p];

param limit {ORIG,DEST} >= 0;   # maximum shipments on routes
param minload > 0;             # minimum nonzero shipment
param maxserve integer > 0;     # maximum destinations served

param vcost {ORIG,DEST,PROD} >= 0; # variable shipment cost on routes
var Trans {ORIG,DEST,PROD} >= 0;   # units to be shipped

param fcost {ORIG,DEST} >= 0;      # fixed cost for using a route
var Ship {i in ORIG, j in DEST} = sum {p in PROD} Trans[i,j,p];

minimize Total_Cost:
   sum {i in ORIG, j in DEST, p in PROD} vcost[i,j,p] * Trans[i,j,p]
 + sum {i in ORIG, j in DEST} if Ship[i,j] >= minload then fcost[i,j];

subject to Supply {i in ORIG, p in PROD}:
   sum {j in DEST} Trans[i,j,p] = supply[i,p];

subject to Demand {j in DEST, p in PROD}:
   sum {i in ORIG} Trans[i,j,p] = demand[j,p];

subject to Ship_Range {i in ORIG, j in DEST}:
   Ship[i,j] = 0  or  minload <= Ship[i,j] <= limit[i,j];

subject to Max_Serve {i in ORIG}:
   count {j in DEST} (Ship[i,j] >= minload) <= maxserve;

When solved using multmip3.mod + multmip3.dat, the optimal value is 235625. When solved using x-multmip3.mod + multmip3.dat, the optimal value is reported in the solver message as 231650, but the objective Total_Cost displays (correctly) as 233150. Also one of the Max_Serve constraints is violated (which accounts for the lower objective value):

ampl: model x-multmip3.mod;
ampl: data multmip3.dat;
ampl: option solver x-gurobi;

ampl: solve;
x-Gurobi 9.5.0: Set parx-Gurobi 9.5.0: optimal solution; objective 231650
177 simplex iterations
1 branching nodes

ampl: display Total_Cost;
Total_Cost = 233150

ampl: display minload, maxserve;
minload = 375
maxserve = 5

ampl: display {i in ORIG} count {j in DEST} (Ship[i,j] >= minload);
count({j in DEST} (Ship[i,j] >= minload)) [*] :=
CLEV  5
GARY  3
PITT  6
;

ampl: display Max_Serve.slack;
Max_Serve.slack [*] :=
CLEV   0
GARY   2
PITT  -1
;

(Note also the glitch in the solver message: x-Gurobi 9.5.0: Set parx-Gurobi 9.5.0: optimal solution; objective 231650.)

[glebbelov]

The reason was a too small eps=1e-6 for strict comparison of continuous variables in MIPConverter. For Ship[i,j] >= minload to count as false it was enough to be 374.999999 which was reported as 375 by Gurobi.

Changed the default comparison eps to 1e-3 and added option cvt:mip:eps for it.

The message glitch: see #158

[4er4er4er]

As I understand the situation, minload is 375, and for some i and j the constraint

Ship[i,j] = 0  or  minload <= Ship[i,j] <= limit[i,j]

is considered by Gurobi to be satisfied by Ship[i,j] = 374.999999. But in the constraint

count {j in DEST} (Ship[i,j] >= minload) <= maxserve

the condition Ship[i,j] >= minload evaluates to False when Ship[i,j] = 374.999999, due to the use of a different tolerance for this comparison. Is that right? If so, I have two questions:

First, is it enough to increase the second tolerance slightly (even to 9e-5) rather than to a much larger value like 1e-3?

Second, why does the optimal value appear as 375 in AMPL rather than as 374.999999? Even when I set option display_precision 0; to request "full" precision is displays, it appears as 375:

ampl: option display_precision 0;
ampl: display Ship;
Ship [*,*] (tr)
:    CLEV  GARY  PITT    :=
DET   625     0   525
FRA   425     0   475
FRE   425   375   375
LAF     0   400   600
LAN   500     0     0
STL   625   625   550
WIN     0     0   375
;

I guess it is being rounded somewhere, but I am curious as to where it is being rounded.

[glebbelov]

Checking Gurobi output (obtained by solving its own LP instance exported), the values are already 375. Probably, when Gurobi checks indicator constraints, it uses the same feasibility tolerance.

9e-5 is almost 1e-4. Are there any solvers with default FeasTol 1e-4 or more? If no, we could set anything as default which seems safe.

Otherwise, we could make cvt:mip:eps a Backend parameter, i.e., solver-dependent, and the Backend would be responsible for setting a custom default value.

[4er4er4er]

I wasn't thinking straight when I said 9e-5. What I intended to ask was: To avoid the problem we're seeing, is it enough to slightly increase the tolerance applied to Ship[i,j] >= minload from the previous value of 1e-6? That could be for example an increase to 2e-6.

It appears that 1e-6 is a very common feasibility tolerance; it's the default in Gurobi, CPLEX, Xpress, Knitro, and MINOS, for example. The user can change this tolerance, however, by setting a "feastol" solver option. Thus our solver interface should (ideally) adapt to the feasibility tolerance setting that the solver is actually going to use.

[4er4er4er]

I am trying to understand what solution you are describing when you refer to "Checking Gurobi output (obtained by solving its own LP instance exported)." I have in mind that there might be two solutions involved:

1. First, at some point, the x-gurobi interface calls Gurobi API functions to solve and to retrieve the optimal values that Gurobi found.
2. Second, after possibly some modifications, the x-gurobi interface returns optimal values to AMPL.

So the first thing I am wondering is, what is the solution (1) returned by the Gurobi API call in this case -- 374.999999 or 375? Then depending on that, I might have some other questions.

[glebbelov]

Gurobi API returns 375.0.

We could have the following mechanism. When cvt:mip:eps<0 (default -1e-6 for example), the tolerance is set to feastol-(cvt:mip:eps). Otherwise it's the given value. Or, just describe the matter in option text.

[glebbelov]

Continued with escrow/#102.