isZero() does not respect config.epsilon

[gwhitney]

This concern came up in the implementation of the polynomial root solver that I am doing as a basic algebraic manipulation benchmark for mathjs/picomath. I was using isZero to test the discriminant of quadratics/cubics, but it involves some cancellation, so it can come out to a minuscule technically non-zero floating point number. I had assumed that isZero was a "relational function" as described in the mathjs docs and hence would consider such a result to be zero. However, that's not currently the case, leading to erroneous results (two distinct roots when there really should be one, or complex roots when they should be real, for example). So should we:

(a) Modify isZero to respect config.epsilon? (Clearly a breaking change...)
(b) Add an equalZero function? (Right now one can use equal(x, 0) but that leads to an unnecessarily complicated double dispatch.)

[gwhitney]

Oh, actually all that equal(x, 0) tests is whether x is less than Number.EPSILON (the smallest number you can add to one that changes it). Not sure why that makes sense, but right now it has absolutely nothing to do with config.epsilon. Since this is a relative comparison, perhaps no x should be equal(x, 0) other than 0 itself, because the difference is always equal to x, and hence not small relative to x. I also don't understand the logic: if two numbers are within Number.EPSILON of each other, they are considered equal. But 1e-19 and 1e-20 are quite distinct (the difference is 90% of the larger) and floating point numbers have no trouble distinguishing them, so considering them to be relatively equal doesn't make any sense. There is a line in the documentation about how equal() can't be used to compare numbers smaller than Number.EPSILON, and it is the case, but only because of the convention that two numbers within Number.EPSILON are always equal regardless of config.epsilon, violating the "relative" aspect of the comparison. I think this violation comes from a desire for it to be the case that if equal(a,b), then equal(a-b, 0) -- but that's not even true right now, and it doesn't actually make sense for relative comparison.

So my thinking at the moment: This violation is unnecessary and the situation is buggy. Change equal to truly be relative, so that the only thing equal to 0 is actual 0. Perhaps add an isNegligible function to mean "is less than config.epsilon in absolute value" -- that may need to be used in some places that relying on the "within Number.EPSILON is equal" fallback.

But I didn't want to file a bug until I could get at least a second opinion. Thoughts?

[josdejong]

Good point. It would be logic indeed if isZero would be consistent with equal and utilizes config.epsilon. I'm open to discuss ways to improve this!

For reference: the current implementation of nearlyEqual used by the relational functions is a mix of an absolute comparison and a relative comparison:

mathjs/src/utils/number.js

Lines 598 to 636 in 77f94c4

	/**
	* Compares two floating point numbers.
	* @param {number} x First value to compare
	* @param {number} y Second value to compare
	* @param {number} [epsilon] The maximum relative difference between x and y
	* If epsilon is undefined or null, the function will
	* test whether x and y are exactly equal.
	* @return {boolean} whether the two numbers are nearly equal
	*/
	export function nearlyEqual (x, y, epsilon) {
	// if epsilon is null or undefined, test whether x and y are exactly equal
	if (epsilon === null || epsilon === undefined) {
	return x === y
	}
	if (x === y) {
	return true
	}
	// NaN
	if (isNaN(x) || isNaN(y)) {
	return false
	}
	// at this point x and y should be finite
	if (isFinite(x) && isFinite(y)) {
	// check numbers are very close, needed when comparing numbers near zero
	const diff = Math.abs(x - y)
	if (diff < DBL_EPSILON) {
	return true
	} else {
	// use relative error
	return diff <= Math.max(Math.abs(x), Math.abs(y)) * epsilon
	}
	}
	// Infinite and Number or negative Infinite and positive Infinite cases
	return false
	}

Unfortunately I didn't put any reference on the formulas and where I found them, I remember I did some reading up on comparing numbers and that this can be challenging and some background reading can help understanding the problem. Calculating the relative difference works very well for any number except zero, since the comparison can only do an absolute comparison there, not a relative one.

So I think we do have two open ideas here:

1. Adjust isZero to utilize epsilon and become consistent with equal. I would say yes lets do that.
2. Can we adjust equal to work better with 0, and make equal(a-b, 0) work. I'm not sure about how to make this work though. In Customizable behavior for nearlyEqual: relative or absolute difference #2752 I mention a workaround , but I think this has problems of it's own.

See also: #2752, #1614

[gwhitney]

the current implementation of nearlyEqual used by the relational functions is a mix of an absolute comparison and a relative comparison

Gotcha. The documentation is not clear on this, it just says it is relative comparison, so at least documentation and reality ought to be synced.

Adjust isZero to utilize epsilon and become consistent with equal. I would say yes lets do that.

Well, but equal(x, 0) makes no use of config.epsilon at the moment. It uses Number.EPSILON (and is not at the moment configurable). We could change isZero to use Number.EPSILON as well to make them consistent, but it is not clear why Number.EPSILON is relevant to comparing to zero, since it has to do with numbers near one, not near zero.

make equal(a-b, 0) work
If someone does

let x = a - b
if (math.equal(x, 0)) ...

and expects that to work like

if (math.equal(a, b)) ...

then I think they are simply out of luck in a world where we are doing relative comparisons, since all scale information is lost by the time you get to equal(x,0).

In light of all of the above, here are some measures we could take that might help:

