Skip to content

Implemented np.allclose in numba/np/arraymath.py #7940

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
sklam merged 16 commits into from
Apr 21, 2022
Merged

Implemented np.allclose in numba/np/arraymath.py #7940

Show file tree
Hide file tree
Changes from 10 commits
Commits
Show all changes
16 commits
Select commit Hold shift + click to select a range
e209bc1
Added np_allclose in numba/np/arraymath.py
czgdp1807 Mar 26, 2022
d16cd63
Added test for np_allclose
czgdp1807 Mar 26, 2022
eee11d0
Added broadcasting support and stronger tests
czgdp1807 Apr 1, 2022
2dc5bb7
Handeled scalars and improved tests
czgdp1807 Apr 2, 2022
ca98b63
Added test data from np.allclose tests
czgdp1807 Apr 2, 2022
6542301
Applied linting
czgdp1807 Apr 2, 2022
f09d8ff
Call a separate function to ensure same rank of arrays
czgdp1807 Apr 2, 2022
50367ef
Applied linting
czgdp1807 Apr 2, 2022
0cc9234
Merge branch 'main' into np_allclose
czgdp1807 Apr 2, 2022
9ee9842
Use np.broadcast_arrays and remove utility functions
czgdp1807 Apr 2, 2022
4d90833
Separate numpy tests for allclose
czgdp1807 Apr 7, 2022
ead6817
Improved tests of np_allclose
czgdp1807 Apr 8, 2022
9c1d44f
Corrected error messages in np_allclose
czgdp1807 Apr 8, 2022
1b24788
Added numpy.allclose to numpysupported.rst
czgdp1807 Apr 9, 2022
a93b0ef
Merge branch 'main' into np_allclose
czgdp1807 Apr 20, 2022
66969ce
Corrected argument number in error message
czgdp1807 Apr 20, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
94 changes: 94 additions & 0 deletions numba/np/arraymath.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -864,6 +864,100 @@ def flat_all(a):
return flat_all


@register_jitable
def _allclose_scalars(a_v, b_v, rtol=1e-05, atol=1e-08, equal_nan=False):
a_v_isnan = np.isnan(a_v)
b_v_isnan = np.isnan(b_v)

# only one of the values is NaN and the
# other is not.
if ( (not a_v_isnan and b_v_isnan) or
(a_v_isnan and not b_v_isnan) ):
return False

# either both of the values are NaN
# or both are numbers
if a_v_isnan and b_v_isnan:
if not equal_nan:
return False
else:
if np.isinf(a_v) or np.isinf(b_v):
return a_v == b_v

if np.abs(a_v - b_v) > atol + rtol * np.abs(b_v * 1.0):
return False

return True


@overload(np.allclose)
@overload_method(types.Array, "allclose")
def np_allclose(a, b, rtol=1e-05, atol=1e-08, equal_nan=False):

if not type_can_asarray(a):
raise TypeError('The first argument "a" must be array-like')

if not type_can_asarray(b):
raise TypeError('The first argument "b" must be array-like')

@sklam sklam Apr 14, 2022

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
raise TypeError('The first argument "b" must be array-like')
raise TypeError('The second argument "b" must be array-like')

@czgdp1807 czgdp1807 Apr 20, 2022

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have addressed the review and the tests are passing. Please let me know if anything else is to be done here. Thanks.


