I like this idea. Following the broadcasting suggestion, it may be nice to try to be reasonably compatible with whatufuncs do generally. I guess that might include an optionalout argument (and correspondingout attribute) or perhaps instead a methodoutput_arrays(out=None) that would create the array with the right dtype and subclass, or check if it is OK if an array is given. Taking this perhaps too far, one could similarly have anaxis argument/method as well (which would just usenormalize_axis on the broadcasted shape under the hood).

p.s. Might it be possible to factor out theufunc machinery (rather that write this from scratch)? Or indeed use something like this inside theufuncs?

I think@mdhaber was hoping to include broadcasting. Not sure, since it is one call away, but also seems fine to add.

That would be nice. It seems like a common enough follow-up step that it would be a convenient addition to this convenience function.

So, is it insane to create a rich object...

We do this frequently in SciPy now, so not insane. Seems unusual for NumPy, but I wouldn't complain.

3. Idea: maybe instead of an extra parameterensure_inexact=True, couldextra_dtype accept classes likenp.inexact, if that makes sense?

4. I'm not used to reading C. Would it be easy to add a docstring or the equivalent sequence of existing NumPy commands so I can check my understanding of what is going on?

5. We had also discussed the possibility of addingaxis to be ignored whenbroadcasting for use in, say, multi-sample reduction statistics. For example, inttest_ind, it is sensible forx.shape == (5,) andy.shape == (3, 4) withaxis=-1; when the arrays are broadcasted, the shapes become(3, 5) and(3, 4). I implemented this already in SciPy, but should I open an issue about it being moved upstream?

Thanks for this!

Yeah, agree that the ufunc machinery should share as much as possible, although was a bit hoping to mainly consolidate Python paths first. The reason was that looking at the array-wrap implementation in ufuncs they are messy and buggy (there is an issue), so fixing that seems like a bit of a messy change easier on its own.

