replace np.isclose for pivot selection warning

[pmeier]

Reference issue

Closes #14081

What does this implement/fix?

This replaces numpy.isclose in

scipy/scipy/optimize/_linprog_simplex.py

Lines 219 to 229 in 0617428

	# The selected pivot should never lead to a pivot value less than the tol.
	if np.isclose(pivval, tol, atol=0, rtol=1e4):
	message = (
	"The pivot operation produces a pivot value of:{0: .1e}, "
	"which is only slightly greater than the specified "
	"tolerance{1: .1e}. This may lead to issues regarding the "
	"numerical stability of the simplex method. "
	"Removing redundant constraints, changing the pivot strategy "
	"via Bland's rule or increasing the tolerance may "
	"help reduce the issue.".format(pivval, tol))
	warn(message, OptimizeWarning, stacklevel=5)

since it relies on the asymmetric definition. Using np.isclose(tol, pivval) leads to a lot of false positives, since rtol is applied to the second argument.

The change communicates the intent better and might also be safer in edge-cases.

[mdhaber]

Thanks @pmeier. Can you provide some example cases of the false positives and add a test that this resolves the problem?

On the other hand, method='simplex' is only really retained for backwards compatibility, so we're not doing much enhancement these days. (Please see the newer HiGHS methods.) If @mckib2 @Kai-Striega are interested in reviewing this, sounds fine, but I might not take a look for a while.

[pmeier]

Can you provide some example cases of the false positives and add a test that this resolves the problem?

There are no false positives yet. There is currently a discussion in numpy about changing the asymmetric definition of np.isclose

abs(actual - expected) <= atol + rtol * abs(expected)

to a symmetric one (following math.isclose)

abs(actual - expected) <= max(atol, rtol * max(abs(actual), abs(expected)))

To judge the impact of such a change I ran a patched numpy version with the above change against the scipy test suite. For example

python runtests.py \
    -t scipy.optimize.tests \
    -- \
    -k "Simplex and test_bounds_mixed"

will give warnings (and in turn test errors) like this:

scipy.optimize.optimize.OptimizeWarning: The pivot operation produces a pivot value of: 2.0e+00, which is only slightly greater than the specified tolerance 1.0e-09. This may lead to issues regarding the numerical stability of the simplex method. Removing redundant constraints, changing the pivot strategy via Bland's rule or increasing the tolerance may help reduce the issue.

Since 2 is far away from being close to 1e-9, I consider that a false positive. This is a result of relying on the asymmetric definition of (unpatched) np.isclose

pivval = 2
tol = 1e-9
atol = 0
rtol = 1e4

np.isclose(pivval, tol, atol=atol, rtol=rtol)  # False
np.isclose(tol, pivval, atol=atol, rtol=rtol)  # True

Since the intent seems to be "check if pivval is not significantly greater than the given tol", IMO it would be better to adapt the code to reflect that

abs(pivval) <= tol * 1e4

Note that this also circumvents checking the difference of pivval and tol, which might lead to other false positives in case they are both in the same magnitude.

---

On the other hand, method='simplex' is only really retained for backwards compatibility, so we're not doing much enhancement these days.

This is more future proofing this than enhancing the functionality.

[mckib2]

Like @mdhaber mentioned, we're in maintenance-only mode for simplex and recommend the newer HiGHS methods. However, I see no issue with accepting simple patches for simplex, so if the CI looks good I'll merge

[mdhaber]

Thanks @pmeier but I doubt np.isclose is going to change. If it does, I suspect there will be a lot of other problems to fix at the same time. Would it be OK with you if we close this until something changes? I suspect they would give downstream projects some time to respond, and we could re-open then if need be.

[pmeier]

@mckib2

However, I see no issue with accepting simple patches for simplex, so if the CI looks good I'll merge

CI seems to be not happy, although I can't find a reason why.

---

@mdhaber

I doubt np.isclose is going to change.

I thought so too, but what I get from the discussion following my comment in numpy/numpy#10161 is that this is actually on the table.

I suspect they would give downstream projects some time to respond

Most definitively. The current state is to create change it in a separate namespace and slowly trickle it down to the core module.

I suspect there will be a lot of other problems to fix at the same time.

It is actually not that bad. If we only change the definition from asymmetric to symmetric, in scipy's test suite there are < 10 test that would fail not counting the errors from the raised warnings I'm trying to fix here. 3 of them should probably not pass in the first place, since they rely on the fact that atol and rtol are added together, which in some cases doubles the effective tolerance.

Would it be OK with you if we close this until something changes?

Normally, I would disagree, since my patch makes it safer, easier to read, and even future proof. But since the CI seems not to be happy about it, I don't know if it is worth the effort to debug. I'll leave it up to you to close this.

[mckib2]

Closing for now, we can re-evaluate if the proposals for np.isclose and friends moves forward