[Numpy-discussion] Allowing broadcasting of code dimensions in generalized ufuncs (original) (raw)

Stephan Hoyer shoyer at gmail.com
Fri Jun 15 14:17:09 EDT 2018

On Mon, Jun 11, 2018 at 11:59 PM Eric Wieser <wieser.eric+numpy at gmail.com> wrote:

I don’t understand your alternative here. If we overload np.matmul using arrayfunction, then it would not use ether of these options for writing the operation in terms of other gufuncs. It would simply look for an arrayfunction attribute, and call that method instead.

Let me explain that suggestion a little more clearly. 1. There’d be a linalg.matmul2d that performs the real matrix case, which would be easy to make as a ufunc right now. 2. matmul and rmatmul would just call np.matmul, as they currently do (for consistency between np.matmul and operator.matmul, needed in python pre- at -operator) 3. np.matmul would be implemented as: @doarrayfunctionoverridesdef matmul(a, b): if a.ndim != 1 and b.ndim != 1: return matmul2d(a, b) elif a.ndim != 1: return matmul2d(a, b[:,None])[...,0] elif b.ndim != 1: return matmul2d(a[None,:], b) else: # this one probably deserves its own ufunf return matmul2d(a[None,:], b[:,None])[0,0] 4. Quantity can just override arrayufunc as with any other ufunc 5. DataArray, knowing the above doesn’t work, would implement something like @matmul.registerarrayfunction(DataArray)def arrayfunction(a, b): if a.ndim != 1 and b.ndim != 1: return matmul2d(a, b) else: # either: # - add/remove dummy dimensions in a dataarray-specific way # - downcast to ndarray and do the dimension juggling there

Advantages of this approach: - Neither the ufunc machinery, nor arrayufunc, nor the inner loop, need to know about optional dimensions. - We get a matmul2d ufunc, that all subclasses support out of the box if they support matmul Eric OK, this sounds pretty reasonable to me -- assuming we manage to figure out the array_function proposal!

There's one additional ingredient we would need to make this work well: some way to guarantee that "ndim" and indexing operations are available without casting to a base numpy array.

For now, np.asanyarray() would probably suffice, but that isn't quite right (e.g., this would fail for np.matrix).

In the long term, I think we need a new coercion protocol for "duck" arrays. Nathaniel Smith and I started writing a NEP on this, but it isn't quite ready yet.

​ -------------- next part -------------- An HTML attachment was scrubbed... URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20180615/c5b747f8/attachment.html>

More information about the NumPy-Discussion mailing list