It seems I touched something in the float -> complex astype. Running

$ spin bench --compare -t bench_ufunc.NDArrayAsType.time_astype

gives me

It seems I touched something in the float -> complex astype.

Maybe the compiler just didn't decide to lift the branch for some reason? I have to say that complex -> real casts, wouldn't be my biggest concern, there is this weird ComplexWarning on it (in some branches at least), for a reason...

The linting failures are from #29197 and PR to fix them in #29210.

The OpenSUSE failure seems to be due to network failures downloading a package.

Somehow this caused test_nditer.py::test_iter_copy_casts to fail for [e-f], [e-d], [e-F], [e-D], [e-g], [F-e], [D-e], [G-e], but only on the debug build, which would suggest I touched something in the half casting.

Edit: it was due to the debug build exposing a missing initialization of context.flags

CI is passing

Hmm. I wonder what we should do with casting float-to-int or float-to-smaller-float? I think the second example should raise, but what about the third?

>>> np.array([1.0, 2.0, 100.0]).astype('int64', casting='same_value')
array([  1,   2, 100])
>>> np.array([1.2, 2.45, 100.0]).astype('int64', casting='same_value')
array([  1,   2, 100])

>>> a = np.array([1.2, 2.45, 3.14156])
>>> a == a.astype(np.half, casting='same_value')
[False, False, False]

I had a list of 5-6 different possible definitions, but for now skipping it. I am tempted to go with your a == a.astype(new_dtype, casting='same_value') must be true definition (one caveat below).
(This implies that a.astype(new_dtype, casting='same_value').astype(a.dtype) round-trips, but is stronger due to == indicating that the two dtype must be promotable int he context of comparing.)

Unless we go all out and say that a same_value cast must also be a same_kind one and disallow float to integer casts entirely here, the float must clearly be an integer.
That allows the use-case of wanting to use a float array as an integer one safely. If this is not desired, the best solution may be an "is_integral" check for the dtype.

• One additional subtlety: The float value 2e60 lacks the mantissa to exactly represent an integral value. Even though a == a.astype(new_dtype, casting='same_value') is true one could consider it unsafe from this perspective.