Either or both of:

1. Add an isNeglible function that just tells you if something is less in absolute value than config.epsilon.
2. Change nearlyEqual to be purely relative, no absolute component, no involvement of Number.EPSILON. In this case, equal(x, 0) will come to match isZero(x) and be true only when x is precisely zero. (That behavior for zero may admittedly be unintuitive for casual users.)

Alternatively:

• Overhaul the approximate comparison system for one that is deliberately and configurably a mix of absolute and relative comparison, with separate relative_epsilon and absolute_epsilon settings. Two numbers would be considered equal if they are equal either because of the relative_epsilon or absolute_epsilon. Setting relative_epsilon to 0 would give pure absolute-difference comparison, as the poster in Customizable behavior for nearlyEqual: relative or absolute difference #2752 seemed to want, and conversely setting absolute_epsilon to 0 would give pure relative-difference comparison, as you likely want if you are dealing with small quantities. Defaulting relative_epsilon to the current default for config.epsilon and absolute_epsilon to Number.EPSILON would more or less exactly recover the current behavior, which does seem to generally work basically OK in a rough-and-ready way for everyday purposes.

In either case perhaps there should be a way to do some comparisons with different epsilon values than the current config? Or is it easy and not computationally expensive to just set config values to something else during a subcomputation and then set them back?

Either route should be combined with carefully documenting the situation.

I guess at the moment I'd lean toward the latter alternative because it involves no new functions, seems to give everyone what they want, and is not a breaking change (except for changing the name config.epsilon to config.relative_epsilon; we could temporarily make config.epsilon an alias for config.relative_epsilon that issues a deprecation warning, and at a future breaking change make it an error to access config.epsilon).

[josdejong]

We definitely have to at least improve the documentation about it 😅 .

Well, but equal(x, 0) makes no use of config.epsilon at the moment

O wow, now I see. I was at first assuming it was also using config.epsilon but it is indeed using Number.EPSILON. Anyway, when removing this absolute comparison with DBL_EPSILON, some of the following unit tests fail. Looking at the FIXME's these tests and logic can use some love anyway:

mathjs/test/unit-tests/utils/number.test.js

Lines 546 to 559 in b050d5e

	it('should test whether two small numbers are nearly equal (always true)', function () {
	const epsilon = 1e-2
	assert.strictEqual(nearlyEqual(1e-200, 0.99e-200, epsilon), true)
	assert.strictEqual(nearlyEqual(1e-200, 10e-200, epsilon), true) // FIXME: why is this true?
	})
	it('should compare with zero', function () {
	const epsilon = 1e-3
	assert.strictEqual(nearlyEqual(0, 0, epsilon), true)
	assert.strictEqual(nearlyEqual(0, -0, epsilon), true)
	assert.strictEqual(nearlyEqual(0, 1.2, epsilon), false)
	assert.strictEqual(nearlyEqual(0, 1e30, epsilon), false)
	assert.strictEqual(nearlyEqual(0, 1e-30, epsilon), true) // FIXME: why is this true?
	})

1.

Yes makes sense to me to introduce two variants of isZero, a "nearly equal" variant and an "exact" variant. I have the feeling that it would be most consistent to make isZero part of the relational functions, using epsilon just like all the other functions like equal and smaller. So, maybe we can do that, an if needed (depending on the decisions in (2) and (3)), create an isStrictlyZero() variant.

2.

The current behavior of comparing to zero is weird anyway, so I think it is a good idea to improve that. It makes sense to me to introduce a relativeEpsilon (the current epsilon) and an absoluteEpsilon. For backward compatibility we should set absoluteEpsilon = Number.EPSILON, but I think in general a better setting would be absoluteEpsilon = 1e-12. I think working with "regular" values and wanting to compare equal(a-b,0) is much more common in real world usage than working with extremely small values like 1.4e-21, so maybe best to cater for that. What do you think?

3.

The function isZero is not the only function that we may want to change to use epsilon. We also have: isPositive, isNegative, and isInteger. I think those functions should ideally use epsilon too, and become a "real" part of the relational functions. Maybe we should also move them from /utils to /function/relational.

[dvd101x]

I did a quick run of changing nearlyEqual and it looks promising. The main difference is that instead of taking epsilon as an argument, it would need atol and rtol it could be with the values previously discussed. It runs but breakers some tests.

Maybe a function should be created as it makes sense also for matrices or units as input.

As a side note solveODE could benefit of that function as currently it only takes atol

[josdejong]

Sounds good!

[gwhitney]

I definitely think that absoluteTolerance and relativeTolerance (or whatever we want the "official" names to be) should become configuration parameters, with default values that minimize the number of breaking unit tests. The ones that do break will need to be examined by hand and decide whether to ship changes before or along with any PR that implements this change. If I understand correctly, epsilon then becomes an unused parameter; or is it still used somewhere else? If not, we may still want to deprecate but allow changes to the epsilon configuration, which would result in roughly corresponding changes to absoluteTolerance and relativeTolerance, for one release cycle to allow clients of mathjs to adjust to the change.

[dvd101x]

Sounds good. I will gladly get the ball rolling over the weekend for all to review the failing tests.

[josdejong]

@glen I think epsilon should become deprecated in the end, though a deprecation period + backward compatibility would be nice (if possible).