Python array API support #1683

Emvlt · 2025-05-26T15:02:59Z

The aim of this pull request is to make ODL multi backend through the Python Array API.

leftaroundabout · 2025-05-28T12:35:42Z

Would it be a lot of work to put the weighting hierarchy in its own PR? I think there are a bunch of questions to be discussed there, and it should be checked thoroughly whether this really works, independent of the Python Array API generalizations.

What certainly needs clarification (if not changes) is that about the methods being "mutually exclusive". It makes some sense from the definition side (we don't want conflicting semantics between e.g. distance and norm). But surely, from the use side we should then have all of them available, as far as mathematically possible? In the most typical use case (say, $L^2$) they can all be derived from the inner product, so it should be sufficient to provide only that, or alternatively the weights (which can then be used to implement first the inner product via array-API and then everything else from that).
But then there are also cases (Banach spaces) where there really is no inner product, and we'd need either weights and and exponent, or a custom norm, which also gives rise to a distance. In principle we could have only a distance and nothing else – a metric space – but I'm not sure if that would ever crop up, given that the TensorSpace class already requires a vector-space structure.

Emvlt · 2025-05-28T12:55:41Z

About the fact that the weighting refactoring should be in its own PR. I chose to add it here as for me, the python array API allows decoupling the backend classes (which will only add a few attributes to the Weighting) and the abstract class. As i am doing it for the TensorSpace and TensorSpaceElement classes i thought it would fit :)
As for the behaviour, I am not adding anything: the classes are refactored but the functionality remain unchanged. I'm happy to open the discussion on the Banach spaces and agree to create a separate PR for that :)

leftaroundabout · 2025-05-28T13:18:04Z

I did not mean to click "ready for review"

leftaroundabout · 2025-05-28T13:19:40Z

the classes are refactored but the functionality remain unchanged

Well, for one thing, I don't think the current state actually works correctly for a simple inner-product space. In the inner branch of the initialization, you're only overriding the __inner attribute, but __norm stays at its default value. Thus calling norm on that weighting will go through the default unit weight and array-norm, which will in general be inconsistent with the custom inner product.

leftaroundabout · 2025-05-28T13:30:00Z

Pretty sure all that could be fixed, but the thing is, it's not so clear-cut to say the functionality remains unchanged (and whether that even makes sense) given that the decision logic is set up in a different way now, with a dict-of-methods instead of inheritance. And this would be best investigated in a separate PR.

Emvlt · 2025-05-28T14:14:05Z

Okay, got you :) I should have checked the behaviour better.

However i don't really understand the point about splitting the PRs. I looked at your pytorch backend PR: the commits are mixed and i think it's still readable. I am also not sure that you ran the tests that would have spotted the errors that we all do while coding. I told you this is a WIP, you told me to make a PR and now you say that it's not thoroughly checked.

I guess we should align in terms of how we push changes to this repo. How about discussing live on Friday?

leftaroundabout · 2025-05-28T15:11:47Z

No critique meant. We're moving in the right direction, I just keep underestimating the amount of individual changes and possible implications. The more we separate concerns, the better we can ensure that each change is safe and won't cause more problems further down the road than it fixes.

leftaroundabout · 2025-06-19T11:41:29Z

odl/space/weighting.py

    if is_real_dtype(x2.dtype):
        return np.vecdot(x1.ravel(), x2.ravel())
    else:
        # This could also be done with `np.vdot`, which has complex conjugation


This comment is now out of date. Fixed in d9575e8

leftaroundabout · 2025-06-27T09:24:51Z

odl/test/space/space_utils_test.py

-    assert all_equal(x, ['1', '2', 'inf'])
+    with pytest.raises(AssertionError):
+       x = vector(inp)
+       # assert isinstance(x, NumpyTensor)


I don't think it makes sense for us to support string dtypes at all. Leaving those commented assertions makes it unclear what's the intention.

leftaroundabout · 2025-07-07T09:23:18Z

odl/space/entry_points.py

