DOC: NEP: array in -> array out

[mdhaber]

Draft NEP. Please review here before we post to NumPy.
@mihaimaruseac @eriknw @sjsrey @matthewfeickert @ksunden @virchan @j-bowhay @dschult @pfackeldey @seberg

[mdhaber]

LMK if you want your name added.

[mdhaber]

One exception is axis=None -> scalar.
If I understood correctly, indexing an array like:

x = np.asarray([1, 2, 3])
x[0]

produces a scalar, too?

axis=None is one thing, but this would make selling the proposal even harder, since it's a half-measure.

[seberg]

since it's a half-measure.

It makes it much easier because it avoids a swaths of BC issues since arr1d[0] will in practice be one of the most common paths to get scalars in user code.
(While it will actually be not that common in major libraries.)

You can call it a half-measure if you like, but I don't see it as such at all. The point is x[0, ...] or x[..., 0] is a simple way to get what you want already.

To home in on what IMO, is the core of all of this: Make scalar returns predictable and avoidable when working with N-D arrays.
Yes, scalars have more annoyances, but I think focusing on the above makes it easier to see that this is the minimal way to achieve it and avoiding the temptation to change things that really don't need changing.

[mdhaber]

Make scalar returns predictable and avoidable

I don't see how to sell that, though. Individual NumPy operations are pretty predictable right now - operations other than array creation functions will always return a scalar instead of a 0d array. There may be exceptions to that rule, but doesn't that sum it up, or is there a ton of inconsistencies that this PR resolves? (Can you provide examples of surprises we get right now that are not summed up by the rule above?)

[seberg]

You are right, the problem isn't really them being unpredictable, the problem is that there is no good way to avoid the 0-D array -> scalar conversion.
The issue is preserving the 0-D arrays to be able to work with N-D arrays without hassle.

[mdhaber]

there is no good way to avoid the 0-D array -> scalar conversion.

OK, so do I understand that somewhere in the code, the results of ufuncs and reducing operations is a 0d array, and then it gets converted (unavoidably) to a scalar. You are reversing that choice except for reducing operations when axis is None.

When we fully index an array (e.g. one integer for a 1d array, 2 integers for a 2d array), is the result naturally a 0-d array or a scalar? Is it like ufuncs in that the result is naturally 0d and it gets converted to a scalar, or is it the other way arround - the result is naturally a scalar and if we wanted it to return a 0d array, it would not be so different from just calling asarray on the result?

[mdhaber]

Should I link https://github.com/seberg/numpy/tree/try-0d-preservation-rebase? (Or will there be a different branch, do you think?)

[mdhaber]

@j-bowhay @pfackeldey @ksunden can you provide links to your branches?

[pfackeldey]

This is the link to Awkward Array with the NumPy branch provided by @seberg: https://github.com/scikit-hep/awkward/tree/pfackeldey/test_NEP57. There are only 2 failures and they are related to pandas. For Awkward Array the migration seems to be simple.

[j-bowhay]

Just catching up on notification now, happy to provide a branch for SciPy and read over this but will probably need a week or two to do so

[j-bowhay]

What's the status here? I am still happy to work on a branch but don't want to invest the time if things are blocked elsewhere

[seberg]

I have worked on it here: numpy#29067 but there are still some issues around polynomials mostly I think, where the NumPy python code needs to properly deal with the changes (which is a bit tedious).

[mdhaber]

It's ok to wait, @j-bowhay. We'll ping you when things are closer.

[mdhaber]

@seberg should I propose something? What would you suggest?

[seberg]

Propose something. I would say:

1. Create a branch and use it to adapt and test downstream.
2. The final NEP acceptance will happen after this with potentially more discussion.

Or in other words, I would just not say much about the precise path of achieving step 2. That makes sense IMO, because without 1 it is even harder to judge how big of a change it is (not that 1 will be a silver bullet).

How exactly you want to write it, I don't care. You could even add a .. note:: early on that a final decision is deferred until after further testing.

