False negatives reported by scipy.optimize.linprog #11326
Description
Sometimes linprog method='interior-point' reports failure despite successful convergence to a solution. That is, in scipy/optimize/_linprog_ip.py, _ip_hsd convergence criteria are met (go = rho_p > tol or rho_d > tol or rho_A > tol is False) and the algorithm terminates, but the "sanity check" for preventing false positives _check_result in scipy/optimize/_linprog_util.py reports failure. The sanity check is not related in an obvious way to the algorithm convergence criteria, and the tolerance used in the sanity check is a semi-arbitrary value that doesn't consider problem scaling, so the failure isn't surprising. I think the solution is for the _check_result to use a relative tolerance as noted in gh-9671. It could also be prevented by adding the sanity check to the algorithm's convergence criteria, but that might lead to other termination issues.
Reproducing code example:
This can be triggered with a minimally-modified version of the example in gh-6139.
from scipy import array, vstack, hstack
from scipy.optimize import linprog
from scipy.linalg import solve
f1 = array([[1., 0., 0.], [-1000., 0., - 1000.]])
g1 = array([5.00000000e+00, -1.00000000e+04])
f2 = -array([[0., 1000000., 1010000000.]]) # added some zeros to f2[2]
g2 = -array([10000000.])
sol = solve(vstack([f1, f2]), hstack([g1, g2-1]))
assert (f1@sol == g1).all()
assert (f2@sol < g2).all()
ret = linprog([1]*3, f2, g2, f1, g1, bounds=[(None, None)]*3)
print(ret)
Error message:
con: array([ 0.00000000e+00, -8.69476935e-07])
fun: -5030.000000000035
message: 'The solution does not satisfy the constraints within the required tolerance of 3.16E-04, yet no errors were raised and there is no certificate of infeasibility or unboundedness. This is known to occur if the `presolve` option is False and the problem is infeasible. This can also occur due to the limited accuracy of the `interior-point` method. Check whether the slack and constraint residuals are acceptable; if not, consider enabling presolve, reducing option `tol`, and/or using method `revised simplex`. If you encounter this message under different circumstances, please submit a bug report.'
nit: 5
slack: array([-0.87733746])
status: 4
success: False
x: array([ 5.00e+00, -5.04e+03, 5.00e+00])
Revised simplex has no trouble:
ret = linprog([1]*3, f2, g2, f1, g1, bounds=[(None, None)]*3, method='revised simplex')
results in:
con: array([0., 0.])
fun: -5030.0
message: 'Optimization terminated successfully.'
nit: 2
slack: array([0.])
status: 0
success: True
x: array([ 5.00e+00, -5.04e+03, 5.00e+00])
Scipy/Numpy/Python version information:
1.3.0 1.17.0 sys.version_info(major=3, minor=7, micro=3, releaselevel='final', serial=0)
But the same issue is likely experienced in more recent versions.