if not isinstance(rtol, types.Float):
raise TypeError('The fourth argument "rtol" must be a '

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

rtol is the third argument, right?

Suggested change
raise TypeError('The fourth argument "rtol" must be a '
raise TypeError('The third argument "rtol" must be a '

'floating point')

if not isinstance(atol, types.Float):
raise TypingError('The fifth argument "atol" must be a '

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
raise TypingError('The fifth argument "atol" must be a '
raise TypingError('The fourth argument "atol" must be a '

'floating point')

if not isinstance(equal_nan, types.Boolean):
raise TypeError('The sixth argument "equal_nan" must be a '

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
raise TypeError('The sixth argument "equal_nan" must be a '
raise TypeError('The fifth argument "equal_nan" must be a '

'boolean')

is_a_scalar = isinstance(a, types.Number)
is_b_scalar = isinstance(b, types.Number)

if is_a_scalar and is_b_scalar:
def np_allclose_impl_scalar_scalar(a, b, rtol=1e-05, atol=1e-08,
equal_nan=False):
return _allclose_scalars(a, b, rtol=rtol, atol=atol,
equal_nan=equal_nan)
return np_allclose_impl_scalar_scalar
elif is_a_scalar and not is_b_scalar:
def np_allclose_impl_scalar_array(a, b, rtol=1e-05, atol=1e-08,
equal_nan=False):
b = np.asarray(b)
for bv in np.nditer(b):
if not _allclose_scalars(a, bv.item(), rtol=rtol, atol=atol,
equal_nan=equal_nan):
return False
return True
return np_allclose_impl_scalar_array
Comment on lines +924 to +933

@guilhermeleobas guilhermeleobas Apr 6, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You can replace this block of code by:

Suggested change
elif is_a_scalar and not is_b_scalar:
def np_allclose_impl_scalar_array(a, b, rtol=1e-05, atol=1e-08,
equal_nan=False):
b = np.asarray(b)
for bv in np.nditer(b):
if not _allclose_scalars(a, bv.item(), rtol=rtol, atol=atol,
equal_nan=equal_nan):
return False
return True
return np_allclose_impl_scalar_array
elif is_a_scalar and not is_b_scalar:
def np_allclose_impl_scalar_array(a, b, rtol=1e-05, atol=1e-08,
equal_nan=False):
return np.allclose(b, a, rtol=rtol, atol=atol, equal_nan=equal_nan)
return np_allclose_impl_scalar_array

@czgdp1807 czgdp1807 Apr 7, 2022

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think this will give right results all the time, because the check depends on the ordering of a and b. For (a, b) the check is abs(a - b) > atol + rtol * abs(b), where as for (b, a) the check is abs(b - a) > atol + rtol * abs(a). See the RHS in both the inequalities, one uses b and the other uses a. So, np.allclose(a, b) might give different result from np.allclose(b, a) for edge cases.

Quoting from "Notes" section of https://numpy.org/doc/stable/reference/generated/numpy.allclose.html

The above equation is not symmetric in a and b, so that allclose(a, b) might be different from allclose(b, a) in some rare cases.

Please let me know if should still apply the above code suggestion.

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's fine. I wasn't aware of this info.

elif not is_a_scalar and is_b_scalar:
def np_allclose_impl_array_scalar(a, b, rtol=1e-05, atol=1e-08,
equal_nan=False):
a = np.asarray(a)
for av in np.nditer(a):
if not _allclose_scalars(av.item(), b, rtol=rtol, atol=atol,
equal_nan=equal_nan):
return False
return True
return np_allclose_impl_array_scalar
elif not is_a_scalar and not is_b_scalar:
def np_allclose_impl_array_array(a, b, rtol=1e-05, atol=1e-08,
equal_nan=False):
a = np.asarray(a)
b = np.asarray(b)
a_a, b_b = np.broadcast_arrays(a, b)

for av, bv in np.nditer((a_a, b_b)):
if not _allclose_scalars(av.item(), bv.item(), rtol=rtol,
atol=atol, equal_nan=equal_nan):
return False

return True

return np_allclose_impl_array_array


@overload(np.any)
@overload_method(types.Array, "any")
def np_any(a):
Expand Down
124 changes: 124 additions & 0 deletions numba/tests/test_np_functions.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -329,6 +329,10 @@ def np_trapz_x_dx(y, x, dx):
return np.trapz(y, x, dx)


def np_allclose(a, b, rtol=1e-05, atol=1e-08, equal_nan=False):
return np.allclose(a, b, rtol, atol, equal_nan)


def np_average(a, axis=None, weights=None):
return np.average(a, axis=axis, weights=weights)

Expand Down Expand Up @@ -3172,6 +3176,126 @@ def test_1D_weights_axis(data, axis, weights):
#test with axis argument
test_1D_weights_axis(data, axis=1, weights=w)

def test_allclose(self):

pyfunc = np_allclose
cfunc = jit(nopython=True)(pyfunc)

min_int = np.iinfo(np.int_).min
a = np.array([min_int], dtype=np.int_)

simple_data = [
(np.asarray([1e10, 1e-7]), np.asarray([1.00001e10, 1e-8])),
(np.asarray([1e10, 1e-8]), np.asarray([1.00001e10, 1e-9])),
(np.asarray([1e10, 1e-8]), np.asarray([1.0001e10, 1e-9])),
(np.asarray([1e10]), np.asarray([1.0001e10, 1e-9])),
(1.0, 1.0),
(np.array([np.inf, 1]), np.array([0, np.inf])),
(a, a)
]

for a, b in simple_data:
py_result = pyfunc(a, b)
c_result = cfunc(a, b)
self.assertEqual(py_result, c_result)

a = np.asarray([1.0, np.nan])
b = np.asarray([1.0, np.nan])
self.assertFalse(cfunc(a, b))
self.assertTrue(cfunc(a, b, equal_nan=True))

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Compare the output of cfunc with pyfunc instead.


b = np.asarray([np.nan, 1.0])
self.assertFalse(cfunc(a, b))

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Compare the output of cfunc with pyfunc instead.


# NumPy test data

@guilhermeleobas guilhermeleobas Apr 6, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can you also include the link to NumPy repository pointing to the tests?

@guilhermeleobas guilhermeleobas Apr 6, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe split the tests into different functions. This makes it easy to debug when there's a test failure.

arr = np.array([100.0, 1000.0])
aran = np.arange(125).astype(dtype=np.float64).reshape((5, 5, 5))
atol = 1e-8
rtol = 1e-5
numpy_data = [
(np.asarray([1, 0]), np.asarray([1, 0])),
(np.asarray([atol]), np.asarray([0.0])),
(np.asarray([1.0]), np.asarray([1 + rtol + atol])),
(arr, arr + arr * rtol),
(arr, arr + arr * rtol + atol * 2),
(aran, aran + aran * rtol),
(np.inf, np.inf),
(np.inf, np.asarray([np.inf]))
]

for (x, y) in numpy_data:
self.assertTrue(cfunc(x, y))

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Instead of asserting to True, compare the outputs of cfunc(x, y) and pyfunc(x, y) are equals.


numpy_data = [
(np.asarray([np.inf, 0]), np.asarray([1.0, np.inf])),
(np.asarray([np.inf, 0]), np.asarray([1.0, 0])),
(np.asarray([np.inf, np.inf]), np.asarray([1.0, np.inf])),
(np.asarray([np.inf, np.inf]), np.asarray([1.0, 0.0])),
(np.asarray([-np.inf, 0.0]), np.asarray([np.inf, 0.0])),
(np.asarray([np.nan, 0.0]), np.asarray([np.nan, 0.0])),
(np.asarray([atol * 2]), np.asarray([0.0])),
(np.asarray([1.0]), np.asarray([1 + rtol + atol * 2])),
(aran, aran + aran * atol + atol * 2),
(np.array([np.inf, 1.0]), np.array([0.0, np.inf]))
]

for (x, y) in numpy_data:
self.assertFalse(cfunc(x, y))

@guilhermeleobas guilhermeleobas Apr 7, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same here. Compare the outputs of cfunc(x, y) and pyfunc(x, y).


class Foo(np.ndarray):
def __new__(cls, *args, **kwargs):
return np.array(*args, **kwargs).view(cls)

a = Foo([1])
self.assertTrue(type(np.allclose(a, a)) is bool)

noise_levels = [1.0, 1e-1, 1e-2, 1e-3, 1e-4, 1e-5, 1e-6, 0.0]
zero_array = np.zeros((25, 4))
a = np.random.ranf((25, 4))
for noise in noise_levels:
for rtol in noise_levels:
for atol in noise_levels:
py_result = pyfunc(zero_array, noise,
atol=atol, rtol=rtol)
c_result = cfunc(zero_array, noise,
atol=atol, rtol=rtol)
self.assertEqual(py_result, c_result)

py_result = pyfunc(noise, zero_array,
atol=atol, rtol=rtol)
c_result = cfunc(noise, zero_array,
atol=atol, rtol=rtol)
self.assertEqual(py_result, c_result)

py_result = pyfunc(np.asarray([noise]), zero_array,
atol=atol, rtol=rtol)
c_result = cfunc(np.asarray([noise]), zero_array,
atol=atol, rtol=rtol)
self.assertEqual(py_result, c_result)

py_result = pyfunc(a, a + noise, atol=atol, rtol=rtol)
c_result = cfunc(a, a + noise, atol=atol, rtol=rtol)
self.assertEqual(py_result, c_result)

py_result = pyfunc(a + noise, a, atol=atol, rtol=rtol)
c_result = cfunc(a + noise, a, atol=atol, rtol=rtol)
self.assertEqual(py_result, c_result)

def test_allcase_exception(self):

@guilhermeleobas guilhermeleobas Apr 6, 2022

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
def test_allcase_exception(self):
def test_allclose_exception(self):

self.disable_leak_check()

pyfunc = np_allclose
cfunc = jit(nopython=True)(pyfunc)

a = np.asarray([1e10, 1e-9, np.nan])
b = np.asarray([1.0001e10, 1e-9])

with self.assertRaises(ValueError) as e:
cfunc(a, b)

self.assertIn(("shape mismatch: objects cannot be broadcast to "
"a single shape"), str(e.exception))

def test_interp_basic(self):
pyfunc = interp
cfunc = jit(nopython=True)(pyfunc)
Expand Down