[mdhaber]

Please replace the sentences after the first for me : )

[seberg]

Minimal change, minimal change, minimal change. And maybe also why not?

If you change arr.sum() to return 0-D arrays you will create a class of places where you return 0-D arrays for, honestly, unclear reasons and you will break user code. Not as much as arr1d[0], but still...
If you change arr1d.sum(-1) I think that is much less, because who would write that in a non N-D context?!

I keep repeating myself, but the point is making the return predictable for N-D arrays, changing arr.sum(), frankly, to me smells like "I don't like scalars", and not "I want to work with N-D arrays in a predictable way".

Now, "I don't like scalars" is a valid opinion, but I think we established well enough that we won't achieve that easily, so I don't see why we should try to make inroads towards it that we don't need to fix the actual core issue (don't get random scalars when working with N-D arrays).

EDIT: To be clear, I am not in love with this rule, also because it is slightly different from the rule for indexing (where adding ... is the magic). But it seems like a rule that should have a very good hit-to-miss ratio in existing code and I suppose I also can't think of a neat rule that would look the same as indexing (relying on out=... or adding ... to the axis seems not so nice).

And yes, this is an opinion that you can reasonably disagree on and the desire for minimal change is at the core of it.

[mihaimaruseac]

The x and y here should both (all 3 cases) be y?

[mdhaber]

Oops; I was messing with the code and forgot to update the text. Thanks!

[sjsrey]

s/same was/same way/

[mdhaber]

Suggested change

	The new functionality will be used in much the same was as old functionality,
	The new functionality will be used in much the same way as old functionality,

[seberg]

Leaving a few comments, I didn't carefully really review most of the text yet.

[seberg]