+TENSOR_SPACE_IMPLS = {
+    'numpy': NumpyTensorSpace
+    }
+AVAILABLE_DEVICES = {


This shouldn't really be needed now that available_devices is in ArrayBackend, right?

I think it would be better if everything that needs to know what devices are available for an impl goes through the corresponding ArrayBackend object instead of fiddling with another global dictionary (which would have to be kept in sync when new backends are added).

I agree that the devices must be only accessed through the backend. I used this variable to build the IMPL_DEVICE_PAIRS fixture as we loaded the backends but that's not clear, actually. I will move that in the pytest_config instead :)

leftaroundabout · 2025-07-07T14:14:28Z

odl/trafos/util/ft_utils.py

-    is_complex_floating_dtype, is_numeric_dtype, is_real_dtype,
-    is_real_floating_dtype, is_string, normalized_axes_tuple,
+    is_complex_dtype, is_numeric_dtype, is_real_dtype,
+    is_floating_dtype, is_string, normalized_axes_tuple,


I'm not sure I like this change of naming. After all, complex numbers (at least the standard complex64 and complex128 are also based on floating point, so it's somewhat misleading to exclude them in something called is_floating_dtype.

Why not just leave those names as they are?

Two things to make you reconsider :-)

as per NumPy's documentation, "There are 5 basic numerical types representing booleans (bool), integers (int), unsigned integers (uint) floating point (float) and complex." https://numpy.org/doc/stable/user/basics.types.html in Numerical data types.

Removing the check for complex inside the is_floating_type did not affect the tests.

as per NumPy's documentation, "There are 5 basic numerical types..."

but those are concrete fixed types (up to bit width). The is_XYZ_dtype don't refer to one particular type-flavour, but each discerns a whole bunch of types. And in the same way both int and float are numeric types (meaning for practical purposes that they support + and * and -), both float and complex are floating types, in the sense that they support / and sqrt etc.. So it might make sense to have also an is_floating_dtype which would include both float and complex, as there are many things that should work for both of those. But at any rate we should still have is_real_floating_dtype that specifically only has float32 and float64.

Okay with that!

Making sure that we are consistent with the use of the ArrayBackend class

leftaroundabout · 2025-07-10T13:19:50Z

Deprecating the test that relied on Numpy mechanisms to perform arithmetic operations between product spaces and lists. In my opinion, this should not be supported. It seems that nowhere in the codebase this behaviour is relied on. What does @leftaroundabout thinks about it?

I suppose you mean test_unary_ops and test_operators? (BTW commenting out a test case isn't really "deprecating" it)

Well, we certainly want -x and x + y and so on to work when x and y are elements of a pspace. This has nothing to do with NumPy, and it should be tested for both NumPy and PyTorch.

Emvlt · 2025-07-10T13:26:12Z

Sure, deprecating is not the right word, you know what i mean though.

perform arithmetic operations between product spaces and lists so not when x and y are elements of a pspace.

leftaroundabout · 2025-07-10T13:34:38Z

perform arithmetic operations between product spaces and lists so not when x and y are elements of a pspace.

Ah, like x + [4,5,6] for x∈ℝ³×ℝ³×ℝ³? Where exactly is/was that tested?

I'm not sure if it should be supported. I tend to say no, it would be better to require it to be written x + my_pspace.element([4,5,6]).

leftaroundabout · 2025-07-10T15:52:55Z

odl/core/array_API_support/array_creation.py

+from .utils import get_array_and_backend
+from numbers import Number
+import numpy as np
+


Do we really need this module? What is the intended use case? What are the semantics of the functions, and are they sensible?

I feel like there is kind of a wheel that keeps being reinvented in ODL, first with the various NumPy-based tricks, then with ArrayOnBackendManager, then access to the Array API, finally the ArrayBackend class.

As soon as you have an array_namespace, you can just look up the functions from that. That's certainly not something the user should need to do for basic elementwise stuff etc., but on the other hand functions such as empty_like are pretty low-level (in the ODL perspective). When on that level, how much do we really win by using functions from yet another module, with yet another semantics to get used to, as opposed to looking into the appropriate array_namespace directly at the use site?

Let's consider x. We don't know if it's an np.ndarray, a torch.Tensor or a LinearSpaceElement.
We need to have the following machinery:

import odl x_arr, backend = get_array_and_backend(x) zeros_like = backend.array_namespace.zeros_like(x_arr)

With the array_creation module, we can have

import odl zeros_like = odl.zeros_like(like_x_arr)

I think that the get_array_and_backend mechanism is a bit verbose for the user (but perfectly fine for us) and, most importantly, does not provide the documentation of the function in the IDE. I know you are a not a user of such features but i personnaly find it more user friendly.

In my opinion, our users overwhelmingly think in terms of arrays (just by personal experience+seeing the commits of non-core developpers) and will expect ODL to be lenient when it comes to the functions they are used to.

Emvlt · 2025-07-11T09:14:33Z

perform arithmetic operations between product spaces and lists so not when x and y are elements of a pspace.

Ah, like x + [4,5,6] for x∈ℝ³×ℝ³×ℝ³? Where exactly is/was that tested?

I'm not sure if it should be supported. I tend to say no, it would be better to require it to be written x + my_pspace.element([4,5,6]).
See commit 2e59ad3 :)

