ENH: Provisional support for NEP 18 (__array_function__ protocol)

[neutrinoceros]

My initial motivation here was to add some unit representation to error messages when comparing two unyt_array instance via functions from numpy.testing, like so:

import numpy as np
import unyt as un

a = [1, 2, 3] * un.cm
b = [1, 2, 3] * un.km
np.testing.assert_array_equal(a, b)

which yields, on master:

...
AssertionError:
Arrays are not equal

Mismatched elements: 3 / 3 (100%)
Max absolute difference: 299997.
Max relative difference: 0.99999
 x: unyt_array([1, 2, 3])
 y: unyt_array([1, 2, 3])

and on this branch

previous version

...
AssertionError:
Arrays are not equal

Mismatched elements: 3 / 3 (100%)
Max absolute difference: 299997. cm
Max relative difference: 0.99999 dimensionless
 x: unyt_array([1, 2, 3] cm)
 y: unyt_array([1, 2, 3] km)

edit:

AssertionError:
Arrays are not equal

Mismatched elements: 3 / 3 (100%)
Max absolute difference: 299997., units='cm'
Max relative difference: 0.99999, units='dimensionless'
 x: unyt_array([1, 2, 3], units='cm')
 y: unyt_array([1, 2, 3], units='km')

Incidentally, it turns out that fixing this necessitated a kick-off implementation of NEP 18, so this work laid the foundation to solve:

• Using np.linalg.inv with unyt_array has no effect on units #69
• np.histogram does not support unyt_array #130
• Inconsistent behaviour between np.arange and np.linspace #50 (most likely out of scope)

More broadly, implementing NEP 18 for any was the topic of #139.
Granted I need to take more time to check that I'm not going against the original intentions here. My current approach is that, since covering the whole numpy public API in one go seems like an gigantic task, I'm implementing unyt_array.__array_function__ with a fallthrough condition: if a special case isn't implemented yet, just fallback to the raw numpy implem (which is currently the behaviour for all functions subject to NEP 18). This way we can add support for more and more functions in a progressive way.
I'm going to set the bar low(ish) for now, and try and fix the already reported cases, as reported above, as a first step.

An important question to address is: what should be done in the general case where we don't have a custom implementation for an array function ?

1. transparently default to the raw numpy implementation without a warning (this is de facto what is done as of unyt 2.8.0, and will remain the case until NEP 18 isn't at least partially implemented).
2. same as 1. but emit a warning (possibly, we could have a whitelist of functions that are known to be perfectly fine without a specialised implementation, for which no warning would be emitted) along the lines of

UserWarning: calling `numpy.FUNC` on a unyt_array. Results may hold incorrect units. A future version of unyt will remove this warning, and possibly change the behaviour of this function to be dimensionally correct.

3. Error out

Option 1 this is the current implementation of this PR because I think it is the less disruptive or noisy one. My personal opinion is that it's probably okay to have incomplete support for NEP 18 for a couple release, as long as it is clearly stated in the release notes.

[neutrinoceros]

I've gone the long way and started implemented a bunch of functions with obvious requirements.
TODO:

inverse-like functions

• np.linalg.inv
• np.linalg.tensorinv
• np.linalg.pinv

product-like functions

• np.dot
• np.vdot
• np.inner
• np.outer
• np.matmul
• np.kron

histogram functions (as requested for #130)

• np.histogram
• np.histogram2d
• np.histogramdd

Pifalls I found along the way

• np.linalg.pinv, np.matmul and likely other functions support "stacks of arrays" as their arguments. In this mode, NEP 18 doesn't help, as far as I can tell, so there is nothing we can do at the moment (as far as I understand), besides mentioning this limitation in docs, and possibly as a custom warning too.
• np.ma.dot isn't as trivial to implement than say np.dot, and the errors I got with my naive attempts were extremely similar o Creating masked array raises AttributeError #169, which was recently fixed by Fix unyt_array.__getitem__ to support numpy masked arrays #170. I think it's probably okay to leave this function out of the scope of the present PR.

[neutrinoceros]

pre-commit.ci autofix

[neutrinoceros]

@jzuhone do you think you'd be able to review this ? Is there someone else I could ask ?

[jzuhone]

@neutrinoceros I can try to look at it next week

[neutrinoceros]

@jzuhone I would like to propose we bundle this PR in the next release alongside the dask layer. IMO, each of these would justify a release on their own, so it's up to you :)

[jzuhone]

@chrishavlin or @matthewturk could you look at this one?

[neutrinoceros]

I rebased the branch but didn't write the history. I added two commit, the first of which corresponds to #288, and the second one adds deprecation warnings as requested.

[neutrinoceros]

pre-commit.ci autofix

[neutrinoceros]

(rebased again now that #288 went in)

[ngoldbaum]

Just a couple minor nits from the changes you made today. Will merge once they're resolved.