Maybe not surprising, I don't like the Array API focus, I would much prefer to only mention it as an additional argument.
First, all of the discussion is much older than it and all of the arguments work just as well without it. (In fact, to some degree better, because scalars aren't as bad at array-api as they at numpy.)

[mdhaber]

What are "all of the arguments"?

It sounds like we have two completely different perspectives about why 0d arrays should be returned instead of scalars. I think we could include both in the NEP so that we have a convincing arguments both with and without Array API, but I am less familiar with why scalars are "bad... at numpy".

[seberg]

Do what you prefer. I would expect all of your problems while doing Array API work to not actually be related to array api adoption itself?
Having to sprinkling in asarray() was always annoying and problematic maybe I am missing the argument why it is much worse now than what it was?

So, I would prefer to start with the issues below, rather than starting with "violating Array API".

The 0-D array->scalar conversion creates difficulties when writing functions that work with N-D inputs that may have 0-D intermediates. I.e. we have to add that np.asarray() call, but even if we do things may misbehave for object arrays.
The object array problem isn't huge for downstream typically, but has no good work-around.

[mdhaber]

OK, I'll need to research this "object array problem" more. I haven't run into it, personally. I assume this is "It changes the dtype of 0d object arrays containing array-likes" from numpy#13105.

maybe I am missing the argument why it is much worse now than what it was?

Probably the fundamental reasons are the same. I've just seen more frequent reminders of the inconsistency.

• Our private array API test assertions (xp_assert_close) have spurred new arguments among SciPy maintainers about whether to return 0d arrays or scalars.
  • stats chose years ago to follow the lead of NumPy and special, and so we find lots of [()] shortly before functions return. But with the array API, all these are getting converted from something like return x[()] to the more verbose return x[()] if x.ndim == 0 else x. Also, the array API prompts us to update our documentation about return types. Having to write "array or scalar" all over the place makes me think how much simpler it would be if we could follow the array in -> array out rule.
• The array API prompts us to rewrite old functions in a vectorized way. if/else branching is often converted to logic based on boolean arrays, and sometimes we are tempted to use boolean index assignment. Scalar immutability causes issues.
• With the array API, we are now more careful about dtypes. We still support NumPy 1.26, though, so I frequently run into dtype issues due to the old scalar dtype promotion rules.

These are just a few examples off the top of my head. These particular things may not be relevant to everyone in the future, so I wouldn't mention them in the NEP. The point is that that the standard prompts code rewrites and reminds us to think about details that might have been ignored before. This inevitably exacerbates the old scalar/0d array problems.

One might say that by the time this NEP would take effect, all array API rewrites will be done, so it won't matter any more. I think the same argument can be made about any bug or undesirable API, though - once you've worked around the problem, it stops bothering you... for a while, at least, until the next time you encounted the same problem. I say we take the cue and fix it once and for all.

[seberg]

Object 0-D is the same as your mutability issue in the end. It leads to some code in NumPy wrapping np.asarray() around it, I have already fixed at least two bugs with out=... that would have been very hard to fix without and this addresses similar issues.
I.e. res = arr0d + arr0d, means that res is the scalar. And in the worst case np.asarray(res) doesn't even round-trip (e.g. because res is itself an array) if:

arr0d = np.array(None)
arr0d[()] = np.arange(10)

But of course typically, we just want to avoid the np.asarray() because it is easy to forget and may need to be sprinkled in multiple times...

For me the issues around this are the main problem I want to solve.

With the array API, we are now more careful about dtypes. We still support NumPy 1.26, though, so I frequently run into dtype issues due to the old scalar dtype promotion rules.

Yeah, although it seems unrelated? Scalars and 0-D always behaved the same in promotion, only Python scalars changed behavior.

stats chose years ago to follow the lead of NumPy and special, and so we find lots of [()] shortly before functions return. But with the array API, all these are getting converted from something like return x[()] to the more verbose return x[()] if x.ndim == 0 else x. Also, the array API prompts us to update our documentation about return types. Having to write "array or scalar" all over the place makes me think how much simpler it would be if we could follow the array in -> array out rule.

This is something where I am honestly not sure what you actually want ☹. If you don't like the arr[()] return, it seems to me that you can choose to stop using it. Is the argument that NumPy leading the way for ufuncs makes that easier?

But if NumPy preserves scalars in -> scalars out, that seems not really easier for on you, unless you are considering something where np.add(np.float64(1), 1) does return an array, but maybe np.float64(1) + 1 doesn't?

So in that case the question is what additional places you want an array return that I bracketed out in the minimal version.
But I suppose something that is short of the maximal change with the env var?

My view has always been that the feasible change to me is the one where we don't implicitly convert 0-D arrays to scalars, but arr0d[()] and arr.sum() are explicit enough to me.

One might say that by the time this NEP would take effect, all array API rewrites will be done, so it won't matter any more.

These type of things never go away, and they do have real, if rare, user impact. But of course the question is how much library author pain is OK if the alternative is for some user code to silently start doing the wrong thing.

[mdhaber]

Is the argument that NumPy leading the way for ufuncs makes that easier?

Yes. We should not have to fight NumPy's behavior to make this happen. We should be able to run the code and not have to either [()] or xp.asarray on the output of every function.

[seberg]

Yeah, but these are two distinct annoyances:

1. The need for you to have to call np.asarray() when working with N-D arrays, since suddenly facing a scalar is very annoying and bug-prone.
2. You want to remove the if arr.ndim == 0: return arr[()] that matches NumPy logic currently. To me, that one is more your choice than mine. Right now you choose to align with NumPy ufuncs so you ensure a scalar result with arr[()].
   In the future you would probably just remove it and always return arrays. That makes sense, but will still diverge from NumPy, which will return scalars for scalar only inputs.

So this NEP solves the first issue fully. The second one is still the downstream libraries problem at least in the sense that you must change your code to make it happen.
And in principle you could have made that choice earlier to say "our functions only work with and return arrays".

Now, transitioning SciPy depending on the NumPy version having implemented this NEP makes sense to me as it makes the transition clearer and cleaner.
But I suppose the fact that libraries are unlikely to fully align with NumPy makes me want this point to be clear, at least between us. Right now, I am not quite sure we are on the same page that this NEP only makes problem 2 go away in the very long term and only by diverging from NumPy by only returning arrays.

(And of course the reason why NumPy doesn't change, is that scalar + scalar -> array is very awkward unless scalars effectively don't exist, and I think the tries showed that this is a too impactful change probably.)

[mdhaber]

Right now, I am not quite sure we are on the same page that this NEP only makes problem 2 go away in the very long term and only by diverging from NumPy by only returning arrays.

I think we're on the same page.

I understand that with this change, I will generally no longer need to use asarray within code written for N-d arrays because (unless axis=None or I fully index a single element) it will no long suddenly turn into a scalar. That is 1.

I understand that removing [()] has always been a choice I could make independently of NumPy.
"transitioning SciPy depending on the NumPy version having implemented this NEP makes sense to me as it makes the transition clearer and cleaner" - yes!
"only by diverging from NumPy by only returning arrays." - I think see why you would write this, but I'm not sure it's necessarily true. IIUC, special uses the ufunc machinery from NumPy, so it would follow NumPy ufunc machinery behavior exactly, whatever that is. These will stay consistent with NumPy.
The functions that have [()] at the end now are mostly reducing functions, so they don't accept scalars or 0d arrays, and therefore there's never a reason for them to return scalars. These will stay consistent with NumPy reducing functions (except maybe that they will not accept 0d objects).
For other functions which do accept scalars and don't use ufunc machinery... I don't know right now what we'd do. In the long run, if SciPy is all-in on the array API, we might choose to only accept arrays (and therefore we'd never return scalars). If we allow scalars, we might choose to explain that this will behave identically to using np.asarray before passing in input. I understand that this would be inconsistent with NumPy, so instead, we might choose to accept scalars and treat them differently from 0d arrays. If internally we use xp.asarray on the input, we'd need to make sure that we keep track of whether the original input was a 0d array or scalar so we could return the appropriate thing at the end.