leftaroundabout · 2025-07-11T13:48:27Z

See commit 2e59ad3 :)

Ah, now I get it! Yes, y_arr = [op(x_) for x_ in x_arr] is much better than the old version.

…ndings.

…rward and backward calls.I found that having two functions to maintain was more error prone than having just one function with if/else logic. Also, I made sure that the cpu calls are compatible with Pytorch by doing explicit cpu conversion.

…spaces are implemented using numpy.

This remains an ongoing TODO.

Instead of having a np.float that isinstace(x, float) evaluates to True, we explicitely convert the step to a float

…ic files are isolated from the core.

…-conversion policy.

…calar. In numpy-2.0, indexing into an array does not give a plain Python number but instead e.g. `np.float64`, which is however still an `isinstance` of `float`. This situation is encountered in some of the ODL solvers.

Changed after directory-structure reorganization.

The only real problem here was that `order` is not supported any more (as it was NumPy-specific, not available in the Array API.

…d dtype. This was problematic particularly for nested product spaces, and caused many failures in the large-scale tests.

…lement of the desired space. This avoids some complications / copying and also ensures a no-op element-generation retains the identity. Change of the identity caused one of the large-scale tests to fail.

Specifically the non-support of calling NumPy ufuncs on ODL objects.

…nder the .

… different PyTorch versions. E.g. 2.7 lacks the `device` and `copy` arguments completely, wheras 2.9 refuses to handle inputs that do not live in the current Cuda device. The version distinction is an ugly hack, but at least it does not rely on exception catching (which is even more unreliable).

It is actually quite worrying that this test _succeeded_ before f67f99b. It seems like the way DLPack-transfer was implemented then caused some elements to come out in NumPy despite `pytorch` being selected as the `impl`.

This version does not use DLPack at all but only handles the relevant NumPy and PyTorch cases manually.

This is slow, but that is already the case for other scenarios covered by the function.

This function is already very forgiving with respect to different types of both input arguments, but there were still some corner cases where it errored.

Emvlt force-pushed the python_array_api_support branch from ab70b6f to c698bfa Compare May 26, 2025 15:09

leftaroundabout marked this pull request as draft May 27, 2025 13:21

leftaroundabout self-assigned this May 27, 2025

leftaroundabout added the type: api change label May 27, 2025

leftaroundabout marked this pull request as ready for review May 28, 2025 13:16

leftaroundabout marked this pull request as draft May 28, 2025 13:18

leftaroundabout mentioned this pull request Jun 10, 2025

Remove the (NumPy-) backend-specific hierarchy of weighting classes #1686

Draft

leftaroundabout reviewed Jun 19, 2025

View reviewed changes

leftaroundabout reviewed Jun 27, 2025

View reviewed changes

leftaroundabout reviewed Jul 7, 2025

View reviewed changes

Emvlt added a commit to Emvlt/odl that referenced this pull request Jul 7, 2025

Change related to odlgroup#1683 (review).

d88cacf

Making sure that we are consistent with the use of the ArrayBackend class

leftaroundabout reviewed Jul 10, 2025

View reviewed changes

Emvlt added 5 commits July 16, 2025 11:57

Tomo test suite compliant with Python-Array API and improved ASTRA bi…

f645a3c

…ndings.

Minor modifications to the skimage_radon backend to ensure that the t…

bcf9fdc

…spaces are implemented using numpy.

Cleanup of the module: removing unused function.

d94d024

Better handling of the error thrown when applying a reduction.

bf73d6f

This remains an ongoing TODO.

leftaroundabout and others added 30 commits October 20, 2025 13:55

Reflect interface change of sampling_function in the doctests.

fb6fadc

Fixing doctests to abide to the new API

b974388

Workaround for duplicate module initialization when running in PyTest.

e146e8c

Minor change on the step dtype in the derivatives module.

bd47b62

Instead of having a np.float that isinstace(x, float) evaluates to True, we explicitely convert the step to a float

Adding a self. in front of the device variable of the to_impl call

44a7258

Creation of a folder for ODL. This makes sure that the backend-specif…

18de1c7

…ic files are isolated from the core.

Make the wavelet operators compatible with ODL-1.0s no-implicit-NumPy…

05bbb9c

…-conversion policy.

Correct import in one of the tests.

08478e7

Changed after directory-structure reorganization.

Re-enable a test on DiscrSpace construction from plain lists.

f3dcfc6

The only real problem here was that `order` is not supported any more (as it was NumPy-specific, not available in the Array API.

Start a new version of the PyTorch-module wrapper for ODL operators.

73c99df

Ongoing rewriting of the torch Operator test

bd721ab

Make the inner product on product spaces independent of a pre-selecte…

c32e8d1

…d dtype. This was problematic particularly for nested product spaces, and caused many failures in the large-scale tests.

Add a forgotten import symbol.

09bc00d

Make equality checks in large-scale test NumPy-independent.

f6b7b13

Short-cut element generation when given something that is already e…

5e75af2

…lement of the desired space. This avoids some complications / copying and also ensures a no-op element-generation retains the identity. Change of the identity caused one of the large-scale tests to fail.

Remove an unnecessary check that prevented boolean operations.

04da9ac

Bring large-scale default-functionals tests in line with ODL-1.0.

a05a39e

Specifically the non-support of calling NumPy ufuncs on ODL objects.

Fix in a largescale test, the ASTRA_VERSION constant is now visible u…

f67f99b

…nder the .

Undo accidental removal of the PyTorch dependency in 18de1c7.

8304df7

Add a version-specific wrapper for from_dlpack for torch-2.7.

2b4e878

Make an assertion conditional that only makes sense on NumPy.

1e72115

It is actually quite worrying that this test _succeeded_ before f67f99b. It seems like the way DLPack-transfer was implemented then caused some elements to come out in NumPy despite `pytorch` being selected as the `impl`.

More simplistic workaround for PyTorch-DLPack inconsistencies.

536c926

This version does not use DLPack at all but only handles the relevant NumPy and PyTorch cases manually.

Allow all_equal to directly work on pairs of PyTorch tensors.

a9be53b

This is slow, but that is already the case for other scenarios covered by the function.

Enable different storage backends/devices in slow ray-trafo tests.

e707c12

Generalize a dtype case distinction beyond NumPy.

1e61123

Typo in the name of a test.

cd3202c

Further relax type restrictions on all_almost_equal.

a074314

This function is already very forgiving with respect to different types of both input arguments, but there were still some corner cases where it errored.

Generalize dtype handling in slow Tensor Space tests.

ccc312d

Python array API support #1683

Are you sure you want to change the base?

Python array API support #1683

Uh oh!

Conversation

Emvlt commented May 26, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

leftaroundabout commented May 28, 2025

Uh oh!

Emvlt commented May 28, 2025

Uh oh!

leftaroundabout commented May 28, 2025

Uh oh!

leftaroundabout commented May 28, 2025

Uh oh!

leftaroundabout commented May 28, 2025

Uh oh!

Emvlt commented May 28, 2025

Uh oh!

leftaroundabout commented May 28, 2025

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

leftaroundabout Jul 10, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

leftaroundabout commented Jul 10, 2025

Uh oh!

Emvlt commented Jul 10, 2025 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

leftaroundabout commented Jul 10, 2025

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Emvlt commented Jul 11, 2025

Uh oh!

leftaroundabout commented Jul 11, 2025

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

Emvlt commented May 26, 2025 •

edited

Loading

leftaroundabout Jul 10, 2025 •

edited

Loading

Emvlt commented Jul 10, 2025 •

edited

Loading