I think a method for output array prep could make sense, adding it to the initial call seems like it would overload meaning (we don't know the common dtype is right for one thing).

Idea: maybe instead of an extra parameter ensure_inexact=True, could extra_dtype accept classes like np.inexact, if that makes sense?

Right, that is how it is implemented. Justdtype=np.inexact doesn't quite work yet and I thought theresult_type(..., 0.) pattern is so common, it might actually make sense to have it explicitly.

5. [discussed] adding axis to be ignored when broadcasting for use in, say, multi-sample reduction statistics [...]

Ah, right! Yeah, maybe make its own issue, it seems to me we should with integrating that intonp.broadcast_arrays first. Admittedly, I think that may require reimplementing the whole function (maybe for better).

@seberg from your description it is not entirely clear to me if you mean for this to be a public or private function? With the namespace we discussed somewhere else,np.lib.array_utils I believe, we should get a place to put utilities like this - but I can't tell from the description yet if this would be useful more widely.

Since it was mentioned above, we recently came across the need fornormalize_axis_index in SciPy - that's also a useful utility, which is not public right now.

The first issue I want to solve would mean private is fine, but yes suspect we should just make it public and if we are unsure about some methods/attributes or so give it an_ (which would mean putting it intonp.lib.array_utils namespace, or different name if anyone suggests a better one soon).
(And@mhvk's interest makes me lean more towards public at last if I can get some of what@mdhaber wants too.)

Agreed that the output makes more sense as a separate method, rather than on input, especially since in most cases something has to be generated.

Also, 👍 to makingnormalize_axis_index a publicly accessible tool (we already import it...); probably all ofshape_base can go into the same collection too.

p.s. A better name thancreation_helper would still be good. Will try to think a bit more.

p.p.s. You'll not be surprised to hear that I think that by default subclasses should be passed through... But it may be good to think a bit more generally what to do with array mimics - the benefit of gettingresult_type depends only on having adtype, not on being anndarray (subclass) (and similarly for a possible broadcast shape). Though easiest perhaps for that is just to let this be overridden by__array_function__.

I just ran into this little gem , which might have been avoided with use of a helper like this.

np.float16([2])<np.float32(2.0009766)# array([False])np.float16([2])<np.float32([2.0009766])# array([ True])

You might think that helps, but it probably wouldn't :). This is NEP 50 stuff.

If before comparison the float16 is promoted to the result type or if the float32 is broadcasted to the result shape, the comparison is done correctly. So I think the helper would help with those things.

@mdhaber no,np.result_type() would do the samedown promotion in the above case until you do NEP 50 .

I always pass the dtypes intoresult_type instead of the arrays.

import numpy as npx = np.float16([2])y = np.float32(2.0009766)print(x < y)x1, y1 = np.broadcast_arrays(x , y)print(x1 < y1)dtype = np.result_type(x.dtype, y.dtype)print(x.astype(dtype) < y.astype(y))

[False][ True][ True]

Right, if you pass the dtypes you get the right thing, but that is the same as adopting NEP 50.

I was trying to get back to this a bit, but having difficulties to quite formulate what we need/want here. I would prefer to not go the feature-creep route for starters.

Things that I would like:

• Allow converting multiple arrays at once more conveniently (easy).
• The result should havecommon_dtype() as it is pretty useful especially with NEP 50.
• The result should have awrap() to apply__array_wrap__() when it makes sense.¹
• We could have more options:
  • shape (braodcast shape)
  • Abroadcast() method to modify the "arrays"` (probably needs a local implementation for starters, which then might be adoptable elsewhere). Broadcastingcould be a kwarg also.

There are some options, which would probably fit as kwargs² :

• subok=True (to set it toFalse)
• broadcastingcould be a kwarg.

Another API, I worry about more: it would be nice to provide a clean way to preserve Python better. This is to allow later normal promotion to use NEP 50 rules.
Due to theisclose and similar use-case, it also seems to me that in many cases it would make sense for broadcasting tonot modify scalars (or 0-D objects).
For that, it would make sense to either have a kwargpreserve_scalars=True or an attributearrays_or_scalars?

Maybe one solution is to make it:

conv=array_converter(*arrays)dt=conv.common_dtype()# always applicable, may do initial conversion to arraysconv.shape# could do this to get the broadcast shape, but could also defer.conv.was_pyscalar# maybe.  Or an enum with `0-D, scalar, python scalar`.arrays=conv.arrays(subok=True,allow_scalars=True,broadcast=True,dtype=dt)# Wrap result (if wanted) (see footnote 2) [^2]returnconv.wrap(result)

Any thoughts? Implementing any of this isn't super hard, but getting more clarity on where to start without a feature creap, is for me.

I think, I like that last thought as API (initialization takes no kwargs, but fetching the arrays is a function and not an attribute/property). But it isn't quite clear to me. Especially how to do scalar preservation.

p.s. On theasanyarray, beyond the rationale that subclasses should behave the same way and thus not be converted (of course, they don't always do...), it is certainly nice if generally they get passed on, since then__array_ufunc__ coverage often is good enough for functions to work. ForQuantity at least, there are still a lot of functions where we do not have to override; seehttps://github.com/astropy/astropy/blob/main/astropy/units/quantity_helper/function_helpers.py (theSUBCLASS_SAFE_FUNCTIONS set; as you'll notice, there are also a lot that do need overrides, which means a very long list of rather boilerplate code that is little fun to write, so all the better if some functions don't need it!)

p.s. Probably time for some squashing of commits too...

I don't want to jinx this because of it, but the[()] logic for subclasses doesn't work for currentmemmap (which incorrectly always does it).
Now that is maybe even never used right now, but it makes me wonder if we need another approach, which might mean passingscalar=True in__array_wrap__!?

I'm a bit confused: why would it not work formemmap? It doesn't hurt to do[()] twice, does it? (I do agree we can and perhaps should adjustmemmap.__array_wrap__ to not create scalars itself, but that seems independent of this PR).

That said, you are right that in principle it would be better to give__array_wrap__ a chance to react to whether or not a scalar is expected, instead of using a heuristic as you have now. And I guess 2.0 is the right time for such an API change. Though on the other hand,__array_wrap__ still is somewhat out-of-date API and generally classes have to be able to deal with[()] anyway. So, perhaps we do not need to worry?

see#25278 which was reported with impressive timing for this PR

I'm still not quite sold onreturn_scalar, but I'm also realizing I've lost what little overview I had of the problems that are caused by the various steps. This makes me wonder if it would be better to have a PR first that creates a helper that simply keeps the status quo, i.e., makes no attempt to be clever on wrapping, just reproducing the__array_wrap__ calls that are present now.

@mhvk and@ngoldbaum, I merged main (i.e. the array-wrap helper) and applied all the remaining doc suggestions (thanks, they were great!). Then squashed it down to three commits.

Ithink this should be pretty much ready now, since IIRC it was relatively settled. But happy to do another round if things come up of course.