So if all that makes sense, I think we're on the same page w.r.t. 2, although I'm not certain that means we will diverge from NumPy by only returning arrays. There are a few possibilities.

[jorenham]

What are "all of the arguments"?

It would also help a lot for static typing. Because, currently, the simplest numpy-esque function is currently pretty difficult to properly type:

@overload
def noop[...](a: ndarray[Shape0D, dtype[ScalarT]]) -> ScalarT: ...
@overload
def noop[ShapeT: AtLeast1D, ...](a: ndarray[ShapeT, DTypeT]) -> ndarray[ShapeT, DTypeT]: ...
@overload
def noop(a: ArrayLike) -> TheShapeIsNotKnownSoEitherArrayOrScalar: ...

If it always returned an array, this could be simplified a lot. It would also help the users, because unknown "array-likes" will retun "array of unknown dtype" instead of "array of unknown dtype or some unknown scalar". The latter often requires users to manually cast it to the "actual" type. That's why currently the numpy-stubs will often return Any in such situations, even though Any is horribly type-unsafe.

tldr; This NEP would help a lot for static typing.

[seberg]

This is confusing, the below examples are almost right, but speaking too generally about reductions isn't?
I have been very clear that in my branch both:

• arr.sum()
• arr1d[0]`
  return scalars.

You could try more for either points (I don't really believe in it from a BC stand-point, because those are the two cases where it makes a lot of sense to have code that continues relying on scalars).

[mdhaber]

I have been very clear that in my branch both:

Must have missed that. I saw:

The simple rule is, if any input is already an array, the output should also be an array. The one extension is that axis=None for reductions should return scalars.

So I suppose there is a second exception, and the easiest way to summarize it is that indexing behavior is unchanged.

[seberg]

🤷 fair, I wanted to refer to ufuncs/functions there. It does continue with:

Since that doesn't work at all, it adds an arr.to_scalar() method, that does the same as arr[()] did previously.

In the context of doing the bigger change with NUMPY_DISLIKE_SCALARS=1

---

If you feel strong about changing indexing, we could still try it out. As you can tell, I am focused on the minimal useful change.

[mdhaber]

I see. And I was hoping to go a little further and follow the array API more closely. We're coming from different perspectives. I think it's fine now that we understand that.

[mdhaber]

As for whether we should try it out - yeah, I think there's value in seeing how much additional trouble it causes. I think the best decision would come from understanding the tradeoff between switching behaviors and downstream pain. Right now, we only know that NUMPY_DISLIKE_SCALARS=1 - which is extremely aggressive, and more so than I have any motivation to propose - looks quite a bit more painful than NUMPY_DISLIKE_SCALARS=0. You made an educated guess that leaving axis=None as a special case would be much less painful - but I don't think we have any data. Maybe the transition from the painful NUMPY_DISLIKE_SCALARS=1 is very sharp and happens mostly because scalar-scalar math suddenly results in 0d arrays or something. It would be good to know!

[seberg]

Yeah, TBH, I just don't expect arr.sum() to be a huge gain either, so I want to err on the safer side.
The next question is how much you want to try arr[()], which is an orthogonal knob (and there may be more?).

[seberg]

The last one doesn't change. In fact, indexing doesn't change at all.

[seberg]

Please replace all isscalar with isinstance(..., np.ndarray). isscalar is very awkward API and much less precise anyway.

[mdhaber]

Well, I can't replace isscalar(x) with literally isinstance(x, np.ndarray) because these would basically be opposites. Do you mean change it to not isinstance(..., np.ndarray)?
I've often found isscalar useful and less awkward than trying to acheive the same thing as isinstance(... np.ndarray) because in the past SciPy has treated lists and DataFrames and such just like ndarrays, but they are not instances of ndarray. Why is it so bad?

[seberg]

Because it has a very narrow meaning of "scalar" that is just not quite right in this context (honestly it is never quite right, but does a decent enough job for numeric types).
For example isscalar(None) is False even though clearly it's not an array and np.asarray(None)[()] round-trips, so it clearly is a scalar (unless you define it very narrowly).

So yes, of course I mean not isinstance, because I think in this context is clear and has no caveats that isscalar doesn't have a precise/correct definition.

[seberg]

I think this is OK for the hashability but needs to be expanded on the mutability. While less common, it does need to be pointed out explicitly because it can lead to incorrect results.

(Unlike hashability, which should mostly lead to a relatively clear error.)

[mdhaber]

I thought about this, but I don't really understand how existing code could run into this problem.

If you try to mutate a scalar now, it raises an error. Presumably this is not what was intended, so this will not be part of a typical code path.

Do you mean that people might be relying on mutating a scalar to raise an error, but with the change it would works as expected, and that would be undesirable?

I was assuming this would hold true:

Going from error to non-error (nan or whatever) should always be fine, the reverse is typically not fine.

[seberg]

If you try to mutate a scalar now, it raises an error.

It doesn't raise for in-place operators. I actually have problems of thinking of an example (which is good!). But a silly one:

def my_metric(x):
    x += 1
    return x

arr = np.random.random()

s = arr.sum(axis=0)
metric = my_metric(s)  # modifies s

metric_s_values.append((s, metric))

The point is, having a 0-D array where code is clearly written for scalars (using in-place operators) is bad news. I could imagine it also for a simplistic integration_step function?

Is that common? Hopefully not, and to me if arr.sum() and arr[0] keep returning scalars it seems like it shouldn't be common to suddenly pass 0-D arrays to a helper function that uses in-place arithmetic.

[mdhaber]

Ah right. I realized that at the conference and recently forgot.

[seberg]

Propose something. I would say:

1. Create a branch and use it to adapt and test downstream.
2. The final NEP acceptance will happen after this with potentially more discussion.

Or in other words, I would just not say much about the precise path of achieving step 2. That makes sense IMO, because without 1 it is even harder to judge how big of a change it is (not that 1 will be a silver bullet).

How exactly you want to write it, I don't care. You could even add a .. note:: early on that a final decision is deferred until after further testing.

[seberg]

I think this somewhat leaves out a huge part of it which is the non-array API part? The whole "just add methods" works very well for array-api numeric types, but badly for general scalars like object/string dtypes.

[mdhaber]

I think this somewhat leaves out a huge part of it which is the non-array API part

Could you be more specific about what the "non-array API part" is? From context, I think you are referring to an important argument in favor of the change, but I don't know what argument that is.

[mdhaber]

Could you reply to this part @seberg?

[seberg]

I think our discussion above covered this, but I'll try to see if there is more (it may make sense to look at NumPy bugs around this, but SciPy fixes seem actually more interesting for the NEP).
I.e. I mainly mean the fact that sometimes you have to sprinkle in np.asarray() because of the magic 0-D -> scalar conversion.

[seberg]

Minimal change, minimal change, minimal change. And maybe also why not?

If you change arr.sum() to return 0-D arrays you will create a class of places where you return 0-D arrays for, honestly, unclear reasons and you will break user code. Not as much as arr1d[0], but still...
If you change arr1d.sum(-1) I think that is much less, because who would write that in a non N-D context?!

I keep repeating myself, but the point is making the return predictable for N-D arrays, changing arr.sum(), frankly, to me smells like "I don't like scalars", and not "I want to work with N-D arrays in a predictable way".

Now, "I don't like scalars" is a valid opinion, but I think we established well enough that we won't achieve that easily, so I don't see why we should try to make inroads towards it that we don't need to fix the actual core issue (don't get random scalars when working with N-D arrays).

EDIT: To be clear, I am not in love with this rule, also because it is slightly different from the rule for indexing (where adding ... is the magic). But it seems like a rule that should have a very good hit-to-miss ratio in existing code and I suppose I also can't think of a neat rule that would look the same as indexing (relying on out=... or adding ... to the axis seems not so nice).

And yes, this is an opinion that you can reasonably disagree on and the desire for minimal change is at the core of it.

[seberg]

One worry I still have -- not sure how big it is -- is the fact that it isn't necessarily trivial for downstream libraries to follow NumPy behavior.
I.e. numpy does scalar in -> scalar out, but for downstream that is hard to impossible to copy.
Unless downstream libraries want to be -> array out always on newer NumPy versions, they will API exposed and use that API likely.

[mdhaber]

but for downstream that is hard to impossible to copy.

I don't see how it is impossible. Even a decorator could be written to detect input type (scalar or array) and ensure that the output types are consistent. It might require an explicit conversion, but aside from the tens to hundreds of nanoseconds, a user would not notice the difference.

Unless downstream libraries want to be -> array out always on newer NumPy versions

They might. This is something SciPy will certainly be discussing for SciPy 2.0.

[seberg]

I don't see how it is impossible. Even a decorator could be written to detect input type (scalar or array) and ensure that the output types are consistent

Yes, but we will need public API from NumPy to do this correctly, and that machinery will again not fit well into your Array API code.

[jorenham]

This suggests that numpy scalars are not compatible with what the array-api defines as an "array object". But, is that really the case? Because on most situations, you can get away with treating a scalar as if it's a full-fledged array, even if you index it as one:

>>> duck = np.array(1)
>>> ducknt = np.int_(-1)
>>> duck[()]
np.int64(1)
>>> ducknt[()]
np.int64(-1)

I noticed that it also has an __array_namespace__. But maybe I'm missing something?

[mdhaber]

I think most of this discussion has already played out in scipy/scipy#22947 (comment). So yes, we can play devil's advocate and say that ndarray and generic work together to implement a weird array protocol where arrays might be immutable when they are 0d. But the standard itself calls out NumPy explicitly in several places for producing scalars instead of arrays, so with the standard as written, it's not a very fun game.

[jorenham]

So array-api objects are mutable, and numpy scalars are no. Got it.

[seberg]

Hehe, even trickier. Array API is either mutable or unmutable. But the abstract class of numpy array|scalar would have undefined mutability.

[jorenham]

What are "all of the arguments"?

It would also help a lot for static typing. Because, currently, the simplest numpy-esque function is currently pretty difficult to properly type:

@overload
def noop[...](a: ndarray[Shape0D, dtype[ScalarT]]) -> ScalarT: ...
@overload
def noop[ShapeT: AtLeast1D, ...](a: ndarray[ShapeT, DTypeT]) -> ndarray[ShapeT, DTypeT]: ...
@overload
def noop(a: ArrayLike) -> TheShapeIsNotKnownSoEitherArrayOrScalar: ...

If it always returned an array, this could be simplified a lot. It would also help the users, because unknown "array-likes" will retun "array of unknown dtype" instead of "array of unknown dtype or some unknown scalar". The latter often requires users to manually cast it to the "actual" type. That's why currently the numpy-stubs will often return Any in such situations, even though Any is horribly type-unsafe.

tldr; This NEP would help a lot for static typing.

[jorenham]

[nitpick] The term "rank" (in the "no. of dimensions" sense) isn't used very often in the numpy docs. Maybe because in linear algebra (and in the UK) "rank" has a different meaning?

Either way, I suppose it would be more numpy-esque to e.g. refer to "zero-rank" as "0-D", "zero-dimensional", or casual "0d".

[mdhaber]

NEP 27 proposed the addition of "Zero-rank arrays" to NumPy. That seemed like a fairly appropriate source to borrow the term from.

So it's fine if the majority of authors would prefer a different term, but I don't immediately see a reason to change it.

[jorenham]

I would've said the same thing if I would've been around when reviewing NEP 27 🤷🏻.

[mdhaber]

It would also help a lot for static typing.

Would you draft a paragraph about that for inclusion in the document?

[jorenham]

Would you draft a paragraph about that for inclusion in the document?

It would also help a lot for static typing.

Would you draft a paragraph about that for inclusion in the document?

Once the dust has settled down a bit here (e.g. when the comments are resolved), then sure.

[mhvk]

Only had a glance through this, but what seems missing most is the discussion of str and object dtype. I think those are important, though, as for those, unlike for numbers, the difference between a scalar and a 0-D array are large. (This slightly echoes @seberg's point about focussing also explicitly on numpy problems -- the array API does not deal with anything but numbers...)

I also think it would help to try to be as explicit as possible about the specific problems, as well as about the pieces of numpy where currently scalars are returned but which should be changed. Is it just,

1. Ufuncs involving 0-D arrays.
2. Reductions with axis=None
3. Indexing.
4. Meaning of indexing with an empty tuple.

Of course, there are also implicit effects on other functions if one changes these.

[mdhaber]

@mhvk thanks for your review. Would you be willing to submit a PR against this branch, or even push changes directly here? I wanted to facilitate this effort, but you and Sebastian have a much better idea of what will make a compelling argument, so I defer to you in that regard.

To answer the question about what is changing, let me modify your list slightly:

1. Ufuncs involving 0-D arrays.
2. Reductions with scalar/tuple axis.
3. Reductions with axis=None.
4. Indexing with anything but an empty tuple.
5. Indexing with an empty tuple.

Sebastian has proposed that 1 and 2 would change to return 0d arrays instead of scalars, and those are the baseline changes in the NEP. I'd like to see 3 and 4 changed, too, but we need to understand the impact of those two changes a bit better before we decide what to write about them in the NEP. I don't think anyone has pushed to change 5 since it is frequently used with the explicit purpose of converting 0d arrays to scalars. Tough to change that now.

[mdhaber]

@mhvk are you able to make those changes so this is more convincing to you and Sebastian?

[mhvk]

@mdhaber - yes, probably, though am currently on holidays so it won't be very soon...

[mdhaber]

No problem @mhvk. I'll ping you in about a month as a friendly reminder, if you'll be back then?