It may actually be worth bringing it up on the mailing list and meeting and sketching at least the main alternatives. (I don't usually expect much input, but sketching the alternatives tends to clarify things a bit either way.)

I'll try to look over the recent changes here soon, but ping me here or on slack if I forget!

I think it makes sense to do the most conservative thing and require "ability to accurately round-trip" for 'same_value' in this first implementation.

I think it makes sense to do the most conservative thing and require "ability to accurately round-trip" for 'same_value' in this first implementation.

I suppose we have at least a brief window where we could still broaden the definition slightly.
There is one use-case where a different definition might be nice, and that is if we were to consider moving towards this cast mode as a default one output parameters/assignments:

arr += something
arr[...] = other_array

That is currently same-kind I think. In that context I think it is best to forgive all float inaccuracies (including over and underflows).

But, I don't want to get lost too much in this. Unless we are sure that we are only interested in one of these two use-cases for floats (or float to int), it seems fine to just focus on one now and add the other if we want to later.

(I suppose I wouldn't mind asking around if there is a preference about the insufficient mantissa for float -> int casts, but it is a detail. -- but more of a sanity check because I don't have a strong gut feeling.)

I ran a benchmark where I disabled the internal dispatching in each call to the function, and the results showed that this is the cause of the slowdown:

 static GCC_CAST_OPT_LEVEL int
 @prefix@_cast_@name1@_to_@name2@(
         PyArrayMethod_Context *context, char *const *args,
         const npy_intp *dimensions, const npy_intp *strides,
         NpyAuxData *data)
 {
-#if !@is_bool2@
+#if 0 && !@is_bool2@
     int same_value_casting = ((context->flags & NPY_SAME_VALUE_CONTEXT_FLAG) == NPY_SAME_VALUE_CONTEXT_FLAG);
     if (same_value_casting) {
         return @prefix@_cast_@name1@_to_@name2@_same_value(context, args, dimensions, strides, data);
     } else {
 #else
     {
 #endif
         return @prefix@_cast_@name1@_to_@name2@_no_same_value(context, args, dimensions, strides, data);
 }}

So that would indicate that a path to restore the degradation in performance would be to avoid the same_value flag check when casting is safe.

There is still a slowdown, even when disabling the same_value check and dispatching. Playing with a version of the code in the compiler explorer, I think the NPY_GCC_OPT_3 attribute is insufficient to efficiently inline things to get rid of function calls. I need to decorate the half conversion utilities as well, and also inline the To BitCast(const From &from)

I also noticed the benchmarks are using itertools.combinations() and not itertools.product, so we are only benchmarking half the casting type combinations.

Edit: we cannot inline things like the float16 conversion routines, since they come from libnpymath.a which we link to.

Even if I replace lowlevel_strided_loops.c.src with the one from main, I still see a slowdown. But if I also remove the flags field from PyArrayMethod_Context, I get no slowdown. Could the size of that struct be influencing runtime?

Running benchmarks locally on a Ryzen AMD machine (the AsType benchmarks now run all possible pairings of int/float/complex dtypes) there is still some slowdowns and some speedups. I think this is acceptable?

| Change   | Before [73f8bd0f] <test_casting~1>   | After [31ee1aa4] <value-based>   |   Ratio | Benchmark (Parameter)                                            |
|----------|--------------------------------------|----------------------------------|---------|------------------------------------------------------------------|
| +        | 18.6±0.2μs                           | 21.0±0.2μs                       |    1.12 | bench_ufunc.NDArrayAsType.time_astype(('float16', 'complex64'))  |
| +        | 2.02±0.01μs                          | 2.24±0.01μs                      |    1.11 | bench_ufunc.NDArrayAsType.time_astype(('int16', 'int32'))        |
| +        | 4.47±0.04μs                          | 4.76±0.05μs                      |    1.07 | bench_ufunc.NDArrayAsType.time_astype(('complex128', 'float64')) |
| +        | 3.05±0.05μs                          | 3.24±0.05μs                      |    1.06 | bench_ufunc.NDArrayAsType.time_astype(('complex64', 'float32'))  |
| -        | 2.64±0.05μs                          | 2.51±0.05μs                      |    0.95 | bench_ufunc.NDArrayAsType.time_astype(('int64', 'int16'))        |
| -        | 24.5±0.5μs                           | 22.7±0.3μs                       |    0.93 | bench_ufunc.NDArrayAsType.time_astype(('complex128', 'float16')) |
| -        | 16.2±0.2μs                           | 15.0±0.2μs                       |    0.93 | bench_ufunc.NDArrayAsType.time_astype(('float16', 'float64'))    |
| -        | 8.06±0.07μs                          | 7.19±0.04μs                      |    0.89 | bench_ufunc.NDArrayAsType.time_astype(('complex128', 'int64'))   |
| -        | 18.7±0.1μs                           | 16.7±0.1μs                       |    0.89 | bench_ufunc.NDArrayAsType.time_astype(('float16', 'int64'))      |
| -        | 18.7±0.2μs                           | 16.3±0.3μs                       |    0.88 | bench_ufunc.NDArrayAsType.time_astype(('float16', 'int32'))      |
| -        | 18.9±0.2μs                           | 16.2±0.2μs                       |    0.86 | bench_ufunc.NDArrayAsType.time_astype(('float16', 'int16'))      |
| -        | 2.54±0.01μs                          | 2.13±0.02μs                      |    0.84 | bench_ufunc.NDArrayAsType.time_astype(('int32', 'float32'))      |

Diff from mypy_primer, showing the effect of this PR on type check results on a corpus of open source code:

spark (https://github.com/apache/spark)
- python/pyspark/ml/functions.py:244: note:     def [_ScalarT: generic[Any]] vstack(tup: Sequence[_SupportsArray[dtype[_ScalarT]] | _NestedSequence[_SupportsArray[dtype[_ScalarT]]]], *, dtype: None = ..., casting: Literal['no', 'equiv', 'safe', 'same_kind', 'unsafe'] = ...) -> ndarray[tuple[Any, ...], dtype[_ScalarT]]
+ python/pyspark/ml/functions.py:244: note:     def [_ScalarT: generic[Any]] vstack(tup: Sequence[_SupportsArray[dtype[_ScalarT]] | _NestedSequence[_SupportsArray[dtype[_ScalarT]]]], *, dtype: None = ..., casting: Literal['no', 'equiv', 'safe', 'same_kind', 'same_value', 'unsafe'] = ...) -> ndarray[tuple[Any, ...], dtype[_ScalarT]]
- python/pyspark/ml/functions.py:244: note:     def [_ScalarT: generic[Any]] vstack(tup: Sequence[Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str]], *, dtype: type[_ScalarT] | dtype[_ScalarT] | _SupportsDType[dtype[_ScalarT]], casting: Literal['no', 'equiv', 'safe', 'same_kind', 'unsafe'] = ...) -> ndarray[tuple[Any, ...], dtype[_ScalarT]]
+ python/pyspark/ml/functions.py:244: note:     def [_ScalarT: generic[Any]] vstack(tup: Sequence[Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str]], *, dtype: type[_ScalarT] | dtype[_ScalarT] | _SupportsDType[dtype[_ScalarT]], casting: Literal['no', 'equiv', 'safe', 'same_kind', 'same_value', 'unsafe'] = ...) -> ndarray[tuple[Any, ...], dtype[_ScalarT]]
- python/pyspark/ml/functions.py:244: note:     def vstack(tup: Sequence[Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str]], *, dtype: type[Any] | dtype[Any] | _SupportsDType[dtype[Any]] | tuple[Any, Any] | list[Any] | _DTypeDict | str | None = ..., casting: Literal['no', 'equiv', 'safe', 'same_kind', 'unsafe'] = ...) -> ndarray[tuple[Any, ...], dtype[Any]]
+ python/pyspark/ml/functions.py:244: note:     def vstack(tup: Sequence[Buffer | _SupportsArray[dtype[Any]] | _NestedSequence[_SupportsArray[dtype[Any]]] | complex | bytes | str | _NestedSequence[complex | bytes | str]], *, dtype: type[Any] | dtype[Any] | _SupportsDType[dtype[Any]] | tuple[Any, Any] | list[Any] | _DTypeDict | str | None = ..., casting: Literal['no', 'equiv', 'safe', 'same_kind', 'same_value', 'unsafe'] = ...) -> ndarray[tuple[Any, ...], dtype[Any]]

Extending the benchmarks for all combinations of dtypes in astype: the benchmarks ran in 3m41sec instead of approximately 3m30secs. I think the extra 10 secs is worth the more extensive benchmarking.

Any more thoughts here?

Don't really have the time for another review, but I was already happy, and you addressed all my more immediate comments. Please do raise follow-up issues about things not done, such as the various casting flags and how they are used and passed on. (@seberg might this also be relevant for how we are trying to implement sorting flags -- a general mechanism to pass on flags may be useful... but maybe I'm missing the boat completely here, don't really have the time to look properly now...)

a general mechanism to pass on flags may be useful... but maybe I'm missing the boat completely here, don't really have the time to look properly now

I think it is almost the same! But here, it seems right -- or at least OK -- to apply this flag to (almost) all ArrayMethods in the future (even ufuncs in principle, that right now are allowed to absorb casts even if there is no such ufunc inside or outside of NumPy, I am sure).

While for sorts it could look basically the same, but the flag and flag space must be specific to the sort gufunc, I think.