BUG: polyval with builtin python types for x follows NEP50

[luxedo]

Fixes #26484

The error occurs because NX.asanyarray and NX.zeros_like that casts python's float to np.float64 instead of np.float16. This fix adds tests for builtin Python types to avoid those casts.

[eendebakpt]

@luxedo Could you also add test cases for the input types that gave incorrect types? They can be added to

https://github.com/numpy/numpy/blob/main/numpy/lib/tests/test_polynomial.py

[eendebakpt]

This one is tricky to get right. What about v = np.polyval(np.array([], dtype=np.float32), 5.0)? With this PR it returns 0 with type int, but before it returned a numpy scalar.

[luxedo]

Indeed! Should it return the same dtype as the input array np.float32(0.0)? I added a suggestion to test for an empty array and return early.

There's still another case that it fails, is when passing python iterables:

>>> np.polyval(np.array([1], dtype=np.float32), [5, 1]).dtype
dtype('float64')

Is there any standardized way of dealing with this?

[eendebakpt]

Indeed! Should it return the same dtype as the input array np.float32(0.0)? I added a suggestion to test for an empty array and return early.

There's still another case that it fails, is when passing python iterables:

>>> np.polyval(np.array([1], dtype=np.float32), [5, 1]).dtype
dtype('float64')

Is there any standardized way of dealing with this?

NEP 50 is about promotion of scalars, so I think the type float64 is ok. For example:

np.array([1,2], dtype=np.float32)+[1,2] # dtype: float64

[luxedo]

Thank you for those hints! I just updated the PR

[seberg]

Having special paths for size == 0 has a code smell, and not worth optimizing for (unless there is a good reason).

To me this looks wrong for two reasons:

1. It matches any shape, even high dimensional arrays, and flattens them.
2. More importantly, why is int64 even the correct result?

Am I missing a reason for why this is needed at all?

---

EDIT: I'll have to concentrate on the rest a bit some time. The case of empty p may be a bit tricky to deal with indeed.

[eendebakpt]

The return type is not int64, as there is the conversion using astype(p.dtype). But you are correct the handing of shapes is not correct. The following case fails:

v = polyval(np.array([], dtype=np.float32), [5.0, 2.0])

The argument p could be higher dimensional, but that gives the same results (there is no broadcasting)

@luxedo You could try to replace

    return NX.int64(0).astype(p.dtype)

with

    return np.zeros_like(x).astype(p.dtype)

and add the example as a test case?

[seberg]

Thanks, I'll have to think about it. To some degree it would be nice if we could just do the identical "operation" in all cases. Like using (x * y + pv[:0]) to get the shape/dtype of the initial array and only then do the loop.